00001 /* 00002 ------------------------------------------------------------------- 00003 00004 Copyright (C) 2006, 2007, 2008, 2009, Andrew W. Steiner 00005 00006 This file is part of O2scl. 00007 00008 O2scl is free software; you can redistribute it and/or modify 00009 it under the terms of the GNU General Public License as published by 00010 the Free Software Foundation; either version 3 of the License, or 00011 (at your option) any later version. 00012 00013 O2scl is distributed in the hope that it will be useful, 00014 but WITHOUT ANY WARRANTY; without even the implied warranty of 00015 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00016 GNU General Public License for more details. 00017 00018 You should have received a copy of the GNU General Public License 00019 along with O2scl. If not, see <http://www.gnu.org/licenses/>. 00020 00021 ------------------------------------------------------------------- 00022 */ 00023 #ifndef O2SCL_ODE_IV_SOLVE_H 00024 #define O2SCL_ODE_IV_SOLVE_H 00025 00026 #include <string> 00027 #include <o2scl/collection.h> 00028 #include <o2scl/adapt_step.h> 00029 #include <o2scl/gsl_astep.h> 00030 00031 #ifndef DOXYGENP 00032 namespace o2scl { 00033 #endif 00034 00035 /** 00036 \brief Solve an initial-value ODE problems given an 00037 adaptive ODE stepper 00038 00039 \future Add error information 00040 */ 00041 template<class param_t, class func_t=ode_funct<param_t>, 00042 class vec_t=ovector_base, class alloc_vec_t=ovector, 00043 class alloc_t=ovector_alloc, class mat_row_t=omatrix_row> 00044 class ode_iv_solve { 00045 00046 #ifndef DOXYGEN_INTERNAL 00047 00048 protected: 00049 00050 /// Derivative 00051 alloc_vec_t dydx; 00052 00053 /// Memory allocator 00054 alloc_t ao; 00055 00056 /// The adaptive stepper 00057 adapt_step<param_t,func_t,vec_t,alloc_vec_t,alloc_t> *astp; 00058 00059 /// Print out iteration information 00060 virtual int print_iter(double x, size_t nv, vec_t &y) { 00061 std::cout << type() << " x: " << x << " y: "; 00062 for(size_t i=0;i<nv;i++) std::cout << y[i] << " "; 00063 std::cout << std::endl; 00064 if (verbose>1) { 00065 char ch; 00066 std::cin >> ch; 00067 } 00068 return 0; 00069 } 00070 00071 #endif 00072 00073 public: 00074 00075 ode_iv_solve() { 00076 verbose=0; 00077 ntrial=1000; 00078 nsteps_out=10; 00079 astp=&gsl_astp; 00080 exit_on_fail=true; 00081 } 00082 00083 virtual ~ode_iv_solve() {} 00084 00085 /// Set the adaptive stepper to use 00086 int set_adapt_step(adapt_step<param_t,func_t,vec_t, 00087 alloc_vec_t,alloc_t> &as) { 00088 astp=&as; 00089 return 0; 00090 } 00091 00092 /** 00093 \brief Solve the initial-value problem and output a table 00094 00095 Initially, \c xsol should be a vector of size \c nsol, and \c 00096 ysol should be a two-dimensional array (i.e. omatrix_view) of 00097 size \c [nsol][n]. On exit, \c nsol will will be the size of 00098 the solution table, less than or equal to the original value 00099 of \c nsol. 00100 00101 If \ref verbose is greater than zero, The solution 00102 at each internal point will be written to \c std::cout. 00103 If \ref verbose is greater than one, a character 00104 will be required after each point. 00105 00106 If the given value of \c h is small enough, the solution 00107 may generate more points than the space initially allocated 00108 and the full solution will not be generated. 00109 00110 \future Consider modifying so that this can handle tables 00111 which are too small by removing half the rows and doubling 00112 the stepsize. 00113 00114 */ 00115 template<class mat_t> 00116 int solve_table(double x0, double x1, double h, size_t n, 00117 vec_t &ystart, size_t &nsol, vec_t &xsol, 00118 mat_t &ysol, param_t &pa, func_t &derivs) { 00119 00120 int ret=0, nmax=nsol, i, first_ret=0; 00121 size_t j; 00122 00123 xsol[0]=x0; 00124 for(j=0;j<n;j++) ysol[0][j]=ystart[j]; 00125 if (verbose>0) { 00126 print_iter(xsol[0],n,ystart); 00127 if (verbose>1) { 00128 char ch; 00129 std::cin >> ch; 00130 } 00131 } 00132 00133 alloc_vec_t yerrl; 00134 ao.allocate(yerrl,n); 00135 ao.allocate(dydx,n); 00136 derivs(x0,n,ystart,dydx,pa); 00137 00138 for(i=1;i<nmax && xsol[i-1]<x1;i++) { 00139 00140 xsol[i]=xsol[i-1]; 00141 for(j=0;j<n;j++) ysol[i][j]=ysol[i-1][j]; 00142 00143 mat_row_t ar(ysol,i); 00144 00145 ret=astp->astep_derivs(xsol[i],h,x1,n,ar,dydx,yerrl,pa,derivs); 00146 if (ret!=0) { 00147 if (exit_on_fail) { 00148 nsol=i+1; 00149 O2SCL_ERR2_RET("Adaptive stepper failed in ", 00150 "ode_iv_solve::solve_table()",ret); 00151 } else if (first_ret!=0) { 00152 first_ret=ret; 00153 } 00154 } 00155 00156 if (verbose>0) { 00157 print_iter(xsol[i],n,ar); 00158 if (verbose>1) { 00159 char ch; 00160 std::cin >> ch; 00161 } 00162 } 00163 } 00164 00165 ao.free(dydx); 00166 ao.free(yerrl); 00167 00168 nsol=i; 00169 00170 if (i==nmax && xsol[i-1]<x1) { 00171 O2SCL_ERR2_RET("Ran out of space in ", 00172 "ode_iv_solve::solve_table().",gsl_etable); 00173 } 00174 return first_ret; 00175 } 00176 00177 /** \brief Solve the initial-value problem from \c x0 00178 to \c x1 over a grid 00179 00180 Initially, \c xsol should be an array of size \c nsol, and 00181 \c ysol should be a \c matrix of size \c [nsol][n]. This 00182 function never takes a step larger than the grid size. 00183 00184 If \ref verbose is greater than zero, The at each grid point 00185 will be written to \c std::cout. If \ref verbose is greater 00186 than one, a character will be required after each point. 00187 */ 00188 template<class mat_t> 00189 int solve_grid(double x0, double x1, double h, size_t n, 00190 vec_t &ystart, size_t nsol, vec_t &xsol, mat_t &ysol, 00191 param_t &pa, func_t &derivs) { 00192 00193 double x=x0, xnext; 00194 int ret=0, first_ret=0; 00195 size_t j; 00196 00197 xsol[0]=x0; 00198 for(j=0;j<n;j++) ysol[0][j]=ystart[j]; 00199 if (verbose>0) print_iter(xsol[0],n,ystart); 00200 00201 alloc_vec_t yerrl; 00202 ao.allocate(dydx,n); 00203 ao.allocate(yerrl,n); 00204 00205 if (x0<x1) { 00206 00207 derivs(x0,n,ystart,dydx,pa); 00208 00209 for(size_t i=1;i<nsol && ret==0;i++) { 00210 00211 xsol[i]=xsol[i-1]; 00212 for(j=0;j<n;j++) ysol[i][j]=ysol[i-1][j]; 00213 00214 xnext=x0+(x1-x0)*((double)i)/((double)(nsol-1)); 00215 00216 mat_row_t ar(ysol,i); 00217 while(xsol[i]<xnext && ret==0) { 00218 00219 ret=astp->astep_derivs(xsol[i],h,xnext,n,ar,dydx,yerrl,pa,derivs); 00220 if (ret!=0) { 00221 if (exit_on_fail) { 00222 O2SCL_ERR2_RET("Adaptive stepper failed in ", 00223 "ode_iv_solve::solve_grid()",ret); 00224 } else if (first_ret!=0) { 00225 first_ret=ret; 00226 } 00227 } 00228 } 00229 if (verbose>0) print_iter(xsol[i],n,ar); 00230 00231 } 00232 00233 } else { 00234 00235 derivs(x0,n,ystart,dydx,pa); 00236 00237 for(size_t i=1;i<nsol && ret==0;i++) { 00238 00239 xsol[i]=xsol[i-1]; 00240 for(j=0;j<n;j++) ysol[i][j]=ysol[i-1][j]; 00241 00242 xnext=x0+(x1-x0)*((double)i)/((double)(nsol-1)); 00243 00244 mat_row_t ar(ysol,i); 00245 while(xsol[i]>xnext && ret==0) { 00246 ret=astp->astep_derivs(xsol[i],h,xnext,n,ar,dydx,yerrl,pa,derivs); 00247 if (ret!=0) { 00248 if (exit_on_fail) { 00249 O2SCL_ERR2_RET("Adaptive stepper failed in ", 00250 "ode_iv_solve::solve_grid()",ret); 00251 } else if (first_ret!=0) { 00252 first_ret=ret; 00253 } 00254 } 00255 } 00256 if (verbose>0) print_iter(xsol[i],n,ar); 00257 00258 } 00259 00260 } 00261 00262 ao.free(dydx); 00263 ao.free(yerrl); 00264 00265 return first_ret; 00266 } 00267 00268 /** \brief Solve the initial-value problem from \c x0 00269 to \c x1 over a grid storing derivatives 00270 00271 Initially, \c xsol should be an array of size \c nsol, and 00272 \c ysol should be a \c omatrix of size \c [nsol][n]. This 00273 function never takes a step larger than the grid size. 00274 00275 If \ref verbose is greater than zero, The solution 00276 at each grid point will be written to \c std::cout. 00277 If \ref verbose is greater than one, a character 00278 will be required after each point. 00279 00280 \todo Add error information 00281 */ 00282 template<class mat_t> 00283 int solve_grid_derivs(double x0, double x1, double h, size_t n, 00284 vec_t &ystart, size_t nsol, vec_t &xsol, 00285 mat_t &ysol, mat_t &dydx_sol, 00286 param_t &pa, func_t &derivs) { 00287 00288 double x=x0, xnext; 00289 int ret=0, first_ret=0; 00290 size_t j; 00291 00292 xsol[0]=x0; 00293 for(j=0;j<n;j++) ysol[0][j]=ystart[j]; 00294 if (verbose>0) print_iter(xsol[0],n,ystart); 00295 00296 alloc_vec_t yerrl; 00297 ao.allocate(yerrl,n); 00298 00299 if (x0<x1) { 00300 00301 mat_row_t mr0(dydx_sol,0); 00302 derivs(x0,n,ystart,mr0,pa); 00303 00304 for(size_t i=1;i<nsol && ret==0;i++) { 00305 00306 xsol[i]=xsol[i-1]; 00307 for(j=0;j<n;j++) { 00308 ysol[i][j]=ysol[i-1][j]; 00309 dydx_sol[i][j]=dydx_sol[i-1][j]; 00310 } 00311 00312 xnext=x0+(x1-x0)*((double)i)/((double)(nsol-1)); 00313 00314 mat_row_t ar(ysol,i); 00315 mat_row_t ar2(dydx_sol,i); 00316 00317 while(xsol[i]<xnext && ret==0) { 00318 00319 ret=astp->astep_derivs(xsol[i],h,xnext,n,ar,ar2,yerrl,pa,derivs); 00320 if (ret!=0) { 00321 if (exit_on_fail) { 00322 ao.free(yerrl); 00323 O2SCL_ERR2_RET("Adaptive stepper failed in ", 00324 "ode_iv_solve::solve_grid()",ret); 00325 } else if (first_ret!=0) { 00326 first_ret=ret; 00327 } 00328 } 00329 } 00330 if (verbose>0) print_iter(xsol[i],n,ar); 00331 00332 } 00333 00334 } else { 00335 00336 mat_row_t mr0(dydx_sol,0); 00337 derivs(x0,n,ystart,dydx,pa); 00338 00339 for(size_t i=1;i<nsol && ret==0;i++) { 00340 00341 xsol[i]=xsol[i-1]; 00342 for(j=0;j<n;j++) { 00343 ysol[i][j]=ysol[i-1][j]; 00344 dydx_sol[i][j]=dydx_sol[i-1][j]; 00345 } 00346 00347 xnext=x0+(x1-x0)*((double)i)/((double)(nsol-1)); 00348 00349 mat_row_t ar(ysol,i); 00350 mat_row_t ar2(dydx_sol,i); 00351 00352 while(xsol[i]>xnext && ret==0) { 00353 ret=astp->astep_derivs(xsol[i],h,xnext,n,ar,ar2,yerrl,pa,derivs); 00354 if (ret!=0) { 00355 if (exit_on_fail) { 00356 ao.free(yerrl); 00357 O2SCL_ERR2_RET("Adaptive stepper failed in ", 00358 "ode_iv_solve::solve_grid()",ret); 00359 } else if (first_ret!=0) { 00360 first_ret=ret; 00361 } 00362 } 00363 } 00364 if (verbose>0) print_iter(xsol[i],n,ar); 00365 00366 } 00367 00368 } 00369 00370 ao.free(yerrl); 00371 00372 return first_ret; 00373 } 00374 00375 #ifdef O2SCL_NEVER_DEFINED 00376 00377 /** \brief Solve the initial-value problem from \c x0 00378 to \c x1 over a grid storing derivatives 00379 00380 Initially, \c xsol should be an array of size \c nsol, and 00381 \c ysol should be a \c omatrix of size \c [nsol][n]. This 00382 function never takes a step larger than the grid size. 00383 00384 If \ref verbose is greater than zero, The solution 00385 at each grid point will be written to \c std::cout. 00386 If \ref verbose is greater than one, a character 00387 will be required after each point. 00388 00389 \todo Add error information 00390 */ 00391 template<class mat_t> 00392 int solve_grid_derivs_err(double x0, double x1, double h, size_t n, 00393 vec_t &ystart, size_t nsol, vec_t &xsol, 00394 mat_t &ysol, mat_t &dydx_sol, mat_t &errsol, 00395 param_t &pa, func_t &derivs) { 00396 00397 double x=x0, xnext; 00398 int ret=0, first_ret=0; 00399 size_t j; 00400 00401 xsol[0]=x0; 00402 for(j=0;j<n;j++) ysol[0][j]=ystart[j]; 00403 if (verbose>0) print_iter(xsol[0],n,ystart); 00404 00405 if (x0<x1) { 00406 00407 mat_row_t mr0(dydx_sol,0); 00408 derivs(x0,n,ystart,mr0,pa); 00409 00410 for(size_t i=1;i<nsol && ret==0;i++) { 00411 00412 xsol[i]=xsol[i-1]; 00413 for(j=0;j<n;j++) { 00414 ysol[i][j]=ysol[i-1][j]; 00415 dydx_sol[i][j]=dydx_sol[i-1][j]; 00416 } 00417 00418 xnext=x0+(x1-x0)*((double)i)/((double)(nsol-1)); 00419 00420 mat_row_t ar(ysol,i); 00421 mat_row_t ar2(dydx_sol,i); 00422 mat_row_t aerr(errsol,i); 00423 00424 while(xsol[i]<xnext && ret==0) { 00425 00426 ret=astp->astep_derivs(xsol[i],h,xnext,n,ar,ar2,aerr,pa,derivs); 00427 if (ret!=0) { 00428 if (exit_on_fail) { 00429 O2SCL_ERR2_RET("Adaptive stepper failed in ", 00430 "ode_iv_solve::solve_grid()",ret); 00431 } else if (first_ret!=0) { 00432 first_ret=ret; 00433 } 00434 } 00435 } 00436 if (verbose>0) print_iter(xsol[i],n,ar); 00437 00438 } 00439 00440 } else { 00441 00442 mat_row_t mr0(dydx_sol,0); 00443 derivs(x0,n,ystart,dydx,pa); 00444 00445 for(size_t i=1;i<nsol && ret==0;i++) { 00446 00447 xsol[i]=xsol[i-1]; 00448 for(j=0;j<n;j++) { 00449 ysol[i][j]=ysol[i-1][j]; 00450 dydx_sol[i][j]=dydx_sol[i-1][j]; 00451 } 00452 00453 xnext=x0+(x1-x0)*((double)i)/((double)(nsol-1)); 00454 00455 mat_row_t ar(ysol,i); 00456 mat_row_t ar2(dydx_sol,i); 00457 mat_row_t aerr(errsol,i); 00458 00459 while(xsol[i]>xnext && ret==0) { 00460 ret=astp->astep_derivs(xsol[i],h,xnext,n,ar,ar2,aerr,pa,derivs); 00461 if (ret!=0) { 00462 if (exit_on_fail) { 00463 O2SCL_ERR2_RET("Adaptive stepper failed in ", 00464 "ode_iv_solve::solve_grid()",ret); 00465 } else if (first_ret!=0) { 00466 first_ret=ret; 00467 } 00468 } 00469 } 00470 if (verbose>0) print_iter(xsol[i],n,ar); 00471 00472 } 00473 00474 } 00475 00476 return first_ret; 00477 } 00478 00479 #endif 00480 00481 /** \brief Solve the initial-value problem 00482 to get the final value 00483 00484 If \ref verbose is greater than zero, The solution at less 00485 than or approximately equal to \ref nsteps_out points will be 00486 written to \c std::cout. If \ref verbose is greater than one, 00487 a character will be required after each selected point. 00488 00489 The solution fails if more than \ref ntrial steps are 00490 required. 00491 */ 00492 int solve_final_value(double x0, double x1, double h, size_t n, 00493 vec_t &ystart, vec_t ¥d, 00494 param_t &pa, func_t &derivs) { 00495 00496 ao.allocate(dydx,n); 00497 int ret=solve_final_value_derivs(x0,x1,h,n,ystart,yend,dydx, 00498 dydx,pa,derivs); 00499 ao.free(dydx); 00500 return ret; 00501 } 00502 00503 /** \brief Solve the initial-value problem 00504 to get the final value and derivative 00505 00506 If \ref verbose is greater than zero, The solution at less 00507 than or approximately equal to \ref nsteps_out points will be 00508 written to \c std::cout. If \ref verbose is greater than one, 00509 a character will be required after each selected point. 00510 00511 The solution fails if more than \ref ntrial steps are 00512 required. 00513 00514 \todo Add error information 00515 */ 00516 int solve_final_value_derivs(double x0, double x1, double h, size_t n, 00517 vec_t &ystart, vec_t ¥d, 00518 vec_t &dydx_start, vec_t &dydx_end, 00519 param_t &pa, func_t &derivs) { 00520 00521 if ((x1>x0 && h<=0.0) || (x0>x1 && h>=0.0)) { 00522 O2SCL_ERR_RET("Interval ordering and sign of h don't match.", 00523 gsl_einval); 00524 } 00525 00526 double x=x0, xverb=0.0, dxverb=0.0; 00527 int ret=0, it=0, first_ret=0; 00528 if (verbose>0) { 00529 print_iter(x0,n,ystart); 00530 if (verbose>1) { 00531 char ch; 00532 std::cin >> ch; 00533 } 00534 if (x1>x0) { 00535 dxverb=(x1-x0)/((double)nsteps_out); 00536 xverb=x0+dxverb; 00537 } else { 00538 dxverb=(x0-x1)/((double)nsteps_out); 00539 xverb=x0-dxverb; 00540 } 00541 } 00542 00543 for(size_t i=0;i<n;i++) yend[i]=ystart[i]; 00544 00545 alloc_vec_t yerrl; 00546 ao.allocate(yerrl,n); 00547 00548 if (x1>x0) { 00549 00550 derivs(x,n,yend,dydx_start,pa); 00551 for(size_t i=0;i<n;i++) dydx_end[i]=dydx_start[i]; 00552 00553 while(x<x1 && ret==0) { 00554 ret=astp->astep_derivs(x,h,x1,n,yend,dydx_end,yerrl,pa,derivs); 00555 if (ret!=0) { 00556 if (exit_on_fail) { 00557 ao.free(yerrl); 00558 O2SCL_ERR2_RET("Adaptive stepper failed in ", 00559 "ode_iv_solve::solve_grid()",ret); 00560 } else if (first_ret!=0) { 00561 first_ret=ret; 00562 } 00563 } 00564 if (verbose>0 && x>xverb) { 00565 print_iter(x,n,yend); 00566 if (verbose>1) { 00567 char ch; 00568 std::cin >> ch; 00569 } 00570 xverb+=dxverb; 00571 } 00572 it++; 00573 if (it>ntrial) { 00574 ao.free(yerrl); 00575 O2SCL_ERR2_RET("Exceeded max number of iterations in ", 00576 "ode_iv_solve::solve_final_value().", 00577 gsl_emaxiter); 00578 } 00579 } 00580 00581 } else { 00582 00583 derivs(x,n,yend,dydx_start,pa); 00584 for(size_t i=0;i<n;i++) dydx_end[i]=dydx_start[i]; 00585 00586 while(x>x1 && ret==0) { 00587 ret=astp->astep_derivs(x,h,x1,n,yend,dydx_end,yerrl,pa,derivs); 00588 if (ret!=0) { 00589 if (exit_on_fail) { 00590 ao.free(yerrl); 00591 O2SCL_ERR2_RET 00592 ("Adaptive stepper failed in ", 00593 "ode_iv_solve::solve_grid()",ret); 00594 } else if (first_ret!=0) { 00595 first_ret=ret; 00596 } 00597 } 00598 if (verbose>0 && x<xverb) { 00599 print_iter(x,n,yend); 00600 if (verbose>1) { 00601 char ch; 00602 std::cin >> ch; 00603 } 00604 xverb-=dxverb; 00605 } 00606 it++; 00607 if (it>ntrial) { 00608 ao.free(yerrl); 00609 O2SCL_ERR2_RET("Exceeded max number of iterations in ", 00610 "ode_iv_solve::solve_final_value().", 00611 gsl_emaxiter); 00612 } 00613 } 00614 } 00615 00616 if (verbose>0) { 00617 print_iter(x,n,yend); 00618 if (verbose>1) { 00619 char ch; 00620 std::cin >> ch; 00621 } 00622 } 00623 00624 ao.free(yerrl); 00625 return first_ret; 00626 } 00627 00628 /// Set output level 00629 int verbose; 00630 00631 /** 00632 \brief Number of output points if \ref verbose is greater 00633 than zero (default 10) 00634 */ 00635 size_t nsteps_out; 00636 00637 /// Maximum number of steps for solve_final_value() (default 1000) 00638 int ntrial; 00639 00640 /** 00641 \brief If true, stop the solution if the adaptive stepper fails 00642 00643 If this is false, then failures in the adaptive stepper are 00644 ignored. 00645 */ 00646 bool exit_on_fail; 00647 00648 /// The default adaptive stepper 00649 gsl_astep<param_t,func_t,vec_t,alloc_vec_t,alloc_t> gsl_astp; 00650 00651 /// Return the type, \c "ode_iv_solve". 00652 virtual const char *type() { return "ode_iv_solve"; } 00653 00654 }; 00655 00656 #ifdef O2SCL_NEVER_DEFINED 00657 00658 /** 00659 \brief Solve an initial-value ODE problems given an 00660 adaptive ODE stepper 00661 */ 00662 template<class param_t,func_t=ode_funct<param_t> > 00663 class ode_iv_solve_x : public ode_iv_solve<param_t,func_t, 00664 ovector_base,ovector,ovector_alloc,omatrix_row> { 00665 00666 public: 00667 00668 int solve_x_derivs(double x0, double x1, double h, size_t n, 00669 ovector_base &ystart, size_t &nsol, 00670 ovector &xsol, std::vector<ovector> &ysol, 00671 std::vector<ovector> &dydx_sol, 00672 param_t &pa, func_t &derivs) { 00673 00674 // Temporary ovectors 00675 ovector y, dy; 00676 00677 // Proceed with solution 00678 double x=x0, xnext; 00679 int ret=0, first_ret=0; 00680 size_t j; 00681 00682 xsol.push_back(x); 00683 for(j=0;j<n;j++) y[j]=ystart[j]; 00684 ysol.push_back(y); 00685 00686 alloc_vec_t yerrl; 00687 ao.allocate(yerrl,n); 00688 00689 if (x0<x1) { 00690 00691 derivs(x,n,y,dy,pa); 00692 dydx_sol.push_back(dy); 00693 00694 while (x<x1 && ret==0) { 00695 00696 ret=astp->astep_derivs(x,h,x1,n,y,dy,yerrl,pa,derivs); 00697 if (ret!=0) { 00698 if (exit_on_fail) { 00699 ao.free(yerrl); 00700 O2SCL_ERR2_RET("Adaptive stepper failed in ", 00701 "ode_iv_solve::solve_grid()",ret); 00702 } else if (first_ret!=0) { 00703 first_ret=ret; 00704 } 00705 } 00706 if (verbose>0) print_iter(xsol[i],n,ar); 00707 00708 xsol.push_back(x); 00709 ysol.push_back(y); 00710 ysol.push_back(dy); 00711 00712 00713 } 00714 00715 } 00716 00717 ao.free(yerrl); 00718 00719 return first_ret; 00720 } 00721 00722 #endif 00723 00724 00725 #ifndef DOXYGENP 00726 } 00727 #endif 00728 00729 #endif
Documentation generated with Doxygen and provided under the GNU Free Documentation License. See License Information for details.
Project hosting provided by
,
O2scl Sourceforge Project Page