MercuryDPM  0.10
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Private Member Functions | Private Attributes | Friends | List of all members
NORMALIZED_POLYNOMIAL< T > Class Template Reference

This class is used to define polynomial axisymmetric coarse-graining functions. More...

#include <Polynomial.h>

Public Member Functions

 NORMALIZED_POLYNOMIAL ()
 Basic constructor; note that this does not determine the particular polynomial; one needs to call set_polynomial to define the coefficients. More...
 
void set_polynomial (std::vector< Mdouble > new_coefficients, unsigned int new_dim)
 Use this function to set the polynomial coefficients $c_i$. More...
 
void set_polynomial (Mdouble *new_coefficients, unsigned int num_coeff, unsigned int new_dim)
 Some as set_polynomial, but avoids the use of a vector. More...
 
void set_name (const char *new_name)
 Use this function to change the name of the polynomial. More...
 
std::string get_name ()
 Returns name of the polynomial. More...
 
Mdouble evaluate (Mdouble r)
 Returns the value of the polynomial, $p(r)=\sum_{i=0}^N c_i r^{N-i}$. More...
 
Mdouble evaluateGradient (Mdouble r)
 Returns the gradient of the polynomial, $\partial_\alpha p(x,y,z)=\sum_{i=0}^N c_{i,\alpha} r^{N-i},\ \alpha=x,y,z$. More...
 
Mdouble evaluateIntegral (Mdouble a, Mdouble b, Mdouble t)
 Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric function. More...
 
int get_Order (void)
 Returns the order of the polynomial. More...
 
template<>
double evaluate (double r)
 
template<>
double evaluate (double r)
 
template<>
double evaluate (double r)
 
template<>
double evaluate (double r)
 
template<>
double evaluate (double r)
 
template<>
double evaluate (double r)
 
template<>
double evaluateGradient (double r)
 
template<>
double evaluateGradient (double r)
 
template<>
double evaluateGradient (double r)
 
template<>
double evaluateGradient (double r)
 
template<>
double evaluateGradient (double r)
 
template<>
double evaluateGradient (double r)
 
template<>
double evaluateIntegral (double a, double b, double t)
 
template<>
double evaluateIntegral (double a, double b, double t)
 
template<>
double evaluateIntegral (double a, double b, double t)
 
template<>
double evaluateIntegral (double a, double b, double t)
 
template<>
double evaluateIntegral (double a, double b, double t)
 
template<>
double evaluateIntegral (double a, double b, double t)
 

Private Member Functions

void finish_set_polynomial ()
 Normalizes the polynomial coefficients $c_i$ such that the integral over the unit sphere of the axisymmetric function $p(r)$ is unity. More...
 
Mdouble get_volume ()
 Returns the integral over the unit sphere of the axisymmetric function $p(r)$. More...
 
Mdouble evaluate_1D (Mdouble r)
 Returns the value of the polynomial averaged over 2 dimensions. More...
 
Mdouble evaluate_2D (Mdouble r)
 Returns the value of the polynomial averaged over 1 dimension. More...
 
void set_average ()
 Sets averaged_coefficients. More...
 
void set_average_1D ()
 Sets averaged_coefficients $\bar{c}_i$ for StatType=X,Y,Z such that $\sum_{i=0}^N \bar{c}_i x^{N-i} = \int\int_{|\vec{x}|\leq 1} p(|\vec{x}|) dy dz$. More...
 
void set_average_2D ()
 For StatType=XY,XZ,XZ, averaged_coefficients is not used since $\bar{p}(r)$ can be evaluated as a function of $c_i$. More...
 
Mdouble evaluateGradient_1D (Mdouble r)
 Returns the value of the gradient averaged over 2 dimensions. More...
 
Mdouble evaluateGradient_2D (Mdouble r)
 Returns the value of the gradient averaged over 1 dimensions. More...
 
Mdouble evaluateIntegral_1D (Mdouble a, Mdouble b, Mdouble t)
 Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric function averaged over 2 dimensions. More...
 
Mdouble evaluateIntegral_2D (Mdouble a, Mdouble b, Mdouble t)
 Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric function averaged over 1 dimensions. More...
 
Mdouble operator[] (int i) const
 Access to the coefficients. More...
 
template<>
void set_average ()
 
template<>
void set_average ()
 
template<>
void set_average ()
 
template<>
void set_average ()
 
template<>
void set_average ()
 
template<>
void set_average ()
 

Private Attributes

std::string name
 Contains the name of the polynomial which will be displayed as CGtype by the statistical code. More...
 
unsigned int dim
 The system dimension. More...
 
std::vector< Mdouble > coefficients
 Stores the coefficients $c_i$. More...
 
std::vector< Mdouble > averaged_coefficients
 Stores some coefficients used in evaluate and evaluateIntegral for StatTypes different from XYZ. More...
 

Friends

std::ostream & operator<< (std::ostream &os, const NORMALIZED_POLYNOMIAL &P)
 Returns a text description of the polynomial. More...
 

Detailed Description

