Using the sundials solver package from D - my first experience with importC

[Robert Kovacs]

1. Wrap the necessary headers to `sundials.c`:

```d
#define _Float16 float // don't care why we don't find _Float16 
on Apple Silicon

#include "cvode/cvode.h"
#include "nvector/nvector_serial.h"
#include "sunmatrix/sunmatrix_dense.h"
#include "sunlinsol/sunlinsol_dense.h"
```

2. Port the cvRoberts_dns.c example to D:

```d
import sundials;
import std.stdio;

extern (C) static int f(sunrealtype t,
                         N_Vector y,
                         N_Vector ydot,
                         void* user_data)
{
   sunrealtype y1, y2, y3, yd1, yd3;

   y1 = NV_DATA_S(y)[0];
   y2 = NV_DATA_S(y)[1];
   y3 = NV_DATA_S(y)[2];

   yd1 = NV_DATA_S(ydot)[0] = -0.04 * y1 + 1.0e4 * y2 * y3;
   yd3 = NV_DATA_S(ydot)[2] = 3.0e7 * y2 * y2;
   NV_DATA_S(ydot)[1]       = -yd1 - yd3;

   return (0);
}

extern (C) static int jac(sunrealtype t,
                           N_Vector y,
                           N_Vector fy,
                           SUNMatrix J,
                           void* user_data,
                           N_Vector tmp1,
                           N_Vector tmp2,
                           N_Vector tmp3)
{
     double y2 = NV_DATA_S(y)[1];
     double y3 = NV_DATA_S(y)[2];

     auto J_content = cast(SUNMatrixContent_Dense) J.content;

     J_content.cols[0][0] = -0.04;
     J_content.cols[1][0] = 1.0e4 * y3;
     J_content.cols[2][0] = 1.0e4 * y2;

     J_content.cols[0][1] = 0.04;
     J_content.cols[1][1] = -1.0e4 * y3 - 6.0e7 * y2;
     J_content.cols[2][1] = -1.0e4 * y2;

     J_content.cols[0][2] = 0.0;
     J_content.cols[1][2] = 6.0e7 * y2;
     J_content.cols[2][2] = 0.0;

     return (0);
}

extern(C) static int g(sunrealtype t,
                        N_Vector y,
                        double* gout,
                        void* user_data)
{
     double y1 = NV_DATA_S(y)[0];
     double y3 = NV_DATA_S(y)[2];

     gout[0] = y1 - 0.0001;
     gout[1] = y3 - 0.01;

     return 0;
}

void main()
{
     SUNContext sunctx;
     int retval = SUNContext_Create(SUN_COMM_NULL, &sunctx);
     writeln("SUNContext created with return value ", retval);

     N_Vector y = N_VNew_Serial(3, sunctx);
     NV_DATA_S(y)[0] = 1.0;
     NV_DATA_S(y)[1] = 0.0;
     NV_DATA_S(y)[2] = 0.0;
     writeln("State vector created");

     void* cvode_mem = CVodeCreate(CV_BDF, sunctx);
     retval = CVodeInit(cvode_mem, &f, 0.0, y);
     writeln("cvode_mem initialized, return value is ", retval);


     SUNMatrix A = SUNDenseMatrix(3, 3, sunctx);

     SUNLinearSolver LS = SUNLinSol_Dense(y, A, sunctx);
     retval = CVodeSetLinearSolver(cvode_mem, LS, A);
     writeln("Linear solver set, return value is ", retval);

     retval = CVodeSetJacFn(cvode_mem, &jac);
     writeln("Jacobian function set, return value is ", retval);

     N_Vector abstol = N_VNew_Serial(3, sunctx);
     NV_DATA_S(abstol)[0] = 1.0e-8;
     NV_DATA_S(abstol)[1] = 1.0e-14;
     NV_DATA_S(abstol)[2] = 1.0e-6;
     writeln("Absolute tolerance vector created");
     retval = CVodeSVtolerances(cvode_mem, 1.0e-4, abstol);

     retval = CVodeRootInit(cvode_mem, 2, &g);
     writeln("Root finding initialized, return value is ", retval);

     double t = 0.0;
     double tout = 0.4;

     int[2] roots_found;
     int n_out = 12;
     int i_out = 0;

     double t_mult = 10.0;

     while(true)
     {
         retval = CVode(cvode_mem, tout, y, &t, CV_NORMAL);
         writeln(t, "\t", NV_DATA_S(y)[0], "\t", NV_DATA_S(y)[1], 
"\t", NV_DATA_S(y)[2]);

         if (retval == CV_ROOT_RETURN)
         {
             retval = CVodeGetRootInfo(cvode_mem, roots_found.ptr);
             writeln("Root found at t = ", t, " with g1 = ", 
roots_found[0], " and g2 = ", roots_found[1]);
         }

         if (retval == CV_SUCCESS)
         {
             i_out++;
             tout *= t_mult;
         }

         if (i_out == n_out) break;
     }
}

```

(Pls pardon my French, this has been done in cca. 15 minutes)

3. Build:

```sh
ldc2 sundials_example.d sundials.c -L-L/usr/local/lib 
-L-lsundials_core -L-lsundials_nvecserial -L-lsundials_cvode
```

4. Enjoy!

This. Is. Amazing.

Thank you for all involved in importC development!