ode_iv_solve.h
Go to the documentation of this file.
1 /*
2  -------------------------------------------------------------------
3 
4  Copyright (C) 2006-2017, Andrew W. Steiner
5 
6  This file is part of O2scl.
7 
8  O2scl is free software; you can redistribute it and/or modify
9  it under the terms of the GNU General Public License as published by
10  the Free Software Foundation; either version 3 of the License, or
11  (at your option) any later version.
12 
13  O2scl is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with O2scl. If not, see <http://www.gnu.org/licenses/>.
20 
21  -------------------------------------------------------------------
22 */
23 #ifndef O2SCL_ODE_IV_SOLVE_H
24 #define O2SCL_ODE_IV_SOLVE_H
25 
26 /** \file ode_iv_solve.h
27  \brief File defining \ref o2scl::ode_iv_solve
28 */
29 
30 #include <string>
31 
32 #include <boost/numeric/ublas/vector.hpp>
33 #include <boost/numeric/ublas/matrix.hpp>
34 #include <boost/numeric/ublas/matrix_proxy.hpp>
35 
36 #include <o2scl/astep.h>
37 #include <o2scl/astep_gsl.h>
38 
39 #ifndef DOXYGEN_NO_O2NS
40 namespace o2scl {
41 #endif
42 
43  /** \brief Solve an initial-value ODE problems given an adaptive ODE
44  stepper
45 
46  This class is experimental.
47 
48  This class gives several functions which solve an initial
49  value ODE problem. The functions \ref solve_final_value()
50  give only the final value of the functions at the end
51  of the ODE integration and are relatively fast.
52 
53  The function solve_store() stores the solution of the ODE over
54  the full range into a set of vectors and matrices which are
55  allocated and specified by the user. This function is designed
56  to give exactly the same results (though this cannot be
57  guaranteed) as solve_final_value() and additionally records some
58  or all of the results from the adaptive steps which were taken.
59 
60  The functions solve_grid() works as in solve_store() except
61  that the solution is stored on a grid of points in the
62  independent variable specified by the user, at the cost of
63  taking extra steps to ensure that function values,
64  derivatives, and errors are computed at each grid point.
65 
66  All of these functions automatically evaluate the derivatives
67  from the specified function at the initial point and
68  user-specified initial derivatives are ignored. The total
69  number of steps taken is limited by \ref ntrial and \ref
70  nsteps stores the number of steps taken by the most recent
71  solution. The variable \ref nsteps_out is the maximum number
72  of points in the interval for which verbose output will be
73  given when \ref o2scl::ode_iv_solve::verbose is greater than zero.
74 
75  There is an example for the usage of this class in
76  <tt>examples/ex_ode.cpp</tt> documented in the \ref ex_ode_sect
77  section.
78 
79  Documentation links for default template arguments
80  - \c func_t - \ref ode_funct11
81  - \c vec_t - \ref boost::numeric::ublas::vector < double >
82 
83  The default adaptive stepper is an object of type \ref astep_gsl.
84 
85  \future The form of solve_final_value() is very similar to that
86  of astep_base::astep_full(), but not quite the same. Maybe
87  should probably be made to be consistent with each other?
88  */
89  template<class func_t=ode_funct11,
91  class ode_iv_solve {
92 
93  public:
94 
96 
97 #ifndef DOXYGEN_INTERNAL
98 
99  protected:
100 
101  /// \name Vectors for temporary storage
102  //@{
103  vec_t vtemp, vtemp2, vtemp3, vtemp4;
104  //@}
105 
106  /// Desc
107  size_t mem_size;
108 
109  /// The adaptive stepper
111 
112  /// Print out iteration information
113  virtual int print_iter(double x, size_t nv, vec_t &y) {
114  std::cout << type() << " x: " << x << " y: ";
115  for(size_t i=0;i<nv;i++) std::cout << y[i] << " ";
116  std::cout << std::endl;
117  if (verbose>1) {
118  char ch;
119  std::cin >> ch;
120  }
121  return 0;
122  }
123 
124  /// Free allocated memory
125  void free() {
126  if (mem_size>0) {
127  vtemp.clear();
128  vtemp2.clear();
129  vtemp3.clear();
130  vtemp4.clear();
131  }
132  }
133 
134  /** \brief Allocate space for temporary vectors
135  */
136  void allocate(size_t n) {
137  if (n!=mem_size) {
138  free();
139  vtemp.resize(n);
140  vtemp2.resize(n);
141  vtemp3.resize(n);
142  vtemp4.resize(n);
143  mem_size=n;
144  }
145  }
146 
147 #endif
148 
149  public:
150 
151  ode_iv_solve() {
152  verbose=0;
153  ntrial=1000;
154  nsteps_out=10;
155  astp=&gsl_astp;
156  exit_on_fail=true;
157  mem_size=0;
158  }
159 
160  virtual ~ode_iv_solve() {
161  free();
162  }
163 
164  /// \name Main solver functions
165  //@{
166  /** \brief Solve the initial-value problem to get the final value
167 
168  Given the \c n initial values of the functions in \c ystart,
169  this function integrates the ODEs specified in \c derivs over
170  the interval from \c x0 to \c x1 with an initial stepsize of \c
171  h. The final values of the function are given in \c yend and the
172  initial values of \c yend are ignored.
173 
174  If \ref verbose is greater than zero, The solution at less than
175  or approximately equal to \ref nsteps_out points will be written
176  to \c std::cout. If \ref verbose is greater than one, a
177  character will be required after each selected point.
178 
179  \todo Document if yend can be the same as ystart.
180  */
181  int solve_final_value(double x0, double x1, double h, size_t n,
182  vec_t &ystart, vec_t &yend, func_t &derivs) {
183 
184  allocate(n);
185  return solve_final_value(x0,x1,h,n,ystart,yend,vtemp2,
186  vtemp3,derivs);
187  }
188 
189  /** \brief Solve the initial-value problem to get the final value
190  with errors
191 
192  Given the \c n initial values of the functions in \c ystart,
193  this function integrates the ODEs specified in \c derivs over
194  the interval from \c x0 to \c x1 with an initial stepsize of \c
195  h. The final values of the function are given in \c yend and the
196  associated errors are given in \c yerr. The initial values of \c
197  yend and \c yerr are ignored.
198 
199  If \ref verbose is greater than zero, The solution at less
200  than or approximately equal to \ref nsteps_out points will be
201  written to \c std::cout. If \ref verbose is greater than one,
202  a character will be required after each selected point.
203  */
204  int solve_final_value(double x0, double x1, double h, size_t n,
205  vec_t &ystart, vec_t &yend, vec_t &yerr,
206  func_t &derivs) {
207 
208  allocate(n);
209  return solve_final_value(x0,x1,h,n,ystart,yend,yerr,
210  vtemp2,derivs);
211  }
212 
213  /** \brief Solve the initial-value problem to get the final value,
214  derivatives, and errors
215 
216  Given the \c n initial values of the functions in \c ystart,
217  this function integrates the ODEs specified in \c derivs over
218  the interval from \c x0 to \c x1 with an initial stepsize of \c
219  h. The final values of the function are given in \c yend, the
220  derivatives in \c dydx_end, and the associated errors are given
221  in \c yerr. The initial values of \c yend and \c yerr are
222  ignored.
223 
224  This function is designed to be relatively fast,
225  avoiding extra copying of vectors back and forth.
226 
227  If \ref verbose is greater than zero, The solution at less
228  than or approximately equal to \ref nsteps_out points will be
229  written to \c std::cout. If \ref verbose is greater than one,
230  a character will be required after each selected point.
231 
232  This function computes \c dydx_start automatically and the
233  values given by the user are ignored.
234 
235  The solution fails if more than \ref ntrial steps are required.
236  This function will also fail if <tt>x1>x0</tt> and <tt>h<0</tt>
237  or if <tt>x1-x0</tt> and <tt>h>0</tt> do not have the same sign.
238  */
239  int solve_final_value(double x0, double x1, double h, size_t n,
240  vec_t &ystart, vec_t &yend, vec_t &yerr,
241  vec_t &dydx_end, func_t &derivs) {
242 
243  if ((x1>x0 && h<=0.0) || (x0>x1 && h>=0.0)) {
244  std::string str="Interval direction (x1-x0="+o2scl::dtos(x1-x0)+
245  ") does not match step direction (h="+o2scl::dtos(h)+
246  " in ode_iv_solve::solve_final_value().";
247  O2SCL_ERR(str.c_str(),exc_einval);
248  }
249  if (x0==x1) {
250  O2SCL_ERR2("Starting and final endpoints identical in ",
251  "ode_iv_solve::solve_final_value().",exc_einval);
252  }
253 
254  // Allocate for temporary vectors
255  allocate(n);
256 
257  // The variable 'x' is the current independent variable, xnext is
258  // the next value of x, xverb is the next value of x for which
259  // verbose output should be provided, and dxverb is the stepsize
260  // for xverb.
261  double x=x0, xverb=0.0, dxverb=0.0, xnext;
262  int ret=0, first_ret=0;
263 
264  nsteps=0;
265 
266  // Create a reference to make the code easier to read
267  vec_t &dydx_start=vtemp;
268 
269  // Compute stepsize for verbose output
270  if (verbose>0) {
271  print_iter(x0,n,ystart);
272  if (verbose>1) {
273  char ch;
274  std::cin >> ch;
275  }
276  // Ensure that 'dxverb' is positive
277  if (x1>x0) {
278  dxverb=(x1-x0)/((double)nsteps_out);
279  xverb=x0+dxverb;
280  } else {
281  dxverb=(x0-x1)/((double)nsteps_out);
282  xverb=x0-dxverb;
283  }
284  }
285 
286  // Use yend as workspace
287  for(size_t i=0;i<n;i++) yend[i]=ystart[i];
288 
289  // Compute initial derivative
290  derivs(x,n,yend,dydx_start);
291 
292  bool done=false;
293  while (done==false) {
294 
295  // Take a step of the first type
296  ret=astp->astep_full(x,x1,xnext,h,n,ystart,dydx_start,
297  yend,yerr,dydx_end,derivs);
298 
299  if (ret!=0) {
300  if (exit_on_fail) {
301  O2SCL_ERR2("Adaptive stepper failed in ",
302  "ode_iv_solve::solve_final_value()",ret);
303  } else if (first_ret!=0) {
304  first_ret=ret;
305  }
306  }
307 
308  // Print out verbose info
309  if (verbose>0 && xnext>xverb) {
310  print_iter(xnext,n,yend);
311  if (verbose>1) {
312  char ch;
313  std::cin >> ch;
314  }
315  xverb+=dxverb;
316  }
317 
318  // Check number of iterations
319  nsteps++;
320  if (nsteps>ntrial) {
321  std::string str="Too many steps required (ntrial="+itos(ntrial)+
322  ", x="+o2scl::dtos(x)+") in ode_iv_solve::solve_final_value().";
323  O2SCL_ERR(str.c_str(),exc_emaxiter);
324  }
325 
326  if (ret!=0) {
327  done=true;
328  } else {
329  if (x1>x0) {
330  if (xnext>=x1) done=true;
331  } else {
332  if (xnext<=x1) done=true;
333  }
334  }
335 
336  if (done==false) {
337 
338  // Take a step of the second type
339  ret=astp->astep_full(xnext,x1,x,h,n,yend,dydx_end,ystart,yerr,
340  dydx_start,derivs);
341 
342  if (ret!=0) {
343  if (exit_on_fail) {
344  O2SCL_ERR2("Adaptive stepper failed in ",
345  "ode_iv_solve::solve_final_value()",ret);
346  } else if (first_ret!=0) {
347  first_ret=ret;
348  }
349  }
350 
351  // Print out verbose info
352  if (verbose>0 && x>xverb) {
353  print_iter(x,n,ystart);
354  if (verbose>1) {
355  char ch;
356  std::cin >> ch;
357  }
358  xverb+=dxverb;
359  }
360 
361  // Check number of iterations
362  nsteps++;
363  if (nsteps>ntrial) {
364  std::string str="Too many steps required (ntrial="+itos(ntrial)+
365  ", x="+o2scl::dtos(x)+") in ode_iv_solve::solve_final_value().";
366  O2SCL_ERR(str.c_str(),exc_emaxiter);
367  }
368 
369  if (ret!=0) {
370  done=true;
371  } else {
372  if (x1>x0) {
373  if (x>=x1) {
374  done=true;
375  for(size_t j=0;j<n;j++) {
376  yend[j]=ystart[j];
377  dydx_end[j]=dydx_start[j];
378  }
379  }
380  } else {
381  if (x<=x1) {
382  done=true;
383  for(size_t j=0;j<n;j++) {
384  yend[j]=ystart[j];
385  dydx_end[j]=dydx_start[j];
386  }
387  }
388  }
389  }
390 
391  }
392 
393  // End of while loop
394  }
395 
396  // Print out final verbose info
397  if (verbose>0) {
398  print_iter(x1,n,yend);
399  if (verbose>1) {
400  char ch;
401  std::cin >> ch;
402  }
403  }
404 
405  return first_ret;
406  }
407 
408  /** \brief Solve the initial-value problem and store the
409  associated output
410 
411  Initially, \c x_sol should be a vector of size \c n_sol, and \c
412  y_sol, \c dydx_sol, and \c yerr_sol should be matrices with size
413  \c [n_sol][n]. On exit, \c n_sol will will be number of points
414  store, less than or equal to the original value of \c
415  n_sol. This function avoids performing extra calls to the
416  adaptive stepper, and the table will be approximately evenly
417  spaced.
418 
419  This function is also designed to give the exactly the same
420  results as solve_final_value(). This cannot be strictly
421  guaranteed, but will likely hold in most applications.
422 
423  This template function works with any matrix class \c mat_t
424  which can be accessed using successive applications of
425  operator[] and which has an associated class \c mat_row_t
426  which returns a row of a matrix of type \c mat_t as
427  an object with type \c vec_t.
428 
429  If \ref verbose is greater than zero, The solution at each
430  internal point will be written to \c std::cout. If \ref verbose
431  is greater than one, a character will be required after each
432  point.
433 
434  \future It might be possible to remove some extra copying
435  by removing the yerrl and dydx vectors
436  */
437  template<class mat_t>
438  int solve_store(double x0, double x1, double h, size_t n, size_t &n_sol,
439  vec_t &x_sol, mat_t &y_sol, mat_t &yerr_sol,
440  mat_t &dydx_sol, func_t &derivs, size_t istart=0) {
441 
442  int ret=0;
443  int first_ret=0;
444  size_t nmax=n_sol-1;
445  nsteps=0;
446 
447  // Stepsize for next verbose output. Use nsteps_out-1 instead of
448  // nsteps_out since the first point is always output below.
449  double dx_verb=(x1-x0)/((double)(nsteps_out-1));
450  // Stepsize for next table entry
451  double dx_tab=(x1-x0)/((double)(n_sol-istart-1));
452 
453  double x_verb=x0+dx_verb;
454  double x_tab=x0+dx_tab;
455 
456  // Allocate space for errors and derivatives and storage
457  allocate(n);
458 
459  // Create some references just to make the code easier
460  // to read
461  vec_t &ystart=vtemp4;
462  vec_t &dydx_start=vtemp2;
463 
464  // Copy first point to ystart
465  for(size_t j=0;j<n;j++) ystart[j]=y_sol(istart,j);
466 
467  // Output first point
468  if (verbose>0) {
469  print_iter(x0,n,ystart);
470  if (verbose>1) {
471  char ch;
472  std::cin >> ch;
473  }
474  }
475 
476  // Initial derivative evaulation
477  derivs(x0,n,ystart,dydx_start);
478 
479  // Add first derivatives to table
480  x_sol[istart]=x0;
481  for(size_t j=0;j<n;j++) {
482  dydx_sol(istart,j)=dydx_start[j];
483  yerr_sol(istart,j)=0.0;
484  }
485 
486  // Copy first point to table again for first step
487  size_t icurr=istart+1;
488  x_sol[icurr]=x0;
489  for(size_t j=0;j<n;j++) y_sol(icurr,j)=ystart[j];
490 
491  double xnext;
492 
493  bool done=false;
494  while (done==false) {
495 
496  // Create some references just to make the code easier
497  // to read
498  vec_t &yerr=vtemp;
499  vec_t &dydx_out=vtemp3;
500 
501  // Use ystart as temporary storage for the end of the current
502  // adaptive step
503  vec_t yrow(n);
504  for(size_t i=0;i<n;i++) yrow[i]=y_sol(icurr,i);
505  ret=astp->astep_full(x0,x1,xnext,h,n,yrow,dydx_start,ystart,yerr,
506  dydx_out,derivs);
507 
508  if (ret!=0) {
509  if (exit_on_fail) {
510  n_sol=icurr+1;
511  // call error handler
512  O2SCL_ERR2("Adaptive stepper returned non-zero in ",
513  "ode_iv_solve::solve_store().",exc_efailed);
514  } else if (first_ret==0) {
515  first_ret=ret;
516  }
517  }
518 
519  // Update step counter and abort if necessary
520  nsteps++;
521  if (nsteps>ntrial) {
522  std::string str="Too many steps required (ntrial="+itos(ntrial)+
523  ", x="+o2scl::dtos(x0)+") in ode_iv_solve::solve_store().";
524  O2SCL_ERR(str.c_str(),exc_emaxiter);
525  }
526 
527  // If we've made enough progress, do verbose output
528  if (xnext>=x_verb && verbose>0) {
529  print_iter(xnext,n,ystart);
530  if (verbose>1) {
531  char ch;
532  std::cin >> ch;
533  }
534  x_verb+=dx_verb;
535  }
536 
537  // If we're at the end
538  if (xnext>=x1) {
539 
540  // Exit the loop
541  done=true;
542 
543  // Store the final entry in the table
544  x_sol[icurr]=xnext;
545  for(size_t j=0;j<n;j++) {
546  y_sol(icurr,j)=ystart[j];
547  dydx_sol(icurr,j)=dydx_out[j];
548  yerr_sol(icurr,j)=yerr[j];
549  }
550 
551  // Update the table size
552  n_sol=icurr+1;
553 
554  } else {
555 
556  if (xnext>=x_tab && icurr<nmax) {
557 
558  // If we've made enough progress, store an entry in the table
559  x_sol[icurr]=xnext;
560  for(size_t j=0;j<n;j++) {
561  y_sol(icurr,j)=ystart[j];
562  dydx_sol(icurr,j)=dydx_out[j];
563  yerr_sol(icurr,j)=yerr[j];
564 
565  // Also prepare for the next step
566  y_sol(icurr+1,j)=ystart[j];
567  dydx_start[j]=dydx_out[j];
568  }
569  x_tab+=dx_tab;
570 
571  // Update x0 and the current and next indicies
572  x0=xnext;
573  icurr++;
574 
575  } else {
576 
577  // Otherwise, just prepare for the next step
578  // without updating icurr
579  x0=xnext;
580  for(size_t j=0;j<n;j++) {
581  dydx_start[j]=dydx_out[j];
582  y_sol(icurr,j)=ystart[j];
583  }
584 
585  }
586 
587  }
588 
589  // End of loop 'while (done==false)'
590  }
591 
592  return first_ret;
593  }
594 
595  /** \brief Solve the initial-value problem from \c x0 to \c x1 over
596  a grid storing derivatives and errors
597 
598  Initially, \c xsol should be an array of size \c nsol, and \c
599  ysol should be a \c ubmatrix of size \c [nsol][n]. This function
600  never takes a step larger than the grid size.
601 
602  If \ref verbose is greater than zero, The solution at each grid
603  point will be written to \c std::cout. If \ref verbose is
604  greater than one, a character will be required after each point.
605 
606  \future Consider making a version of grid which gives the
607  same answers as solve_final_value(). After each proposed step,
608  it would go back and fill in the grid points if the proposed
609  step was past the next grid point.
610  */
611  template<class mat_t, class mat_row_t>
612  int solve_grid(double h, size_t n, size_t nsol, vec_t &xsol,
613  mat_t &ysol, mat_t &err_sol, mat_t &dydx_sol,
614  func_t &derivs) {
615 
616  double x0=xsol[0];
617  double x1=xsol[nsol-1];
618 
619  double x=x0, xnext;
620  int ret=0, first_ret=0;
621  nsteps=0;
622 
623  // Copy the initial point to the first row
624  xsol[0]=x0;
625 
626  // Create space for the initial point of each step
627  allocate(n);
628 
629  // Create some references just to make the code easier
630  // to read
631  vec_t &ystart=vtemp;
632  vec_t &dydx_start=vtemp2;
633 
634  // Copy ysol[0] to ystart
635  for(size_t j=0;j<n;j++) {
636  ystart[j]=ysol(0,j);
637  }
638 
639  // Verbose output for the first row
640  if (verbose>0) print_iter(xsol[0],n,ystart);
641 
642  // Evalulate the first derivative
643  derivs(x0,n,ystart,dydx_start);
644 
645  // Set the derivatives and the uncertainties for the first row
646  for(size_t j=0;j<n;j++) {
647  dydx_sol(0,j)=dydx_start[j];
648  err_sol(0,j)=0.0;
649  }
650 
651  for(size_t i=1;i<nsol && ret==0;i++) {
652 
653  mat_row_t y_row(ysol,i);
654  mat_row_t dydx_row(dydx_sol,i);
655  mat_row_t yerr_row(err_sol,i);
656 
657  // Step until we reach the next grid point
658  bool done=false;
659  while(done==false && ret==0) {
660 
661  ret=astp->astep_full(x,xsol[i],xnext,h,n,ystart,dydx_start,
662  y_row,yerr_row,dydx_row,derivs);
663  //ret=gsl_astp2.astep_full(x,xsol[i],xnext,h,n,y_row,yerr_row,
664  //y_row,yerr_row,dydx_row,derivs);
665 
666  nsteps++;
667  if (ret!=0) {
668  if (exit_on_fail) {
669  O2SCL_ERR2("Adaptive stepper failed in ",
670  "ode_iv_solve::solve_grid()",ret);
671  } else if (first_ret!=0) {
672  first_ret=ret;
673  }
674  }
675 
676  if (nsteps>ntrial) {
677  std::string str="Too many steps required (ntrial="+itos(ntrial)+
678  ", x="+o2scl::dtos(x)+") in ode_iv_solve::solve_grid().";
679  O2SCL_ERR(str.c_str(),exc_emaxiter);
680  }
681 
682  // Copy the results of the last step to the starting
683  // point for the next step
684  x=xnext;
685  for(size_t j=0;j<n;j++) {
686  ystart[j]=y_row[j];
687  dydx_start[j]=dydx_row[j];
688  }
689 
690  // Decide if we have reached the grid point
691  if (x1>x0) {
692  if (x>=xsol[i]) done=true;
693  } else {
694  if (x<=xsol[i]) done=true;
695  }
696 
697  }
698 
699  //if (verbose>0) print_iter(xsol[i],n,dydx_row);
700 
701  }
702 
703  return first_ret;
704  }
705  //@}
706 
707  /// Set output level
708  int verbose;
709 
710  /** \brief Number of output points for verbose output (default 10)
711 
712  This is used in functions solve_table() and solve_final_value()
713  to control how often steps are output when verbose is greater
714  than zero.
715  */
716  size_t nsteps_out;
717 
718  /// Maximum number of applications of the adaptive stepper (default 1000)
719  size_t ntrial;
720 
721  /// Number of adaptive ste!ps employed
722  size_t nsteps;
723 
724  /// \name The adaptive stepper
725  //@{
726  /// Set the adaptive stepper to use
728  astp=&as;
729  return 0;
730  }
731 
732  /** \brief If true, stop the solution if the adaptive stepper fails
733  (default true)
734 
735  If this is false, then failures in the adaptive stepper are
736  ignored.
737  */
739 
740  /// The default adaptive stepper
742  //@}
743 
744  /// Return the type, \c "ode_iv_solve".
745  virtual const char *type() { return "ode_iv_solve"; }
746 
747  };
748 
749 #ifndef DOXYGEN_NO_O2NS
750 }
751 #endif
752 
753 #endif
Solve an initial-value ODE problems given an adaptive ODE stepper.
Definition: ode_iv_solve.h:91
astep_base< vec_t, vec_t, vec_t, func_t > * astp
The adaptive stepper.
Definition: ode_iv_solve.h:110
The main O<span style=&#39;position: relative; top: 0.3em; font-size: 0.8em&#39;>2</span>scl O$_2$scl names...
Definition: anneal.h:42
astep_gsl< vec_t, vec_t, vec_t, func_t > gsl_astp
The default adaptive stepper.
Definition: ode_iv_solve.h:741
invalid argument supplied by user
Definition: err_hnd.h:59
virtual const char * type()
Return the type, "ode_iv_solve".
Definition: ode_iv_solve.h:745
size_t nsteps
Number of adaptive ste!ps employed.
Definition: ode_iv_solve.h:722
exceeded max number of iterations
Definition: err_hnd.h:73
bool exit_on_fail
If true, stop the solution if the adaptive stepper fails (default true)
Definition: ode_iv_solve.h:738
int verbose
Set output level.
Definition: ode_iv_solve.h:708
generic failure
Definition: err_hnd.h:61
int solve_final_value(double x0, double x1, double h, size_t n, vec_t &ystart, vec_t &yend, vec_t &yerr, func_t &derivs)
Solve the initial-value problem to get the final value with errors.
Definition: ode_iv_solve.h:204
int solve_final_value(double x0, double x1, double h, size_t n, vec_t &ystart, vec_t &yend, vec_t &yerr, vec_t &dydx_end, func_t &derivs)
Solve the initial-value problem to get the final value, derivatives, and errors.
Definition: ode_iv_solve.h:239
int solve_final_value(double x0, double x1, double h, size_t n, vec_t &ystart, vec_t &yend, func_t &derivs)
Solve the initial-value problem to get the final value.
Definition: ode_iv_solve.h:181
#define O2SCL_ERR2(d, d2, n)
Set an error, two-string version.
Definition: err_hnd.h:281
std::string dtos(double x, int prec=6, bool auto_prec=false)
Convert a double to a string.
#define O2SCL_ERR(d, n)
Set an error with message d and code n.
Definition: err_hnd.h:273
virtual int astep_full(double x, double xlimit, double &x_out, double &h, size_t n, vec_y_t &y, vec_dydx_t &dydx, vec_y_t &yout, vec_yerr_t &yerr, vec_dydx_t &dydx_out, func_t &derivs)=0
Make an adaptive integration step of the system derivs with derivatives.
size_t ntrial
Maximum number of applications of the adaptive stepper (default 1000)
Definition: ode_iv_solve.h:719
int set_astep(astep_base< vec_t, vec_t, vec_t, func_t > &as)
Set the adaptive stepper to use.
Definition: ode_iv_solve.h:727
int solve_store(double x0, double x1, double h, size_t n, size_t &n_sol, vec_t &x_sol, mat_t &y_sol, mat_t &yerr_sol, mat_t &dydx_sol, func_t &derivs, size_t istart=0)
Solve the initial-value problem and store the associated output.
Definition: ode_iv_solve.h:438
void allocate(size_t n)
Allocate space for temporary vectors.
Definition: ode_iv_solve.h:136
virtual int print_iter(double x, size_t nv, vec_t &y)
Print out iteration information.
Definition: ode_iv_solve.h:113
void free()
Free allocated memory.
Definition: ode_iv_solve.h:125
int solve_grid(double h, size_t n, size_t nsol, vec_t &xsol, mat_t &ysol, mat_t &err_sol, mat_t &dydx_sol, func_t &derivs)
Solve the initial-value problem from x0 to x1 over a grid storing derivatives and errors...
Definition: ode_iv_solve.h:612
size_t nsteps_out
Number of output points for verbose output (default 10)
Definition: ode_iv_solve.h:716
size_t mem_size
Desc.
Definition: ode_iv_solve.h:107
std::string itos(int x)
Convert an integer to a string.
static const double x1[5]
Definition: inte_qng_gsl.h:48
std::function< int(double, size_t, const boost::numeric::ublas::vector< double > &, boost::numeric::ublas::vector< double > &)> ode_funct11
Ordinary differential equation function.
Definition: ode_funct.h:46

Documentation generated with Doxygen. Provided under the GNU Free Documentation License (see License Information).