template<StatType T>
class NORMALIZED_POLYNOMIAL< T >

This class is used to define polynomial axisymmetric coarse-graining functions.

This class stores a polynomial, $p(r)=\sum_{i=0}^N c_i r^{N-i}$, which is normalized such that the integral over the unit sphere of the axisymmetric function $p(|\vec{x}|)$ is unity.

Use set_polynomial to define the polynomial. Use evaluate to evaluate the polynomial. Use evaluateGradient to evaluate the polynomial's gradient.

Calculations can be found in src/docs/Polynomials.nb

This is used to define polynomial axisymmetric coarse-graining functions (see StatisticsVector).

Note: not everything is implemented yet: only dim=3 is working, no gradients are computed.

Definition at line 52 of file Polynomial.h.

Constructor & Destructor Documentation

template<StatType T>
NORMALIZED_POLYNOMIAL< T >::NORMALIZED_POLYNOMIAL ( )
inline

Basic constructor; note that this does not determine the particular polynomial; one needs to call set_polynomial to define the coefficients.

Definition at line 76 of file Polynomial.h.

References NORMALIZED_POLYNOMIAL< T >::coefficients, NORMALIZED_POLYNOMIAL< T >::dim, and NORMALIZED_POLYNOMIAL< T >::set_name().

76  {
77  set_name("Polynomial");
78  coefficients.resize(0);
79  dim = 0;
80  }
NORMALIZED_POLYNOMIAL::dim
unsigned int dim
The system dimension.
Definition: Polynomial.h:63
NORMALIZED_POLYNOMIAL::set_name
void set_name(const char *new_name)
Use this function to change the name of the polynomial.
Definition: Polynomial.h:90
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66

Member Function Documentation

template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::evaluate ( Mdouble  r)

Returns the value of the polynomial, $p(r)=\sum_{i=0}^N c_i r^{N-i}$.

For averaged StatType this function is templated. If averaging statistics are used, then an averaged function is stored as well; for averaging a over certain dimensions is stored as well.

For averaging over two dimensions, $(y_{max}-y_{min})\cdot (z_{max}-z_{min})\cdot \bar{p}(x)=\int_{\vec{x}\leq1} p(|\vec{x}|) dy\,dz = \sum_{i=0}^N \bar{c}_i r^{N-i}$.

For averaging over one dimensions, $(z_{max}-z_{min})\cdot \bar{p}(x,y)=\int_{\vec{x}\leq1} p(|\vec{x}|) dz = \sum_{i=0}^N \bar{c}_i r^{N-i}$.

Todo:
I originally forgot to put a return statement here but the compiler did not give me a warning. Can we give a compiler option to throw warnings if there are no return statements?

Definition at line 109 of file Polynomial.hcc.

109  {
110  unsigned int N = coefficients.size();
111  double value = coefficients[0];
112  for (unsigned int i=1; i<N; i++)
113  value = value*r+coefficients[i];
115  return value;
116 }
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<>
double NORMALIZED_POLYNOMIAL< X >::evaluate ( double  r)

Definition at line 293 of file Polynomial.hcc.

293 {return evaluate_1D(r);}
NORMALIZED_POLYNOMIAL::evaluate_1D
Mdouble evaluate_1D(Mdouble r)
Returns the value of the polynomial averaged over 2 dimensions.
Definition: Polynomial.hcc:129
template<>
double NORMALIZED_POLYNOMIAL< Y >::evaluate ( double  r)

Definition at line 294 of file Polynomial.hcc.

294 {return evaluate_1D(r);}
NORMALIZED_POLYNOMIAL::evaluate_1D
Mdouble evaluate_1D(Mdouble r)
Returns the value of the polynomial averaged over 2 dimensions.
Definition: Polynomial.hcc:129
template<>
double NORMALIZED_POLYNOMIAL< Z >::evaluate ( double  r)

Definition at line 295 of file Polynomial.hcc.

295 {return evaluate_1D(r);}
NORMALIZED_POLYNOMIAL::evaluate_1D
Mdouble evaluate_1D(Mdouble r)
Returns the value of the polynomial averaged over 2 dimensions.
Definition: Polynomial.hcc:129
template<>
double NORMALIZED_POLYNOMIAL< XY >::evaluate ( double  r)

Definition at line 296 of file Polynomial.hcc.

296 {return evaluate_2D(r);}
NORMALIZED_POLYNOMIAL::evaluate_2D
Mdouble evaluate_2D(Mdouble r)
Returns the value of the polynomial averaged over 1 dimension.
Definition: Polynomial.hcc:157
template<>
double NORMALIZED_POLYNOMIAL< XZ >::evaluate ( double  r)

Definition at line 297 of file Polynomial.hcc.

297 {return evaluate_2D(r);}
NORMALIZED_POLYNOMIAL::evaluate_2D
Mdouble evaluate_2D(Mdouble r)
Returns the value of the polynomial averaged over 1 dimension.
Definition: Polynomial.hcc:157
template<>
double NORMALIZED_POLYNOMIAL< YZ >::evaluate ( double  r)

