www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
MonotoneCubicInterpolation Class Reference

This class interpolates values given a set of data pairs and an abscissa. More...

#include <MonotoneCubicInterpolation.h>

Public Member Functions

 MonotoneCubicInterpolation ()
 Empty constructor. More...
 
 MonotoneCubicInterpolation (const std::vector< Real > &x, const std::vector< Real > &y)
 Constructor, Takes two vectors of points for which to apply the fit. More...
 
virtual ~MonotoneCubicInterpolation ()=default
 
virtual void setData (const std::vector< Real > &x, const std::vector< Real > &y)
 Method generally used when MonotoneCubicInterpolation object was created using the empty constructor. More...
 
virtual Real sample (const Real &x) const
 This function will take an independent variable input and will return the dependent variable based on the generated fit. More...
 
virtual Real sampleDerivative (const Real &x) const
 This function will take an independent variable input and will return the derivative of the dependent variable with respect to the independent variable based on the generated fit. More...
 
virtual Real sample2ndDerivative (const Real &x) const
 This function will take an independent variable input and will return the second derivative of the dependent variable with respect to the independent variable based on the generated fit. More...
 
virtual void dumpCSV (std::string filename, const std::vector< Real > &xnew)
 This function takes an array of independent variable values and writes a CSV file with values corresponding to y, y', and y''. More...
 
virtual unsigned int getSampleSize ()
 This method returns the length of the independent variable vector. More...
 

Protected Member Functions

virtual void errorCheck ()
 
Real sign (const Real &x) const
 
Real phi (const Real &t) const
 
Real psi (const Real &t) const
 
Real phiPrime (const Real &t) const
 
Real psiPrime (const Real &t) const
 
Real phiDoublePrime (const Real &t) const
 
Real psiDoublePrime (const Real &t) const
 
