Go to the documentation of this file.
23 #ifndef O2SCL_DERIV_EQI_H
24 #define O2SCL_DERIV_EQI_H
32 #include <boost/numeric/ublas/vector.hpp>
34 #include <o2scl/deriv.h>
35 #include <o2scl/err_hnd.h>
37 #ifndef DOXYGEN_NO_O2NS
142 template<
class func_t=
funct,
170 }
else if (npoints==3) {
173 }
else if (npoints==4) {
180 if (npoints<=1 || npoints>5) {
181 O2SCL_ERR(
"Invalid # of points in set_npoints(). Using default",
195 }
else if (npoints==4) {
200 if (npoints<=2 || npoints>5) {
201 O2SCL_ERR(
"Invalid # of points in set_npoints2(). Using default",
210 double &dfdx,
double &err) {
212 dfdx=(this->*
cp)(x,p,func)/
h;
220 double &dfdx,
double &err) {
222 dfdx=(this->*
c2p)(x,p,func)/
h/
h;
230 double &dfdx,
double &err) {
232 dfdx=(this->*
c3p)(x,
h,p,func)/
h;
248 size_t nx,
const vec_t &y) {
249 size_t ix=(size_t)((x-x0)/dx);
250 return (this->*
cap)(x,x0,dx,nx,y,ix)/dx;
264 size_t nx,
const vec_t &y)
266 size_t ix=(size_t)((x-x0)/dx);
267 return (this->*
c2ap)(x,x0,dx,nx,y,ix)/dx;
281 size_t nx,
const vec_t &y)
283 size_t ix=(size_t)((x-x0)/dx);
284 return (this->*
c3ap)(x,x0,dx,nx,y,ix)/dx;
296 dydx[0]=(-25.0/12.0*y[0]+4.0*y[1]-3.0*y[2]+4.0/3.0*y[3]-0.25*y[4])/dx;
297 dydx[1]=(-0.25*y[0]-5.0/6.0*y[1]+1.5*y[2]-0.5*y[3]+1.0/12.0*y[4])/dx;
298 for(
size_t i=2;i<nv-2;i++) {
299 dydx[i]=(1.0/12.0*y[i-2]-2.0/3.0*y[i-1]+2.0/3.0*y[i+1]-
302 dydx[nv-2]=(-1.0/12.0*y[nv-5]+0.5*y[nv-4]-1.5*y[nv-3]+
303 5.0/6.0*y[nv-2]+0.25*y[nv-1])/dx;
304 dydx[nv-1]=(0.25*y[nv-5]-4.0/3.0*y[nv-4]+3.0*y[nv-3]-
305 4.0*y[nv-2]+25.0/12.0*y[nv-1])/dx;
310 virtual const char *
type() {
return "deriv_eqi"; }
312 #ifndef DOXYGEN_NO_O2NS
323 (
double x,
funct &func,
double &dfdx,
double &err) {
328 double derivp2(
double x,
double p, func_t &func) {
329 return (func(x+
h)-func(x));
333 double derivp3(
double x,
double p, func_t &func) {
335 return ((-0.5)*func(x-
h)+(0.5)*func(x+
h));
337 return ((p-0.5)*func(x-
h)-2.0*p*func(x)+(p+0.5)*func(x+
h));
344 return (-(3.0*p2-6.0*p+2.0)/6.0*func(x-
h)+
345 (1.5*p2-2.0*p-0.5)*func(x)-
346 (1.5*p2-p-1.0)*func(x+
h)+
347 (3.0*p2-1.0)/6.0*func(x+2.0*
h));
353 double p2=p*p, p3=p*p*p;
355 return ((1.0)/12.0*func(x-2.0*
h)-
357 (-4.0)/6.0*func(x+
h)+
358 (-1.0)/12.0*func(x+2.0*
h));
360 return ((2.0*p3-3.0*p2-p+1.0)/12.0*func(x-2.0*
h)-
361 (4.0*p3-3.0*p2-8.0*p+4.0)/6.0*func(x-
h)+
363 (4.0*p3+3.0*p2-8.0*p-4.0)/6.0*func(x+
h)+
364 (2.0*p3+3.0*p2-p-1.0)/12.0*func(x+2.0*
h));
370 const vec_t &y,
size_t ix) {
372 if (ix>0 && ix<nx-1) {
374 return ((p-0.5)*y[ix-1]-2.0*p*y[ix]+(p+0.5)*y[ix+1]);
377 return ((p-0.5)*y[0]-2.0*p*y[1]+(p+0.5)*y[2]);
380 return ((p-0.5)*y[nx-3]-2.0*p*y[nx-2]+(p+0.5)*y[nx-1]);
385 const vec_t &y,
size_t ix) {
387 if (ix>0 && ix<nx-2) {
390 return (-(3.0*p2-6.0*p+2.0)/6.0*y[ix-1]+
391 (1.5*p2-2.0*p-0.5)*y[ix]-
392 (1.5*p2-p-1.0)*y[ix+1]+
393 (3.0*p2-1.0)/6.0*y[ix+2]);
397 return (-(3.0*p2-6.0*p+2.0)/6.0*y[0]+
398 (1.5*p2-2.0*p-0.5)*y[1]-
400 (3.0*p2-1.0)/6.0*y[3]);
404 return (-(3.0*p2-6.0*p+2.0)/6.0*y[nx-4]+
405 (1.5*p2-2.0*p-0.5)*y[nx-3]-
406 (1.5*p2-p-1.0)*y[nx-2]+
407 (3.0*p2-1.0)/6.0*y[nx-1]);
412 double dx,
size_t nx,
413 const vec_t &y,
size_t ix) {
415 if (ix>1 && ix<nx-2) {
418 return ((2.0*p3-3.0*p2-p+1.0)/12.0*y[ix-2]-
419 (4.0*p3-3.0*p2-8.0*p+4.0)/6.0*y[ix-1]+
421 (4.0*p3+3.0*p2-8.0*p-4.0)/6.0*y[ix+1]+
422 (2.0*p3+3.0*p2-p-1.0)/12.0*y[ix+2]);
426 return ((2.0*p3-3.0*p2-p+1.0)/12.0*y[0]-
427 (4.0*p3-3.0*p2-8.0*p+4.0)/6.0*y[1]+
429 (4.0*p3+3.0*p2-8.0*p-4.0)/6.0*y[3]+
430 (2.0*p3+3.0*p2-p-1.0)/12.0*y[4]);
434 return ((2.0*p3-3.0*p2-p+1.0)/12.0*y[nx-5]-
435 (4.0*p3-3.0*p2-8.0*p+4.0)/6.0*y[nx-4]+
437 (4.0*p3+3.0*p2-8.0*p-4.0)/6.0*y[nx-2]+
438 (2.0*p3+3.0*p2-p-1.0)/12.0*y[nx-1]);
442 double deriv2p3(
double x,
double p, func_t &func) {
443 return (func(x-
h)-2*func(x)+func(x+
h));
447 double deriv2p4(
double x,
double p, func_t &func) {
448 return ((1.0-2.0*p)*func(x-
h)-(1.0-6.0*p)*func(x)
449 -(1.0-6.0*p)*func(x+
h)+(1.0+2.0*p)*func(x+2.0*
h))/2.0;
453 double deriv2p5(
double x,
double p, func_t &func) {
454 return ((1.0-2.0*p)*(1.0-2.0*p)*func(x-2.0*
h)
455 +(8.0*p-16.0*p*p)*func(x-
h)
456 -(2.0-24.0*p*p)*func(x)
457 -(8.0*p+16.0*p*p)*func(x+
h)
458 +(1.0+2.0*p)*(1.0+2.0*p)*func(x+2.0*
h))/4.0;
467 double dx,
size_t nx,
468 const vec_t &y,
size_t ix);
476 double dx,
size_t nx,
477 const vec_t &y,
size_t ix);
485 double dx,
size_t nx,
486 const vec_t &y,
size_t ix);
492 #ifndef DOXYGEN_NO_O2NS
double(deriv_eqi::* c2p)(double x, double p, func_t &func)
Pointer to the second derivative function.
double derivp4(double x, double p, func_t &func)
Four-point first derivative.
virtual int deriv2_err(double x, func_t &func, double &dfdx, double &err)
Calculate the second derivative of func w.r.t. x.
double xoff
Offset (default 0.0)
int set_npoints2(int npoints)
Set the number of points to use for second derivatives (default 5)
The main O<span style='position: relative; top: 0.3em; font-size: 0.8em'>2</span>scl O$_2$scl names...
int set_npoints(int npoints)
Set the number of points to use for first derivatives (default 5)
Numerical differentiation base [abstract base].
double deriv2_vector(double x, double x0, double dx, size_t nx, const vec_t &y)
Calculate the second derivative at x given an array.
double derivp3(double x, double p, func_t &func)
Three-point first derivative.
double(deriv_eqi::* c3ap)(double x, double x0, double dx, size_t nx, const vec_t &y, size_t ix)
Pointer to the third derivative for arrays function.
double deriv2p5(double x, double p, func_t &func)
Five-point second derivative.
double deriv_vector(double x, double x0, double dx, size_t nx, const vec_t &y)
Calculate the derivative at x given an array.
virtual int deriv3_err(double x, func_t &func, double &dfdx, double &err)
Calculate the third derivative of func w.r.t. x.
double(deriv_eqi::* cap)(double x, double x0, double dx, size_t nx, const vec_t &y, size_t ix)
Pointer to the first derivative for arrays function.
double(deriv_eqi::* c3p)(double x, double h, double p, func_t &func)
Pointer to the third derivative function.
virtual int deriv_err_int(double x, funct &func, double &dfdx, double &err)
Calculate the first derivative of func w.r.t. x and the uncertainty.
@ exc_einval
invalid argument supplied by user
double deriv_vector5(double x, double x0, double dx, size_t nx, const vec_t &y, size_t ix)
Five-point first derivative for arrays.
double deriv2p3(double x, double p, func_t &func)
Three-point second derivative.
Derivatives for equally-spaced abscissas.
double derivp5(double x, double p, func_t &func)
Five-point first derivative.
double deriv_vector4(double x, double x0, double dx, size_t nx, const vec_t &y, size_t ix)
Four-point first derivative for arrays.
int deriv_vector(size_t nv, double dx, const vec_t &y, vec_t &dydx)
Calculate the derivative of an entire array.
double h
Stepsize (Default ).
virtual int deriv_err(double x, func_t &func, double &dfdx, double &err)
Calculate the first derivative of func w.r.t. x.
double derivp2(double x, double p, func_t &func)
Two-point first derivative.
virtual const char * type()
Return string denoting type ("deriv_eqi")
double deriv2p4(double x, double p, func_t &func)
Four-point second derivative.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
std::function< double(double)> funct
One-dimensional function typedef in src/base/funct.h.
double(deriv_eqi::* cp)(double x, double p, func_t &func)
Pointer to the first derivative function.
double deriv_vector3(double x, double x0, double dx, size_t nx, const vec_t &y, size_t ix)
Three-point first derivative for arrays.
double deriv3_vector(double x, double x0, double dx, size_t nx, const vec_t &y)
Calculate the third derivative at x given an array.
double(deriv_eqi::* c2ap)(double x, double x0, double dx, size_t nx, const vec_t &y, size_t ix)
Pointer to the second derivative for arrays function.
Documentation generated with Doxygen. Provided under the
GNU Free Documentation License (see License Information).