Definition at line 298 of file Polynomial.hcc.

298 {return evaluate_2D(r);}
NORMALIZED_POLYNOMIAL::evaluate_2D
Mdouble evaluate_2D(Mdouble r)
Returns the value of the polynomial averaged over 1 dimension.
Definition: Polynomial.hcc:157
template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::evaluate_1D ( Mdouble  r)
private

Returns the value of the polynomial averaged over 2 dimensions.

For StatType=X, $r=|x|$. See also set_average_1D.

Definition at line 129 of file Polynomial.hcc.

129  {
130  unsigned int N = averaged_coefficients.size();
131  double value = averaged_coefficients[0];
132  for (unsigned int i=1; i<N; i++)
133  value = value*r+averaged_coefficients[i];
134  return value;
135 }
NORMALIZED_POLYNOMIAL::averaged_coefficients
std::vector< Mdouble > averaged_coefficients
Stores some coefficients used in evaluate and evaluateIntegral for StatTypes different from XYZ...
Definition: Polynomial.h:69
template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::evaluate_2D ( Mdouble  r)
private

Returns the value of the polynomial averaged over 1 dimension.

For StatType=XY, $r=\sqrt{x^2+y^2}$.

Definition at line 157 of file Polynomial.hcc.

157  {
158  //~ std::cout << "eval_2D" << std::endl;
159  unsigned int N = coefficients.size();
160  double value1 = 2.*coefficients[N-1];
161  double value2 = 0;
162  double r2 = r*r;
163  switch (N) {
164  case 7:
165  value1 += coefficients[N-7]*2./35.*(5.+r2*(6.+r2*(8*r2+16)));
166  case 6:
167  value1 += coefficients[N-6]*(8.+(10.+15.*r2)*r2)/24.;
168  value2 += coefficients[N-6]*(15.*r2*r2*r2)/24.;
169  case 5:
170  value1 += coefficients[N-5]*.4*(1.+r2*(4./3.+r2*8./3.));
171  case 4:
172  value1 += coefficients[N-4]*(2.+r2*3.)/4.;
173  value2 += coefficients[N-4]*.75*r2*r2;
174  case 3:
175  value1 += coefficients[N-3]*2./3.*(1.+r2*2.);
176  case 2:
177  value1 += coefficients[N-2]*1.;
178  value2 += coefficients[N-2]*r2;
179  case 1:
180  break;
181  default:
182  std::cerr << "Error: no rules set for high-order polynomials" << std::endl; exit(-1);
183  }
184  double root = sqrt(1-r2);
185  double invz = root/r;
186  double arccsch = log(sqrt(1.+invz*invz)+invz);
187  return value1*root+value2*arccsch;
188 }
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::evaluateGradient ( Mdouble  r)

Returns the gradient of the polynomial, $\partial_\alpha p(x,y,z)=\sum_{i=0}^N c_{i,\alpha} r^{N-i},\ \alpha=x,y,z$.

Definition at line 120 of file Polynomial.hcc.

120  {
121  unsigned int N = coefficients.size();
122  double value = (N-1)*coefficients[0];
123  for (unsigned int i=1; i<N-1; i++)
124  value = value*r+(N-1-i)*coefficients[i];
125  return value;
126 }
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<>
double NORMALIZED_POLYNOMIAL< X >::evaluateGradient ( double  r)

Definition at line 300 of file Polynomial.hcc.

300 {return evaluateGradient_1D(r);}
NORMALIZED_POLYNOMIAL::evaluateGradient_1D
Mdouble evaluateGradient_1D(Mdouble r)
Returns the value of the gradient averaged over 2 dimensions.
Definition: Polynomial.hcc:138
template<>
double NORMALIZED_POLYNOMIAL< Y >::evaluateGradient ( double  r)

Definition at line 301 of file Polynomial.hcc.

301 {return evaluateGradient_1D(r);}
NORMALIZED_POLYNOMIAL::evaluateGradient_1D
Mdouble evaluateGradient_1D(Mdouble r)
Returns the value of the gradient averaged over 2 dimensions.
Definition: Polynomial.hcc:138
template<>
double NORMALIZED_POLYNOMIAL< Z >::evaluateGradient ( double  r)

Definition at line 302 of file Polynomial.hcc.

302 {return evaluateGradient_1D(r);}
NORMALIZED_POLYNOMIAL::evaluateGradient_1D
Mdouble evaluateGradient_1D(Mdouble r)
Returns the value of the gradient averaged over 2 dimensions.
Definition: Polynomial.hcc:138
template<>
double NORMALIZED_POLYNOMIAL< XY >::evaluateGradient ( double  r)

Definition at line 303 of file Polynomial.hcc.

303 {return evaluateGradient_2D(r);}
NORMALIZED_POLYNOMIAL::evaluateGradient_2D
Mdouble evaluateGradient_2D(Mdouble r)
Returns the value of the gradient averaged over 1 dimensions.
Definition: Polynomial.hcc:150
template<>
double NORMALIZED_POLYNOMIAL< XZ >::evaluateGradient ( double  r)

