diff options
Diffstat (limited to 'src/c/differential_calculus')
| -rw-r--r-- | src/c/differential_calculus/includes/ode.h | 6 | ||||
| -rw-r--r-- | src/c/differential_calculus/interfaces/int_ode.h | 31 | ||||
| -rw-r--r-- | src/c/differential_calculus/ode/dodea.c | 32 | ||||
| -rw-r--r-- | src/c/differential_calculus/ode/dodes.c | 44 |
4 files changed, 83 insertions, 30 deletions
diff --git a/src/c/differential_calculus/includes/ode.h b/src/c/differential_calculus/includes/ode.h index 0a870201..79962434 100644 --- a/src/c/differential_calculus/includes/ode.h +++ b/src/c/differential_calculus/includes/ode.h @@ -14,11 +14,11 @@ #define __ODE_H__ double dodes(double initial_value, double start_time, double end_time, \ - char *ode_function, double nequs, double eps_abs, double eps_rel, \ - double step_size); + int (*ode_function), char *solver_type, double nequs, double eps_abs, double eps_rel, \ + double step_size, int *params); void dodea(double *initial_value, double start_time, double end_time, \ - char *ode_function, double nequs, double eps_abs, double eps_rel, \ + int (*ode_function), char *solver_type, double nequs, double eps_abs, double eps_rel, \ double step_size, int *params, double *out); #endif /*__ODE_H__*/
\ No newline at end of file diff --git a/src/c/differential_calculus/interfaces/int_ode.h b/src/c/differential_calculus/interfaces/int_ode.h index 1cb64497..28f43999 100644 --- a/src/c/differential_calculus/interfaces/int_ode.h +++ b/src/c/differential_calculus/interfaces/int_ode.h @@ -17,11 +17,34 @@ extern "C" { #endif -#define d0d0d0g2oded0(in1, in2, in3, in4) dodes(in1, in2, in3, in4, 1, 1.0e-2,\ - 1.0e-2, 1.0e-6) +#define d0d0d0fn0oded0(in1, in2, in3, func_name) dodes(in1, in2, in3, \ + func_name, "rkf",1, 1.0e-2, 1.0e-2, 1.0e-6, NULL) -#define d2d0d0f0oded2(in1, size1, in2, in3, func_name, out) dodea(in1, in2, in3, func_name, \ - size1[1], 1.0e-2, 1.0e-2, 1.0e-6, size1, out) +#define d2d0d0fn0oded2(in1, size1, in2, in3, func_name, out) dodea(in1, \ + in2, in3, func_name, "rkf",size1[1], 1.0e-2, 1.0e-2, \ + 1.0e-6, size1, out) + +#define d0d0d2fn0oded2(in1, in2, in3, size3, func_name, out) dodea(in1, in2, \ + in3, func_name, "rkf", 1, 1.0e-2, 1.0e-2, 1.0e-6, size3, out) + +#define d2d0d2fn0oded2(in1, size1, in2, in3, size3, func_name, out) dodea(in1, \ + in2, in3, func_name, "rkf",size1[1], 1.0e-2, 1.0e-2, \ + 1.0e-6, size1, out) + +#define g2d0d0d0fn0oded0(solvertype, typesize, in1, in2, in3, func_name) dodes(in1, in2, in3, \ + func_name, solvertype, 1, 1.0e-2, 1.0e-2, 1.0e-6, NULL) + +#define g2d2d0d0fn0oded2(solvertype, typesize, in1, size1, in2, in3, func_name, out) \ + dodea(in1, in2, in3, func_name, solvertype, size1[1], \ + 1.0e-2, 1.0e-2, 1.0e-6, size1, out) + +#define g2d0d0d2fn0oded2(solvertype, typesize, in1, in2, in3, size3, func_name, out) \ + dodea(in1, in2, in3, func_name, solvertype,1, 1.0e-2, 1.0e-2, \ + 1.0e-6, size3, out) + +#define g2d2d0d2fn0oded2(solvertype, typesize, in1, size1, in2, in3, size3, func_name, out) dodea(in1, \ + in2, in3, func_name, solvertype, size1[1], 1.0e-2, 1.0e-2, \ + 1.0e-6, size1, out) #ifdef __cplusplus } /* extern "C" */ diff --git a/src/c/differential_calculus/ode/dodea.c b/src/c/differential_calculus/ode/dodea.c index ff3cbdec..1cb07fab 100644 --- a/src/c/differential_calculus/ode/dodea.c +++ b/src/c/differential_calculus/ode/dodea.c @@ -10,7 +10,7 @@ Email: toolbox@scilab.in */ -//Function for solving ODEs using GSL library +/*Function for solving ODEs using GSL library*/ #include "ode.h" #include "types.h" @@ -20,21 +20,39 @@ void dodea(double *initial_value, double start_time, double end_time, \ - char *ode_function, double nequs, double eps_abs, double eps_rel, \ - double step_size, int *params, double *out) + int (*ode_function), char *solver_type, double nequs, double eps_abs, \ + double eps_rel, double step_size, int *params, double *out) { double t = start_time; - //Initialise output to initial state + gsl_odeiv2_step_type *step_type; + + /*Initialise output to initial state*/ int counter = 0; for (counter = 0; counter<nequs;counter++) { out[counter] = initial_value[counter]; } - //Setup ODE related parameters - gsl_odeiv2_system sys = {ode_function, NULL, 2, params}; + /*Setup ODE related parameters*/ + gsl_odeiv2_system sys = {ode_function, NULL, nequs, params}; - gsl_odeiv2_step *s = gsl_odeiv2_step_alloc (gsl_odeiv2_step_rkf45, nequs); + /*Select step solver*/ + if (solver_type == "adams") + step_type = gsl_odeiv2_step_msadams; + if (solver_type == "stiff") + step_type = gsl_odeiv2_step_msbdf; + if (solver_type == "rk") + step_type = gsl_odeiv2_step_rk4; + if (solver_type == "rkf") + step_type = gsl_odeiv2_step_rkf45; + if (solver_type == "root") + step_type = gsl_odeiv2_step_rkck; + if (solver_type == "discrete") + step_type = gsl_odeiv2_step_rk8pd; + else + step_type = gsl_odeiv2_step_rkf45; + + gsl_odeiv2_step *s = gsl_odeiv2_step_alloc (step_type, nequs); gsl_odeiv2_control *c = gsl_odeiv2_control_y_new (eps_abs, eps_rel); gsl_odeiv2_evolve *e = gsl_odeiv2_evolve_alloc (nequs); diff --git a/src/c/differential_calculus/ode/dodes.c b/src/c/differential_calculus/ode/dodes.c index ee6718bb..adef1ba7 100644 --- a/src/c/differential_calculus/ode/dodes.c +++ b/src/c/differential_calculus/ode/dodes.c @@ -1,13 +1,16 @@ -// Copyright (C) 2016 - IIT Bombay - FOSSEE -// -// This file must be used under the terms of the CeCILL. -// This source file is licensed as described in the file COPYING, which -// you should have received as part of this distribution. The terms -// are also available at -// http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt -// Author: Siddhesh Wani -// Organization: FOSSEE, IIT Bombay -// Email: toolbox@scilab.in +/* Copyright (C) 2016 - IIT Bombay - FOSSEE + + This file must be used under the terms of the CeCILL. + This source file is licensed as described in the file COPYING, which + you should have received as part of this distribution. The terms + are also available at + http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt + Author: Siddhesh Wani + Organization: FOSSEE, IIT Bombay + Email: toolbox@scilab.in +*/ + +/*Function for solving ODEs using GSL library*/ #include "ode.h" #include "types.h" @@ -17,20 +20,29 @@ double dodes(double initial_value, double start_time, double end_time, \ - char *ode_function, double nequs, double eps_abs, double eps_rel, \ - double step_size) + int (*ode_function), char *solver_type, double nequs, double eps_abs, \ + double eps_rel, double step_size, int *params) { - double out = 0; + double out = 0, t = 0; //int status; + out = initial_value; + t = start_time; //Setup ODE related parameters - gsl_odeiv2_system sys = {ode_function, NULL, 2, NULL}; + gsl_odeiv2_system sys = {ode_function, NULL, nequs, NULL}; gsl_odeiv2_step *s = gsl_odeiv2_step_alloc (gsl_odeiv2_step_rkf45, nequs); gsl_odeiv2_control *c = gsl_odeiv2_control_y_new (eps_abs, eps_rel); gsl_odeiv2_evolve *e = gsl_odeiv2_evolve_alloc (nequs); - gsl_odeiv2_evolve_apply_fixed_step (e, c, s, &sys, &start_time, step_size, &out); - + while(t < end_time) + { + gsl_odeiv2_evolve_apply_fixed_step (e, c, s, &sys, &t, step_size, &out); + } + + gsl_odeiv2_evolve_free (e); + gsl_odeiv2_control_free (c); + gsl_odeiv2_step_free (s); + return out; } |