Changeset 16225 for branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Problem.cs

Timestamp:

10/10/18 11:22:43 (6 years ago)

Author:

gkronber

Message:

#2925: first working version of CVODES integration

File:

• branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Problem.cs (modified) (3 diffs)

branches/2925_AutoDiffForDynamicalModels/HeuristicLab.Problems.DynamicalSystemsModelling/3.3/Problem.cs

@@ -24 +24 @@
 using System.Diagnostics;
 using System.Linq;
+using System.Runtime.InteropServices;
 using HeuristicLab.Analysis;
 using HeuristicLab.Collections;
@@ -39 +40 @@
 
 namespace HeuristicLab.Problems.DynamicalSystemsModelling {
+
   public class Vector {
     public readonly static Vector Zero = new Vector(new double[0]);
@@ -126 +128 @@
   [StorableClass]
   public sealed class Problem : SingleObjectiveBasicProblem<MultiEncoding>, IRegressionProblem, IProblemInstanceConsumer<IRegressionProblemData>, IProblemInstanceExporter<IRegressionProblemData> {
+
+    // CVODES types
+    public enum MultistepMethod : int { CV_ADAMS = 1, CV_BDF = 2 };
+    public enum NonlinearSolverIteration : int { CV_NEWTON = 1, CV_FUNCTIONAL = 2 };
+
+
+    /* itask */
+    public const int CV_NORMAL = 1;
+    public const int CV_ONE_STEP = 2;
+
+    /* ism */
+    public const int CV_SIMULTANEOUS = 1;
+    public const int CV_STAGGERED = 2;
+    public const int CV_STAGGERED1 = 3;
+
+    /* DQtype */
+    public const int CV_CENTERED = 1;
+    public const int CV_FORWARD = 2;
+
+    /* interp */
+    public const int CV_HERMITE = 1;
+    public const int CV_POLYNOMIAL = 2;
+
+    /* 
+     * ----------------------------------------
+     * CVODES return flags
+     * ----------------------------------------
+     */
+
+    public const int CV_SUCCESS = 0;
+    public const int CV_TSTOP_RETURN = 1;
+    public const int CV_ROOT_RETURN = 2;
+
+    public const int CV_WARNING = 99;
+
+    public const int CV_TOO_MUCH_WORK = -1;
+    public const int CV_TOO_MUCH_ACC = -2;
+    public const int CV_ERR_FAILURE = -3;
+    public const int CV_CONV_FAILURE = -4;
+
+    public const int CV_LINIT_FAIL = -5;
+    public const int CV_LSETUP_FAIL = -6;
+    public const int CV_LSOLVE_FAIL = -7;
+    public const int CV_RHSFUNC_FAIL = -8;
+    public const int CV_FIRST_RHSFUNC_ERR = -9;
+    public const int CV_REPTD_RHSFUNC_ERR = -10;
+    public const int CV_UNREC_RHSFUNC_ERR = -11;
+    public const int CV_RTFUNC_FAIL = -12;
+    public const int CV_CONSTR_FAIL = -13;
+
+    public const int CV_MEM_FAIL = -20;
+    public const int CV_MEM_NULL = -21;
+    public const int CV_ILL_INPUT = -22;
+    public const int CV_NO_MALLOC = -23;
+    public const int CV_BAD_K = -24;
+    public const int CV_BAD_T = -25;
+    public const int CV_BAD_DKY = -26;
+    public const int CV_TOO_CLOSE = -27;
+
+    public const int CV_NO_QUAD = -30;
+    public const int CV_QRHSFUNC_FAIL = -31;
+    public const int CV_FIRST_QRHSFUNC_ERR = -32;
+    public const int CV_REPTD_QRHSFUNC_ERR = -33;
+    public const int CV_UNREC_QRHSFUNC_ERR = -34;
+
+    public const int CV_NO_SENS = -40;
+    public const int CV_SRHSFUNC_FAIL = -41;
+    public const int CV_FIRST_SRHSFUNC_ERR = -42;
+    public const int CV_REPTD_SRHSFUNC_ERR = -43;
+    public const int CV_UNREC_SRHSFUNC_ERR = -44;
+
+    public const int CV_BAD_IS = -45;
+
+    public const int CV_NO_QUADSENS = -50;
+    public const int CV_QSRHSFUNC_FAIL = -51;
+    public const int CV_FIRST_QSRHSFUNC_ERR = -52;
+    public const int CV_REPTD_QSRHSFUNC_ERR = -53;
+    public const int CV_UNREC_QSRHSFUNC_ERR = -54;
+
+    /* 
+     * ----------------------------------------
+     * CVODEA return flags
+     * ----------------------------------------
+     */
+
+    public const int CV_NO_ADJ = -101;
+    public const int CV_NO_FWD = -102;
+    public const int CV_NO_BCK = -103;
+    public const int CV_BAD_TB0 = -104;
+    public const int CV_REIFWD_FAIL = -105;
+    public const int CV_FWD_FAIL = -106;
+    public const int CV_GETY_BADT = -107;
+
+    [UnmanagedFunctionPointer(CallingConvention.Cdecl)]
+    public delegate int CVRhsFunc(
+        double t, // realtype
+        IntPtr y, // N_Vector
+        IntPtr ydot, // N_Vector
+        IntPtr user_data
+      );
+
+    [UnmanagedFunctionPointer(CallingConvention.Cdecl)]
+    public delegate int CVDlsJacFunc(
+      double t,
+      IntPtr y, // N_Vector
+      IntPtr fy, // N_Vector
+      IntPtr Jac, // SUNMatrix
+      IntPtr user_data,
+      IntPtr tmp1, // N_Vector
+      IntPtr tmp2, // N_Vector
+      IntPtr tmp3 // N_Vector
+      );
+
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "CVodeCreate", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    // returns a void* to the cvodes memory block if successful otherwise NULL
+    public static extern IntPtr CVodeCreate(MultistepMethod lmm, NonlinearSolverIteration iter);
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "CVodeInit", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+
+    public static extern int CVodeInit(
+      IntPtr cvode_mem, // pointer returned by CVodeCreate
+      CVRhsFunc f,
+      double t0, // realtype, the inital value of t
+      IntPtr y0 // N_Vector the initial value of y
+    );
+
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "CVodeSStolerances", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern int CVodeSStolerances(
+      IntPtr cvode_mem, 
+      double reltol, 
+      double abstol
+      );
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "CVDlsSetLinearSolver", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern int CVDlsSetLinearSolver(
+      IntPtr cvode_mem,
+      IntPtr linearSolver, // SUNLinearSolver
+      IntPtr j // SUNMatrix
+      );
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "CVDlsSetJacFn", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern int CVDlsSetJacFn(
+      IntPtr cvode_mem, 
+      CVDlsJacFunc jacFunc
+      );
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "CVode", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern int CVode(
+      IntPtr cvode_mem,
+      double tout, // next time at which a solution is desired
+      IntPtr yout, // N_Vector
+      ref double tret, // the time reached by the solver (output)
+      int itask // flag indicating the job of the solver for the next step. 
+      );
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "CVodeFree", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+
+    public static extern int CVodeFree(IntPtr cvode_mem);
+
+    #region matrix
+    [DllImport("sundials_cvodes.dll", EntryPoint = "SUNDenseMatrix", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern IntPtr SUNDenseMatrix(long m, long n);
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "SUNMatDestroy", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern void SUNMatDestroy(IntPtr A);
+    #endregion
+
+    #region linear solver
+    [DllImport("sundials_cvodes.dll", EntryPoint = "SUNDenseLinearSolver", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern IntPtr SUNDenseLinearSolver(
+      IntPtr y, // N_Vector
+      IntPtr A // SUNMatrix
+      );
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "SUNLinSolInitialize", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern int SUNLinSolInitialize(IntPtr linearSolver);
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "SUNLinSolSetup", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern int SUNLinSolSetup(
+      IntPtr linearSolver, 
+      IntPtr A // SUNMatrix
+      );
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "SUNLinSolSolve", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern int SUNLinSolSolve(
+      IntPtr linearSolver,
+      IntPtr A, // SUNMatrix
+      IntPtr x, // N_Vector
+      IntPtr b, // N_Vector
+      double tol
+      );
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "SUNLinSolFree", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern int SUNLinSolFree(IntPtr linearSolver);
+
+    #endregion
+
+    #region N_Vector
+    [DllImport("sundials_cvodes.dll", EntryPoint = "N_VNew_Serial", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern IntPtr N_VNew_Serial(long vec_length);
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "N_VDestroy_Serial", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern void N_VDestroy_Serial(IntPtr vec);
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "N_VPrint_Serial", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern void N_VPrint_Serial(IntPtr vec);
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "N_VConst_Serial", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern void N_VConst_Serial(double c, IntPtr vec);
+
+    
+    [DllImport("sundials_cvodes.dll", EntryPoint = "N_VL1Norm_Serial", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern double N_VL1Norm_Serial(IntPtr vec);
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "N_VMake_Serial", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public static extern IntPtr N_VMake_Serial(long vec_length, double[] v_data);
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "N_VScale", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    ///  Performs the operation z = c*x
+    public static extern void N_VScale(double s,
+      IntPtr x, // N_Vector
+      IntPtr z // N_Vector
+      );
+
+    [DllImport("sundials_cvodes.dll", EntryPoint = "N_VMake_Serial", ExactSpelling = true, CallingConvention = CallingConvention.Cdecl)]
+    public unsafe static extern double* N_VGetArrayPointer_Serial(IntPtr vec);
+    #endregion
+
+
+
+
 
     #region parameter names