Definition at line 304 of file Polynomial.hcc.

304 {return evaluateGradient_2D(r);}
NORMALIZED_POLYNOMIAL::evaluateGradient_2D
Mdouble evaluateGradient_2D(Mdouble r)
Returns the value of the gradient averaged over 1 dimensions.
Definition: Polynomial.hcc:150
template<>
double NORMALIZED_POLYNOMIAL< YZ >::evaluateGradient ( double  r)

Definition at line 305 of file Polynomial.hcc.

305 {return evaluateGradient_2D(r);}
NORMALIZED_POLYNOMIAL::evaluateGradient_2D
Mdouble evaluateGradient_2D(Mdouble r)
Returns the value of the gradient averaged over 1 dimensions.
Definition: Polynomial.hcc:150
template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::evaluateGradient_1D ( Mdouble  r)
private

Returns the value of the gradient averaged over 2 dimensions.

Definition at line 138 of file Polynomial.hcc.

References constants::pi.

138  {
139  unsigned int N = coefficients.size();
140  double value = coefficients[0];
141  double value2 = coefficients[0];
142  for (unsigned int i=1; i<N-1; i++) {
143  value = value*r+coefficients[i];
144  value2 += coefficients[i];
145  }
146  return 2.*constants::pi*r*(value2-value*r);
147 }
constants::pi
const Mdouble pi
Definition: ExtendedMath.h:54
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::evaluateGradient_2D ( Mdouble  r)
private

Returns the value of the gradient averaged over 1 dimensions.

Todo:

Definition at line 150 of file Polynomial.hcc.

150  {
152  return 0.0;
153 }
template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::evaluateIntegral ( Mdouble  a,
Mdouble  b,
Mdouble  t 
)

Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric function.

averaging.jpeg
circle denotes the cutoff radius of the cg function around P

For averaged StatType this function is templated.

Definition at line 191 of file Polynomial.hcc.

References A.

191  {
192  unsigned int N = coefficients.size();
193  double value1 = coefficients[N-1]*(b-a);
194  double value2 = 0;
195  double t2 = t*t;
196  double a2 = a*a;
197  double b2 = b*b;
198  double roota = sqrt(a*a+t2);
199  double rootb = sqrt(b*b+t2);
200  switch (N) {
201  case 7:
202  value1 += coefficients[N-7]*(b*(b2*b2*b2/7.0+0.6*b2*b2*t2+b2*t2*t2+t2*t2*t2)-a*(a2*a2*a2/7.0+0.6*a2*a2*t2+a2*t2*t2+t2*t2*t2));
203  case 6:
204  value1 += coefficients[N-6]*(b*rootb*(8.*b2*b2+26.*b2*t2+33.*t2*t2)-a*roota*(8.*a2*a2+26.*a2*t2+33.*t2*t2))/48.;
205  value2 += coefficients[N-6]*15./48.*t2*t2*t2;
206  case 5:
207  value1 += coefficients[N-5]*(b*(b2*b2/5.+2./3.*b2*t2+t2*t2)-a*(a2*a2/5.+2./3.*a2*t2+t2*t2));
208  case 4:
209  value1 += coefficients[N-4]*(b*rootb*(2.*b2+5.*t2)-a*roota*(2.*a2+5.*t2))/8.;
210  value2 += coefficients[N-4]*3./8.*t2*t2;
211  case 3:
212  value1 += coefficients[N-3]*(b*(b2/3.+t2)-a*(a2/3.+t2));
213  case 2:
214  value1 += coefficients[N-2]*(b*rootb-a*roota)/2.;
215  value2 += coefficients[N-2]*t2/2.;
216  case 1:
217  break;
218  default:
219  std::cerr << "Error: no rules set for high-order polynomials" << std::endl; exit(-1);
220  }
221  double A = a/t;
222  double B = b/t;
223  double arcsinh = t<1e-12 ? 0.0 : log(B+sqrt(1.+B*B))-log(A+sqrt(1.+A*A));
224  //~ std::cout
225  //~ << "value1 " << value1 << std::endl
226  //~ << "value2 " << value2 << std::endl
227  //~ << "arcsinh " << arcsinh << std::endl
228  //~ << "a " << a << std::endl
229  //~ << "b " << b << std::endl
230  //~ << "t " << t << std::endl;
231  return value1+value2*arcsinh;
232 }
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<>
double NORMALIZED_POLYNOMIAL< X >::evaluateIntegral ( double  a,
double  b,
double  t 
)

Definition at line 307 of file Polynomial.hcc.

307 {return evaluateIntegral_1D(a,b,t);}
NORMALIZED_POLYNOMIAL::evaluateIntegral_1D
Mdouble evaluateIntegral_1D(Mdouble a, Mdouble b, Mdouble t)
Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric func...
Definition: Polynomial.hcc:242
template<>
double NORMALIZED_POLYNOMIAL< Y >::evaluateIntegral ( double  a,
double  b,
double  t 
)