Real h1 (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h2 (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h3 (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h4 (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h1Prime (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h2Prime (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h3Prime (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h4Prime (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h1DoublePrime (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h2DoublePrime (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h3DoublePrime (const Real &xhi, const Real &xlo, const Real &x) const
 
Real h4DoublePrime (const Real &xhi, const Real &xlo, const Real &x) const
 
virtual Real p (const Real &xhi, const Real &xlo, const Real &fhi, const Real &flo, const Real &dhi, const Real &dlo, const Real &x) const
 
virtual Real pPrime (const Real &xhi, const Real &xlo, const Real &fhi, const Real &flo, const Real &dhi, const Real &dlo, const Real &x) const
 
virtual Real pDoublePrime (const Real &xhi, const Real &xlo, const Real &fhi, const Real &flo, const Real &dhi, const Real &dlo, const Real &x) const
 
virtual void initialize_derivs ()
 
virtual void modify_derivs (const Real &alpha, const Real &beta, const Real &delta, Real &yp_lo, Real &yp_hi)
 
virtual void solve ()
 
virtual void findInterval (const Real &x, unsigned int &klo, unsigned int &khi) const
 

Protected Attributes

std::vector< Real > _x
 
std::vector< Real > _y
 
std::vector< Real > _h
 
std::vector< Real > _yp
 
std::vector< Real > _delta
 
std::vector< Real > _alpha
 
std::vector< Real > _beta
 
unsigned int _n_knots
 
unsigned int _n_intervals
 
unsigned int _internal_knots
 

Detailed Description

This class interpolates values given a set of data pairs and an abscissa.

The interpolation is cubic with at least C1 continuity; C2 continuity can be violated in favor of ensuring the interpolation is monotonic. The algorithm used is laid out in Fritsch and Carlson, SIAM J. Numer. Anal. Vol. 17(2) April 1980

Definition at line 26 of file MonotoneCubicInterpolation.h.

Constructor & Destructor Documentation

◆ MonotoneCubicInterpolation() [1/2]

MonotoneCubicInterpolation::MonotoneCubicInterpolation ( )

Empty constructor.

Definition at line 19 of file MonotoneCubicInterpolation.C.

19 {}

◆ MonotoneCubicInterpolation() [2/2]

MonotoneCubicInterpolation::MonotoneCubicInterpolation ( const std::vector< Real > &  x,
const std::vector< Real > &  y 
)

Constructor, Takes two vectors of points for which to apply the fit.

One should be of the independent variable while the other should be of the dependent variable. These values should have a one-to-one correspondence, e.g. the vectors must be of the same size.

Definition at line 21 of file MonotoneCubicInterpolation.C.

23  : _x(x), _y(y)
24 {
25  errorCheck();
26  solve();
27 }
static PetscErrorCode Vec x

◆ ~MonotoneCubicInterpolation()

virtual MonotoneCubicInterpolation::~MonotoneCubicInterpolation ( )
virtualdefault

Member Function Documentation

◆ dumpCSV()

void MonotoneCubicInterpolation::dumpCSV ( std::string  filename,
const std::vector< Real > &  xnew 
)
virtual

This function takes an array of independent variable values and writes a CSV file with values corresponding to y, y', and y''.

This can be used for sanity checks of the interpolation curve.

Definition at line 370 of file MonotoneCubicInterpolation.C.

371 {
372  unsigned int n = xnew.size();
373  std::vector<Real> ynew(n), ypnew(n), yppnew(n);
374 
375  std::ofstream out(filename.c_str());
376  for (unsigned int i = 0; i < n; ++i)
377  {
378  ynew[i] = sample(xnew[i]);
379  ypnew[i] = sampleDerivative(xnew[i]);
380  yppnew[i] = sample2ndDerivative(xnew[i]);
381  out << xnew[i] << ", " << ynew[i] << ", " << ypnew[i] << ", " << yppnew[i] << "\n";
382  }
383  out.close();
384 }
virtual Real sample2ndDerivative(const Real &x) const
This function will take an independent variable input and will return the second derivative of the de...
virtual Real sampleDerivative(const Real &x) const
This function will take an independent variable input and will return the derivative of the dependent...
virtual Real sample(const Real &x) const
This function will take an independent variable input and will return the dependent variable based on...
PetscInt n

◆ errorCheck()

void MonotoneCubicInterpolation::errorCheck ( )
protectedvirtual

Definition at line 39 of file MonotoneCubicInterpolation.C.

Referenced by MonotoneCubicInterpolation(), and setData().

40 {
41  if (_x.size() != _y.size())
42  throw std::domain_error("MonotoneCubicInterpolation: x and y vectors are not the same length");
43 
44  bool error = false;
45  for (unsigned i = 0; !error && i + 1 < _x.size(); ++i)
46  if (_x[i] >= _x[i + 1])
47  error = true;
48 
49  if (error)
50  throw std::domain_error("x-values are not strictly increasing");
51 }

◆ findInterval()

void MonotoneCubicInterpolation::findInterval ( const Real &  x,
unsigned int &  klo,
unsigned int &  khi 
) const
protectedvirtual

Definition at line 325 of file MonotoneCubicInterpolation.C.

Referenced by sample(), sample2ndDerivative(), and sampleDerivative().

328 {
329  klo = 0;
330  khi = _n_knots - 1;
331  while (khi - klo > 1)
332  {
333  unsigned int k = (khi + klo) >> 1;
334  if (_x[k] > x)
335  khi = k;
336  else
337  klo = k;
338  }
339 }
static PetscErrorCode Vec x

◆ getSampleSize()

unsigned int MonotoneCubicInterpolation::getSampleSize ( )
virtual

This method returns the length of the independent variable vector.

Definition at line 387 of file MonotoneCubicInterpolation.C.

388 {
389  return _x.size();
390 }

◆ h1()

Real MonotoneCubicInterpolation::h1 ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 101 of file MonotoneCubicInterpolation.C.

Referenced by p().

102 {
103  Real h = xhi - xlo;
104  Real t = (xhi - x) / h;
105  return phi(t);
106 }
static PetscErrorCode Vec x

◆ h1DoublePrime()

Real MonotoneCubicInterpolation::h1DoublePrime ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 118 of file MonotoneCubicInterpolation.C.

Referenced by pDoublePrime().

119 {
120  Real h = xhi - xlo;
121  Real t = (xhi - x) / h;
122  Real tPrime = -1. / h;
123  return phiDoublePrime(t) * tPrime * tPrime;
124 }
static PetscErrorCode Vec x
Real phiDoublePrime(const Real &t) const

◆ h1Prime()

Real MonotoneCubicInterpolation::h1Prime ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 109 of file MonotoneCubicInterpolation.C.

Referenced by pPrime().

110 {
111  Real h = xhi - xlo;
112  Real t = (xhi - x) / h;
113  Real tPrime = -1. / h;
114  return phiPrime(t) * tPrime;
115 }
static PetscErrorCode Vec x
Real phiPrime(const Real &t) const

◆ h2()

Real MonotoneCubicInterpolation::h2 ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 127 of file MonotoneCubicInterpolation.C.

Referenced by p().

128 {
129  Real h = xhi - xlo;
130  Real t = (x - xlo) / h;
131  return phi(t);
132 }
static PetscErrorCode Vec x

◆ h2DoublePrime()

Real MonotoneCubicInterpolation::h2DoublePrime ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 144 of file MonotoneCubicInterpolation.C.

Referenced by pDoublePrime().

145 {
146  Real h = xhi - xlo;
147  Real t = (x - xlo) / h;
148  Real tPrime = 1. / h;
149  return phiDoublePrime(t) * tPrime * tPrime;
150 }
static PetscErrorCode Vec x
Real phiDoublePrime(const Real &t) const

◆ h2Prime()

Real MonotoneCubicInterpolation::h2Prime ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 135 of file MonotoneCubicInterpolation.C.

Referenced by pPrime().

136 {
137  Real h = xhi - xlo;
138  Real t = (x - xlo) / h;
139  Real tPrime = 1. / h;
140  return phiPrime(t) * tPrime;
141 }
static PetscErrorCode Vec x
Real phiPrime(const Real &t) const

◆ h3()

Real MonotoneCubicInterpolation::h3 ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 153 of file MonotoneCubicInterpolation.C.

Referenced by p().

154 {
155  Real h = xhi - xlo;
156  Real t = (xhi - x) / h;
157  return -h * psi(t);
158 }
static PetscErrorCode Vec x

◆ h3DoublePrime()

Real MonotoneCubicInterpolation::h3DoublePrime ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 170 of file MonotoneCubicInterpolation.C.

Referenced by pDoublePrime().

171 {
172  Real h = xhi - xlo;
173  Real t = (xhi - x) / h;
174  Real tPrime = -1. / h;
175  return psiDoublePrime(t) * tPrime;
176 }
static PetscErrorCode Vec x
Real psiDoublePrime(const Real &t) const

◆ h3Prime()

Real MonotoneCubicInterpolation::h3Prime ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 161 of file MonotoneCubicInterpolation.C.

Referenced by pPrime().

162 {
163  Real h = xhi - xlo;
164  Real t = (xhi - x) / h;
165  Real tPrime = -1. / h;
166  return -h * psiPrime(t) * tPrime; // psiPrime(t)
167 }
Real psiPrime(const Real &t) const
static PetscErrorCode Vec x

◆ h4()

Real MonotoneCubicInterpolation::h4 ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 179 of file MonotoneCubicInterpolation.C.

Referenced by p().

180 {
181  Real h = xhi - xlo;
182  Real t = (x - xlo) / h;
183  return h * psi(t);
184 }
static PetscErrorCode Vec x

◆ h4DoublePrime()

Real MonotoneCubicInterpolation::h4DoublePrime ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 196 of file MonotoneCubicInterpolation.C.

Referenced by pDoublePrime().

197 {
198  Real h = xhi - xlo;
199  Real t = (x - xlo) / h;
200  Real tPrime = 1. / h;
201  return psiDoublePrime(t) * tPrime;
202 }
static PetscErrorCode Vec x
Real psiDoublePrime(const Real &t) const

◆ h4Prime()

Real MonotoneCubicInterpolation::h4Prime ( const Real &  xhi,
const Real &  xlo,
const Real &  x 
) const
protected

Definition at line 187 of file MonotoneCubicInterpolation.C.

Referenced by pPrime().

188 {
189  Real h = xhi - xlo;
190  Real t = (x - xlo) / h;
191  Real tPrime = 1. / h;
192  return h * psiPrime(t) * tPrime; // psiPrime(t)
193 }
Real psiPrime(const Real &t) const
static PetscErrorCode Vec x

◆ initialize_derivs()

void MonotoneCubicInterpolation::initialize_derivs ( )
protectedvirtual

Definition at line 244 of file MonotoneCubicInterpolation.C.

Referenced by solve().

245 {
246  for (unsigned int i = 1; i < _n_knots - 1; ++i)
247  _yp[i] = (std::pow(_h[i - 1], 2) * _y[i + 1] - std::pow(_h[i], 2) * _y[i - 1] -
248  _y[i] * (_h[i - 1] - _h[i]) * (_h[i - 1] + _h[i])) /
249  (_h[i - 1] * _h[i] * (_h[i - 1] * _h[i]));
250 
251  _yp[0] = (-std::pow(_h[0], 2) * _y[2] - _h[1] * _y[0] * (2 * _h[0] + _h[1]) +
252  _y[1] * std::pow(_h[0] + _h[1], 2)) /
253  (_h[0] * _h[1] * (_h[0] + _h[1]));
254 
255  Real hlast = _h[_n_intervals - 1];
256  Real hsecond = _h[_n_intervals - 2];
257  Real ylast = _y[_n_knots - 1];
258  Real ysecond = _y[_n_knots - 2];
259  Real ythird = _y[_n_knots - 3];
260  _yp[_n_knots - 1] = (hsecond * ylast * (hsecond + 2 * hlast) + std::pow(hlast, 2) * ythird -
261  ysecond * std::pow(hsecond + hlast, 2)) /
262  (hsecond * hlast * (hsecond + hlast));
263 }
Real pow(Real x, int e)
Definition: MathUtils.C:211

◆ modify_derivs()

void MonotoneCubicInterpolation::modify_derivs ( const Real &  alpha,
const Real &  beta,
const Real &  delta,
Real &  yp_lo,
Real &  yp_hi 
)
protectedvirtual

Definition at line 266 of file MonotoneCubicInterpolation.C.

Referenced by solve().

268 {
269  Real tau = 3. / std::sqrt(std::pow(alpha, 2) + std::pow(beta, 2));
270  Real alpha_star = alpha * tau;
271  Real beta_star = beta * tau;
272  yp_lo = alpha_star * delta;
273  yp_hi = beta_star * delta;
274 }
Real pow(Real x, int e)
Definition: MathUtils.C:211

◆ p()

Real MonotoneCubicInterpolation::p ( const Real &  xhi,
const Real &  xlo,
const Real &  fhi,
const Real &  flo,
const Real &  dhi,
const Real &  dlo,
const Real &  x 
) const
protectedvirtual

Definition at line 205 of file MonotoneCubicInterpolation.C.

Referenced by sample().

212 {
213  return flo * h1(xhi, xlo, x) + fhi * h2(xhi, xlo, x) + dlo * h3(xhi, xlo, x) +
214  dhi * h4(xhi, xlo, x);
215 }
static PetscErrorCode Vec x
Real h4(const Real &xhi, const Real &xlo, const Real &x) const
Real h2(const Real &xhi, const Real &xlo, const Real &x) const
Real h3(const Real &xhi, const Real &xlo, const Real &x) const
Real h1(const Real &xhi, const Real &xlo, const Real &x) const

◆ pDoublePrime()

Real MonotoneCubicInterpolation::pDoublePrime ( const Real &  xhi,
const Real &  xlo,
const Real &  fhi,
const Real &  flo,
const Real &  dhi,
const Real &  dlo,
const Real &  x 
) const
protectedvirtual

Definition at line 231 of file MonotoneCubicInterpolation.C.

Referenced by sample2ndDerivative().

238 {
239  return flo * h1DoublePrime(xhi, xlo, x) + fhi * h2DoublePrime(xhi, xlo, x) +
240  dlo * h3DoublePrime(xhi, xlo, x) + dhi * h4DoublePrime(xhi, xlo, x);
241 }
Real h3DoublePrime(const Real &xhi, const Real &xlo, const Real &x) const
Real h4DoublePrime(const Real &xhi, const Real &xlo, const Real &x) const
static PetscErrorCode Vec x
Real h1DoublePrime(const Real &xhi, const Real &xlo, const Real &x) const
Real h2DoublePrime(const Real &xhi, const Real &xlo, const Real &x) const

◆ phi()

Real MonotoneCubicInterpolation::phi ( const Real &  t) const
protected

Definition at line 65 of file MonotoneCubicInterpolation.C.

Referenced by h1(), and h2().

66 {
67  return 3. * t * t - 2. * t * t * t;
68 }

◆ phiDoublePrime()

Real MonotoneCubicInterpolation::phiDoublePrime ( const Real &  t) const
protected

Definition at line 77 of file MonotoneCubicInterpolation.C.

Referenced by h1DoublePrime(), and h2DoublePrime().

78 {
79  return 6. - 12. * t;
80 }

◆ phiPrime()

Real MonotoneCubicInterpolation::phiPrime ( const Real &  t) const
protected

Definition at line 71 of file MonotoneCubicInterpolation.C.

Referenced by h1Prime(), and h2Prime().

72 {
73  return 6. * t - 6. * t * t;
74 }

◆ pPrime()

Real MonotoneCubicInterpolation::pPrime ( const Real &  xhi,
const Real &  xlo,
const Real &  fhi,
const Real &  flo,
const Real &  dhi,
const Real &  dlo,
const Real &  x 
) const
protectedvirtual

Definition at line 218 of file MonotoneCubicInterpolation.C.

Referenced by sampleDerivative().

225 {
226  return flo * h1Prime(xhi, xlo, x) + fhi * h2Prime(xhi, xlo, x) + dlo * h3Prime(xhi, xlo, x) +
227  dhi * h4Prime(xhi, xlo, x);
228 }
Real h1Prime(const Real &xhi, const Real &xlo, const Real &x) const
Real h2Prime(const Real &xhi, const Real &xlo, const Real &x) const
Real h3Prime(const Real &xhi, const Real &xlo, const Real &x) const
static PetscErrorCode Vec x
Real h4Prime(const Real &xhi, const Real &xlo, const Real &x) const

◆ psi()

Real MonotoneCubicInterpolation::psi ( const Real &  t) const
protected

Definition at line 83 of file MonotoneCubicInterpolation.C.

Referenced by h3(), and h4().

84 {
85  return t * t * t - t * t;
86 }

◆ psiDoublePrime()

Real MonotoneCubicInterpolation::psiDoublePrime ( const Real &  t) const
protected

Definition at line 95 of file MonotoneCubicInterpolation.C.

Referenced by h3DoublePrime(), and h4DoublePrime().

96 {
97  return 6. * t - 2.;
98 }

◆ psiPrime()

Real MonotoneCubicInterpolation::psiPrime ( const Real &  t) const
protected

Definition at line 89 of file MonotoneCubicInterpolation.C.

Referenced by h3Prime(), and h4Prime().

90 {
91  return 3. * t * t - 2. * t;
92 }

◆ sample()

Real MonotoneCubicInterpolation::sample ( const Real &  x) const
virtual

This function will take an independent variable input and will return the dependent variable based on the generated fit.

Definition at line 342 of file MonotoneCubicInterpolation.C.

Referenced by dumpCSV().

343 {
344  // sanity check (empty MontoneCubicInterpolations are constructable
345  // so we cannot put this into the errorCheck)
346  assert(_x.size() > 0);
347 
348  unsigned int klo, khi;
349  findInterval(x, klo, khi);
350  return p(_x[khi], _x[klo], _y[khi], _y[klo], _yp[khi], _yp[klo], x);
351 }
static PetscErrorCode Vec x
virtual void findInterval(const Real &x, unsigned int &klo, unsigned int &khi) const
virtual Real p(const Real &xhi, const Real &xlo, const Real &fhi, const Real &flo, const Real &dhi, const Real &dlo, const Real &x) const

◆ sample2ndDerivative()

Real MonotoneCubicInterpolation::sample2ndDerivative ( const Real &  x) const
virtual

This function will take an independent variable input and will return the second derivative of the dependent variable with respect to the independent variable based on the generated fit.

Note that this can be discontinous at the knots.

Definition at line 362 of file MonotoneCubicInterpolation.C.

Referenced by dumpCSV().

363 {
364  unsigned int klo, khi;
365  findInterval(x, klo, khi);
366  return pDoublePrime(_x[khi], _x[klo], _y[khi], _y[klo], _yp[khi], _yp[klo], x);
367 }
static PetscErrorCode Vec x
virtual Real pDoublePrime(const Real &xhi, const Real &xlo, const Real &fhi, const Real &flo, const Real &dhi, const Real &dlo, const Real &x) const
virtual void findInterval(const Real &x, unsigned int &klo, unsigned int &khi) const

◆ sampleDerivative()

Real MonotoneCubicInterpolation::sampleDerivative ( const Real &  x) const
virtual

This function will take an independent variable input and will return the derivative of the dependent variable with respect to the independent variable based on the generated fit.

Definition at line 354 of file MonotoneCubicInterpolation.C.

Referenced by dumpCSV().

355 {
356  unsigned int klo, khi;
357  findInterval(x, klo, khi);
358  return pPrime(_x[khi], _x[klo], _y[khi], _y[klo], _yp[khi], _yp[klo], x);
359 }
static PetscErrorCode Vec x
virtual Real pPrime(const Real &xhi, const Real &xlo, const Real &fhi, const Real &flo, const Real &dhi, const Real &dlo, const Real &x) const
virtual void findInterval(const Real &x, unsigned int &klo, unsigned int &khi) const

◆ setData()

void MonotoneCubicInterpolation::setData ( const std::vector< Real > &  x,
const std::vector< Real > &  y 
)
virtual

Method generally used when MonotoneCubicInterpolation object was created using the empty constructor.

Takes two vectors of points for which to apply the fit. One should be of the independent variable while the other should be of the dependent variable. These values should have a one-to-one correspondence, e.g. the vectors must be of the same size.

Definition at line 30 of file MonotoneCubicInterpolation.C.

31 {
32  _x = x;
33  _y = y;
34  errorCheck();
35  solve();
36 }
static PetscErrorCode Vec x

◆ sign()

Real MonotoneCubicInterpolation::sign ( const Real &  x) const
protected

Definition at line 54 of file MonotoneCubicInterpolation.C.

Referenced by solve().

55 {
56  if (x < 0)
57  return -1;
58  else if (x > 0)
59  return 1;
60  else
61  return 0;
62 }
static PetscErrorCode Vec x

◆ solve()

void MonotoneCubicInterpolation::solve ( )
protectedvirtual

Definition at line 277 of file MonotoneCubicInterpolation.C.

Referenced by MonotoneCubicInterpolation(), and setData().

278 {
279  _n_knots = _x.size(), _n_intervals = _x.size() - 1, _internal_knots = _x.size() - 2;
280  _h.resize(_n_intervals);
281  _yp.resize(_n_knots);
282  _delta.resize(_n_intervals);
283  _alpha.resize(_n_intervals);
284  _beta.resize(_n_intervals);
285 
286  for (unsigned int i = 0; i < _n_intervals; ++i)
287  _h[i] = _x[i + 1] - _x[i];
288 
290  for (unsigned int i = 0; i < _n_intervals; ++i)
291  _delta[i] = (_y[i + 1] - _y[i]) / _h[i];
292  if (sign(_delta[0]) != sign(_yp[0]))
293  _yp[0] = 0;
294  if (sign(_delta[_n_intervals - 1]) != sign(_yp[_n_knots - 1]))
295  _yp[_n_knots - 1] = 0;
296  for (unsigned int i = 0; i < _internal_knots; ++i)
297  if (sign(_delta[i + 1]) == 0 || sign(_delta[i]) == 0 || sign(_delta[i + 1]) != sign(_delta[i]))
298  _yp[1 + i] = 0;
299 
300  for (unsigned int i = 0; i < _n_intervals; ++i)
301  {
302  // Test for zero slope
303  if (_yp[i] == 0)
304  _alpha[i] = 0;
305  else if (_delta[i] == 0)
306  _alpha[i] = 4;
307  else
308  _alpha[i] = _yp[i] / _delta[i];
309 
310  // Test for zero slope
311  if (_yp[i + 1] == 0)
312  _beta[i] = 0;
313  else if (_delta[i] == 0)
314  _beta[i] = 4;
315  else
316  _beta[i] = _yp[i + 1] / _delta[i];
317 
318  // check if outside region of monotonicity
319  if (std::pow(_alpha[i], 2) + std::pow(_beta[i], 2) > 9.)
320  modify_derivs(_alpha[i], _beta[i], _delta[i], _yp[i], _yp[i + 1]);
321  }
322 }
virtual void modify_derivs(const Real &alpha, const Real &beta, const Real &delta, Real &yp_lo, Real &yp_hi)
Real pow(Real x, int e)
Definition: MathUtils.C:211
Real sign(const Real &x) const

Member Data Documentation

◆ _alpha

std::vector<Real> MonotoneCubicInterpolation::_alpha
protected

Definition at line 147 of file MonotoneCubicInterpolation.h.

Referenced by solve().

◆ _beta

std::vector<Real> MonotoneCubicInterpolation::_beta
protected

Definition at line 148 of file MonotoneCubicInterpolation.h.

Referenced by solve().

◆ _delta

std::vector<Real> MonotoneCubicInterpolation::_delta
protected

Definition at line 146 of file MonotoneCubicInterpolation.h.

Referenced by solve().

◆ _h

std::vector<Real> MonotoneCubicInterpolation::_h
protected

Definition at line 144 of file MonotoneCubicInterpolation.h.

Referenced by initialize_derivs(), and solve().

◆ _internal_knots

unsigned int MonotoneCubicInterpolation::_internal_knots
protected

Definition at line 152 of file MonotoneCubicInterpolation.h.

Referenced by solve().

◆ _n_intervals

unsigned int MonotoneCubicInterpolation::_n_intervals
protected

Definition at line 151 of file MonotoneCubicInterpolation.h.

Referenced by initialize_derivs(), and solve().

◆ _n_knots

unsigned int MonotoneCubicInterpolation::_n_knots
protected

Definition at line 150 of file MonotoneCubicInterpolation.h.

Referenced by findInterval(), initialize_derivs(), and solve().

◆ _x

std::vector<Real> MonotoneCubicInterpolation::_x
protected

◆ _y

std::vector<Real> MonotoneCubicInterpolation::_y
protected

◆ _yp

std::vector<Real> MonotoneCubicInterpolation::_yp
protected

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