www.mooseframework.org
PolynomialFit.h
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
15 #ifndef POLYNOMIALFIT_H
16 #define POLYNOMIALFIT_H
17 
18 #include <vector>
19 #include <string>
20 
21 #include "Moose.h"
22 
30 {
31 public:
32  /* Constructor, Takes two vectors of points for which to apply the fit. One should be of the
33  * independent variable while the other should be of the dependent variable. These values should
34  * correspond to one and other in the same position. The third parameter is the requested
35  * polynomial
36  * order and the forth parameter tells the class whether or not it should truncate the order if
37  * there
38  * are not enough points for which to apply the polynomial fit.
39  */
40  PolynomialFit(std::vector<Real> X,
41  std::vector<Real> Y,
42  unsigned int order,
43  bool truncate_order = false);
44 
45  virtual ~PolynomialFit() = default;
46 
52  void generate();
53 
58  Real sample(Real x);
59 
64  void dumpSampleFile(std::string base_name,
65  std::string x_label = "X",
66  std::string y_label = "Y",
67  float xmin = 0,
68  float xmax = 0,
69  float ymin = 0,
70  float ymax = 0);
71 
76  unsigned int getSampleSize();
77 
81  const std::vector<Real> & getCoefficients();
82 
83 private:
87  void fillMatrix();
88 
92  void doLeastSquares();
93 
94  std::vector<Real> _x;
95  std::vector<Real> _y;
96  std::vector<Real> _matrix;
97  std::vector<Real> _coeffs;
98  unsigned int _order;
100 
101  static int _file_number;
102 };
103 
104 #endif // POLYNOMIALFIT_H
std::vector< Real > _matrix
Definition: PolynomialFit.h:96
const std::vector< Real > & getCoefficients()
Get a const reference to the coefficients of the least squares fit.
static PetscErrorCode Vec x
PolynomialFit(std::vector< Real > X, std::vector< Real > Y, unsigned int order, bool truncate_order=false)
Definition: PolynomialFit.C:24
static int _file_number
void generate()
This function generates the polynomial fit.
Definition: PolynomialFit.C:48
unsigned int getSampleSize()
This function returns the size of the array holding the points, i.e.
virtual ~PolynomialFit()=default
std::vector< Real > _x
Definition: PolynomialFit.h:94
This class applies the Least Squares algorithm to a set of points to provide a smooth curve for sampl...
Definition: PolynomialFit.h:29
bool _truncate_order
Definition: PolynomialFit.h:99
std::vector< Real > _y
Definition: PolynomialFit.h:95
std::vector< Real > _coeffs
Definition: PolynomialFit.h:97
void doLeastSquares()
This function is the wrapper for the LAPACK dgels function and is called by generate.
Definition: PolynomialFit.C:75
void fillMatrix()
This is a helper function that creates the matrix necessary for the Least Squares algorithm...
Definition: PolynomialFit.C:55
Real sample(Real x)
This function will take an independent variable input and will return the dependent variable based on...
void dumpSampleFile(std::string base_name, std::string x_label="X", std::string y_label="Y", float xmin=0, float xmax=0, float ymin=0, float ymax=0)
This function will dump GNUPLOT input files that can be run to show the data points and function fits...
unsigned int _order
Definition: PolynomialFit.h:98