Definition at line 308 of file Polynomial.hcc.

308 {return evaluateIntegral_1D(a,b,t);}
NORMALIZED_POLYNOMIAL::evaluateIntegral_1D
Mdouble evaluateIntegral_1D(Mdouble a, Mdouble b, Mdouble t)
Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric func...
Definition: Polynomial.hcc:242
template<>
double NORMALIZED_POLYNOMIAL< Z >::evaluateIntegral ( double  a,
double  b,
double  t 
)

Definition at line 309 of file Polynomial.hcc.

309 {return evaluateIntegral_1D(a,b,t);}
NORMALIZED_POLYNOMIAL::evaluateIntegral_1D
Mdouble evaluateIntegral_1D(Mdouble a, Mdouble b, Mdouble t)
Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric func...
Definition: Polynomial.hcc:242
template<>
double NORMALIZED_POLYNOMIAL< XY >::evaluateIntegral ( double  a,
double  b,
double  t 
)

Definition at line 310 of file Polynomial.hcc.

310 {return evaluateIntegral_2D(a,b,t);}
NORMALIZED_POLYNOMIAL::evaluateIntegral_2D
Mdouble evaluateIntegral_2D(Mdouble a, Mdouble b, Mdouble t)
Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric func...
Definition: Polynomial.hcc:236
template<>
double NORMALIZED_POLYNOMIAL< XZ >::evaluateIntegral ( double  a,
double  b,
double  t 
)

Definition at line 311 of file Polynomial.hcc.

311 {return evaluateIntegral_2D(a,b,t);}
NORMALIZED_POLYNOMIAL::evaluateIntegral_2D
Mdouble evaluateIntegral_2D(Mdouble a, Mdouble b, Mdouble t)
Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric func...
Definition: Polynomial.hcc:236
template<>
double NORMALIZED_POLYNOMIAL< YZ >::evaluateIntegral ( double  a,
double  b,
double  t 
)

Definition at line 312 of file Polynomial.hcc.

312 {return evaluateIntegral_2D(a,b,t);}
NORMALIZED_POLYNOMIAL::evaluateIntegral_2D
Mdouble evaluateIntegral_2D(Mdouble a, Mdouble b, Mdouble t)
Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric func...
Definition: Polynomial.hcc:236
template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::evaluateIntegral_1D ( Mdouble  a,
Mdouble  b,
Mdouble  t 
)
private

Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric function averaged over 2 dimensions.

Definition at line 242 of file Polynomial.hcc.

References A.

242  {
243  unsigned int N = averaged_coefficients.size();
244  double value1 = averaged_coefficients[N-1]*(b-a);
245  double value2 = 0;
246  double t2 = t*t;
247  double a2 = a*a;
248  double b2 = b*b;
249  double roota = sqrt(a*a+t2);
250  double rootb = sqrt(b*b+t2);
251  switch (N) {
252  case 7:
253  value1 += averaged_coefficients[N-7]*(b*(b2*b2*b2/7.0+0.6*b2*b2*t2+b2*t2*t2+t2*t2*t2)-a*(a2*a2*a2/7.0+0.6*a2*a2*t2+a2*t2*t2+t2*t2*t2));
254  case 6:
255  value1 += averaged_coefficients[N-6]*(b*rootb*(8.*b2*b2+26.*b2*t2+33.*t2*t2)-a*roota*(8.*a2*a2+26.*a2*t2+33.*t2*t2))/48.;
256  value2 += averaged_coefficients[N-6]*15./48.*t2*t2*t2;
257  case 5:
258  value1 += averaged_coefficients[N-5]*(b*(b2*b2/5.+2./3.*b2*t2+t2*t2)-a*(a2*a2/5.+2./3.*a2*t2+t2*t2));
259  case 4:
260  value1 += averaged_coefficients[N-4]*(b*rootb*(2.*b2+5.*t2)-a*roota*(2.*a2+5.*t2))/8.;
261  value2 += averaged_coefficients[N-4]*3./8.*t2*t2;
262  case 3:
263  value1 += averaged_coefficients[N-3]*(b*(b2/3.+t2)-a*(a2/3.+t2));
264  case 2:
265  value1 += averaged_coefficients[N-2]*(b*rootb-a*roota)/2.;
266  value2 += averaged_coefficients[N-2]*t2/2.;
267  case 1:
268  break;
269  default:
270  std::cerr << "Error: no rules set for high-order polynomials" << std::endl; exit(-1);
271  }
272  double A = a/t;
273  double B = b/t;
274  double arcsinh = t<1e-12 ? 0.0 : log(B+sqrt(1.+B*B))-log(A+sqrt(1.+A*A));
275  return value1+value2*arcsinh;
276 
277  // 3D Heaviside: ( b-a );
278  // 3D Polynomial: evaluateIntegral(max(a,-wn)/w,min(b,wn)/w,tangential.GetLength()/w)*w;
279  // 1D Polynomial: evaluateIntegral(max(a,-wn)/w,min(b,wn)/w,tangential.GetLength()/w)*w;
280  // 1D Heaviside: ((b*pi*(wn2-sqr(b)/3.0)))-((a*pi*(wn2-sqr(a)/3.0)));
281  //~ averaged_coefficients[N+1] += coefficients[i]*2*pi/(N-1-i+2);
282  //~ averaged_coefficients[i] -= coefficients[i]*2*pi/(N-1-i+2);
283 
284 }
NORMALIZED_POLYNOMIAL::averaged_coefficients
std::vector< Mdouble > averaged_coefficients
Stores some coefficients used in evaluate and evaluateIntegral for StatTypes different from XYZ...
Definition: Polynomial.h:69
template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::evaluateIntegral_2D ( Mdouble a  UNUSED,
Mdouble b  UNUSED,
Mdouble t  UNUSED 
)
private

Returns the value of the line integral along the normal P1P2 "from a to b" over the axisymmetric function averaged over 1 dimensions.

Todo:
{Thomas: has yet to be implemented}

Definition at line 236 of file Polynomial.hcc.

236  {
237  std::cerr << "evaluateIntegral_2D is not implemented yet" << std::endl;
238  exit(-1);
239 }
template<StatType T>
void NORMALIZED_POLYNOMIAL< T >::finish_set_polynomial ( )
private

Normalizes the polynomial coefficients $c_i$ such that the integral over the unit sphere of the axisymmetric function $p(r)$ is unity.

In detail, $\int_0^1 f(r) p(r) dr = 1$, with $f(r)=4\pi r^2$ for 3D, $f(r)=2\pi r$ for 2D, $f(r)=2$ for 1D systems.

Also sets averaged_coefficients

Assumes that dim and coefficients are already set.

Definition at line 43 of file Polynomial.hcc.

43  {
44  //normalizes the Polynomial.
45  double volume = get_volume();
46  if (fabs(volume-1.0)>1e-12)
47  std::cout << "Warning: Polynomial has to be normalized by " << volume << std::endl;
48  for (unsigned int i=0; i<coefficients.size(); i++)
49  coefficients[i] /= volume;
50 
51  //sets #averaged_coefficients
52  set_average();
53 
54  //std::couts the Polynomial
55  std::stringstream sstr;
56  sstr << *this;
57  set_name(sstr.str().c_str());
58  std::cout << "p=" << *this << std::endl;
59 }
NORMALIZED_POLYNOMIAL::set_average
void set_average()
Sets averaged_coefficients.
Definition: Polynomial.hcc:86
NORMALIZED_POLYNOMIAL::set_name
void set_name(const char *new_name)
Use this function to change the name of the polynomial.
Definition: Polynomial.h:90
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
NORMALIZED_POLYNOMIAL::get_volume
Mdouble get_volume()
Returns the integral over the unit sphere of the axisymmetric function .
Definition: Polynomial.hcc:62
template<StatType T>
std::string NORMALIZED_POLYNOMIAL< T >::get_name ( )
inline

Returns name of the polynomial.

Definition at line 93 of file Polynomial.h.

References NORMALIZED_POLYNOMIAL< T >::name.

93 {return name;}
NORMALIZED_POLYNOMIAL::name
std::string name
Contains the name of the polynomial which will be displayed as CGtype by the statistical code...
Definition: Polynomial.h:60
template<StatType T>
int NORMALIZED_POLYNOMIAL< T >::get_Order ( void  )
inline

Returns the order of the polynomial.

Definition at line 117 of file Polynomial.h.

References NORMALIZED_POLYNOMIAL< T >::coefficients.

117  {
118  return coefficients.size()-1;
119  }
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<StatType T>
double NORMALIZED_POLYNOMIAL< T >::get_volume ( )
private

Returns the integral over the unit sphere of the axisymmetric function $p(r)$.

In detail, $\int_{|\vec{x}|\leq1} p(|\vec{x}|) d\vec{x} = \int_0^1 f(r) p(r) dr = 1$, with $f(r)=4\pi r^2$ for 3D, $f(r)=2\pi r$ for 2D, $f(r)=2$ for 1D systems.

For $p(r)=\sum_{i=0}^{N-1} c_i r^{N-1-i}$, we obtain $V = \sum_{i=0}^{N-1} 4\pi c_i/(2+N-i)$ for 3D, $V = \sum_{i=0}^{N-1} 2\pi c_i/(1+N-i)$ for 2D, $V = \sum_{i=0}^{N-1} 2 c_i/(N-i)$ for 1D systems.

Definition at line 62 of file Polynomial.hcc.

References constants::pi.

62  {
63  double volume = 0;
64  unsigned int N = coefficients.size();
65  if (dim==3) {
66  for (unsigned int i=0; i<N; i++)
67  volume += coefficients[i]/(2.+N-i);
68  volume *= 4.*constants::pi;
69  } else if (dim==2) {
70  std::cerr << "dim=2 is not working yet" << std::endl; exit(-1);
71  for (unsigned int i=0; i<coefficients.size(); i++)
72  volume += coefficients[i]/(1.+N-i);
73  volume *= 2.*constants::pi;
74  } else if (dim==1) {
75  std::cerr << "dim=1 is not working yet" << std::endl; exit(-1);
76  for (unsigned int i=0; i<coefficients.size(); i++)
77  volume += coefficients[i]/(0.+N-i);
78  volume *= 2.;
79  } else {
80  std::cerr << "Error in get_volume: dim=" << dim << std::endl; exit(-1);
81  }
82  return volume;
83 }
NORMALIZED_POLYNOMIAL::dim
unsigned int dim
The system dimension.
Definition: Polynomial.h:63
constants::pi
const Mdouble pi
Definition: ExtendedMath.h:54
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<StatType T>
Mdouble NORMALIZED_POLYNOMIAL< T >::operator[] ( int  i) const
inlineprivate

Access to the coefficients.

Definition at line 193 of file Polynomial.h.

References NORMALIZED_POLYNOMIAL< T >::coefficients.

193  {
194  return coefficients[i];
195  }
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<StatType T>
void NORMALIZED_POLYNOMIAL< T >::set_average ( )
private

Sets averaged_coefficients.

This function is templated, with the default used only for StatType=XYZ, so it does nothing.

Definition at line 86 of file Polynomial.hcc.

86 {std::cout << "set_average" << std::endl;}
template<>
void NORMALIZED_POLYNOMIAL< X >::set_average ( )
private

Definition at line 286 of file Polynomial.hcc.

286 {set_average_1D();}
NORMALIZED_POLYNOMIAL::set_average_1D
void set_average_1D()
Sets averaged_coefficients for StatType=X,Y,Z such that .
Definition: Polynomial.hcc:89
template<>
void NORMALIZED_POLYNOMIAL< Y >::set_average ( )
private

Definition at line 287 of file Polynomial.hcc.

287 {set_average_1D();}
NORMALIZED_POLYNOMIAL::set_average_1D
void set_average_1D()
Sets averaged_coefficients for StatType=X,Y,Z such that .
Definition: Polynomial.hcc:89
template<>
void NORMALIZED_POLYNOMIAL< Z >::set_average ( )
private

Definition at line 288 of file Polynomial.hcc.

288 {set_average_1D();}
NORMALIZED_POLYNOMIAL::set_average_1D
void set_average_1D()
Sets averaged_coefficients for StatType=X,Y,Z such that .
Definition: Polynomial.hcc:89
template<>
void NORMALIZED_POLYNOMIAL< XY >::set_average ( )
private

Definition at line 289 of file Polynomial.hcc.

289 {set_average_2D();}
NORMALIZED_POLYNOMIAL::set_average_2D
void set_average_2D()
For StatType=XY,XZ,XZ, averaged_coefficients is not used since can be evaluated as a function of ...
Definition: Polynomial.hcc:104
template<>
void NORMALIZED_POLYNOMIAL< XZ >::set_average ( )
private

Definition at line 290 of file Polynomial.hcc.

290 {set_average_2D();}
NORMALIZED_POLYNOMIAL::set_average_2D
void set_average_2D()
For StatType=XY,XZ,XZ, averaged_coefficients is not used since can be evaluated as a function of ...
Definition: Polynomial.hcc:104
template<>
void NORMALIZED_POLYNOMIAL< YZ >::set_average ( )
private

Definition at line 291 of file Polynomial.hcc.

291 {set_average_2D();}
NORMALIZED_POLYNOMIAL::set_average_2D
void set_average_2D()
For StatType=XY,XZ,XZ, averaged_coefficients is not used since can be evaluated as a function of ...
Definition: Polynomial.hcc:104
template<StatType T>
void NORMALIZED_POLYNOMIAL< T >::set_average_1D ( )
private

Sets averaged_coefficients $\bar{c}_i$ for StatType=X,Y,Z such that $\sum_{i=0}^N \bar{c}_i x^{N-i} = \int\int_{|\vec{x}|\leq 1} p(|\vec{x}|) dy dz$.

See evaluate_1D.

Definition at line 89 of file Polynomial.hcc.

References constants::pi.

89  {
90  std::cout << "set_average_1D" << std::endl;
91  unsigned int N = coefficients.size();
92  averaged_coefficients.resize(N+2);
93  for (unsigned int i=0; i<N+2; i++) {
94  averaged_coefficients[i] = 0;
95  }
96  for (unsigned int i=0; i<N; i++) {
97  averaged_coefficients[N+1] += coefficients[i]*2*constants::pi/(N-1-i+2);
98  averaged_coefficients[i] -= coefficients[i]*2*constants::pi/(N-1-i+2);
99  }
100 }
constants::pi
const Mdouble pi
Definition: ExtendedMath.h:54
NORMALIZED_POLYNOMIAL::averaged_coefficients
std::vector< Mdouble > averaged_coefficients
Stores some coefficients used in evaluate and evaluateIntegral for StatTypes different from XYZ...
Definition: Polynomial.h:69
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<StatType T>
void NORMALIZED_POLYNOMIAL< T >::set_average_2D ( )
private

For StatType=XY,XZ,XZ, averaged_coefficients is not used since $\bar{p}(r)$ can be evaluated as a function of $c_i$.

See evaluate_2D.

Definition at line 104 of file Polynomial.hcc.

104  {
105 }
template<StatType T>
void NORMALIZED_POLYNOMIAL< T >::set_name ( const char *  new_name)
inline

Use this function to change the name of the polynomial.

Definition at line 90 of file Polynomial.h.

References NORMALIZED_POLYNOMIAL< T >::name.

Referenced by NORMALIZED_POLYNOMIAL< T >::NORMALIZED_POLYNOMIAL().

90 {name = new_name;}
NORMALIZED_POLYNOMIAL::name
std::string name
Contains the name of the polynomial which will be displayed as CGtype by the statistical code...
Definition: Polynomial.h:60
template<StatType T>
void NORMALIZED_POLYNOMIAL< T >::set_polynomial ( std::vector< Mdouble >  new_coefficients,
unsigned int  new_dim 
)

Use this function to set the polynomial coefficients $c_i$.

This function calls finish_set_polynomial to normalize the coefficients.

Definition at line 27 of file Polynomial.hcc.

27  {
28  coefficients = new_coefficients;
29  dim = new_dim;
30  finish_set_polynomial();
31 }
NORMALIZED_POLYNOMIAL::dim
unsigned int dim
The system dimension.
Definition: Polynomial.h:63
NORMALIZED_POLYNOMIAL::finish_set_polynomial
void finish_set_polynomial()
Normalizes the polynomial coefficients such that the integral over the unit sphere of the axisymmetr...
Definition: Polynomial.hcc:43
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66
template<StatType T>
void NORMALIZED_POLYNOMIAL< T >::set_polynomial ( Mdouble *  new_coefficients,
unsigned int  num_coeff,
unsigned int  new_dim 
)

Some as set_polynomial, but avoids the use of a vector.

Definition at line 34 of file Polynomial.hcc.

34  {
35  coefficients.resize(num_coeff);
36  for (unsigned int i=0; i<num_coeff; i++)
37  coefficients[i] = new_coefficients[i];
38  dim = new_dim;
39  finish_set_polynomial();
40 }
NORMALIZED_POLYNOMIAL::dim
unsigned int dim
The system dimension.
Definition: Polynomial.h:63
NORMALIZED_POLYNOMIAL::finish_set_polynomial
void finish_set_polynomial()
Normalizes the polynomial coefficients such that the integral over the unit sphere of the axisymmetr...
Definition: Polynomial.hcc:43
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66

Friends And Related Function Documentation

template<StatType T>
std::ostream& operator<< ( std::ostream &  os,
const NORMALIZED_POLYNOMIAL< T > &  P 
)
friend

Returns a text description of the polynomial.

Definition at line 122 of file Polynomial.h.

122  {
123  unsigned int N = P.coefficients.size();
124  for (unsigned int i=0; i<N; i++) {
125  if (!P[i]) continue;
126  if (P[i]>=0) os << "+";
127  os << std::setprecision(2) << P[i];
128  if (N-1-i>1) os << "r^" << N-1-i;
129  else if (N-1-i==1) os << "r";
130  }
131  return os;
132  }
NORMALIZED_POLYNOMIAL::coefficients
std::vector< Mdouble > coefficients
Stores the coefficients .
Definition: Polynomial.h:66

Member Data Documentation

template<StatType T>
std::vector<Mdouble> NORMALIZED_POLYNOMIAL< T >::averaged_coefficients
private

Stores some coefficients used in evaluate and evaluateIntegral for StatTypes different from XYZ.

Definition at line 69 of file Polynomial.h.

template<StatType T>
std::vector<Mdouble> NORMALIZED_POLYNOMIAL< T >::coefficients
private

Stores the coefficients $c_i$.

Definition at line 66 of file Polynomial.h.

Referenced by NORMALIZED_POLYNOMIAL< T >::get_Order(), NORMALIZED_POLYNOMIAL< T >::NORMALIZED_POLYNOMIAL(), and NORMALIZED_POLYNOMIAL< T >::operator[]().

template<StatType T>
unsigned int NORMALIZED_POLYNOMIAL< T >::dim
private

The system dimension.

Definition at line 63 of file Polynomial.h.

Referenced by NORMALIZED_POLYNOMIAL< T >::NORMALIZED_POLYNOMIAL().

template<StatType T>
std::string NORMALIZED_POLYNOMIAL< T >::name
private

Contains the name of the polynomial which will be displayed as CGtype by the statistical code.

Definition at line 60 of file Polynomial.h.

Referenced by NORMALIZED_POLYNOMIAL< T >::get_name(), and NORMALIZED_POLYNOMIAL< T >::set_name().

The documentation for this class was generated from the following files: