[sundials-users] implementing custom linear solver in Fortran

[Mettupalayam Sivaselvan]

Hi,

 

I am trying to implement a custom linear solver in Fortran for use with IDA.

 

I am following the example in idaAnalytic_mels.c, and am looking to write Fortran equivalents of the 4 functions,

 

1. static SUNLinearSolver MatrixEmbeddedLS(void* ida_mem, SUNContext ctx);
2. static SUNLinearSolver_Type MatrixEmbeddedLSType(SUNLinearSolver S);
3. static int MatrixEmbeddedLSSolve(SUNLinearSolver S, SUNMatrix A, N_Vector x, N_Vector b, sunrealtype tol);
4. static SUNErrCode MatrixEmbeddedLSFree(SUNLinearSolver S);

 

I don’t anticipate any issues in implementing equivalents of the last three functions in Fortran.

 

However, for the first function, the problem I am running into is that since the SUNContext argument is passed by value, the SUNContext type would need to be defined in Fortran. However, I realized that fsundials_core_mod does not have a definition of SUNContext, since all the Fortran interface routine the use SUNContext pass it by value, and all occurrences of SUNContext in Fortran are type(c_ptr).

 

So I was wondering if the only way to implement a custom linear solver in Fortran was to implement the first function in C, and possibly the other 3 in Fortran.

 

Thank you.

 

Siva

---

To unsubscribe from the SUNDIALS-USERS list: write to: mailto:SUNDIALS-USERS-...@LISTSERV.LLNL.GOV

[Balos, Cody]

Hi Siva,

 

SUNContext is just a typedef to a pointer type so you should be able to just treat it as a type(c_ptr).  

 

Cody

[Mettupalayam Sivaselvan]

Hi Cody,

Thank you for your quick response. I'll give it a try.

One other minor question. FCreateLinSolNewEmpty returns a pointer to SUNLinearSolver, but MatrixEmbeddedLS is to return by value. Could you please advise how best to do this? Would there be a copy involved in MatrixEmbeddedLS?

Thank you again.

Siva

On Oct 1, 2024, at 5:16 PM, Balos, Cody <bal...@llnl.gov> wrote:



[Balos, Cody]

Hi Siva,

 

See this example https://github.com/LLNL/sundials/blob/main/examples/arkode/F2003_custom/fsunlinsol_fortran_mod.f90#L42. It shows how you may return a SUNLinearSolver from a Fortran function.   There is no copy of the underlying struct involved. SUNLinearSolver is also a typedef to a pointer to a struct.

[Mettupalayam Sivaselvan]

Hi Cody,

 

Thank you again. That example has been very helpful. I was a bit confused about when to use type(c_ptr) and SUNDIALS types, but it is clear now.

 

I am now able to run my code with IDA, and all the callbacks (the residual function and those associated with the custom linear solver) being called correctly.

 

I am getting the following error messages at the very first FIDASolve call. I am fairly certain my jacobian calculation and linear solve are correct. But I was wondering if based on the error message, you might have a suggestion on what the problem might be, or if I can compile IDA with options that would generate more descriptive output.

 

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 6.10352e-10, poor iterative algorithm performance. Nonlinear convergence failure rate is 7.000000e+00.

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 9.15527e-10, poor iterative algorithm performance. Nonlinear convergence failure rate is 4.000000e+00.

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 1.52588e-09, poor iterative algorithm performance. Nonlinear convergence failure rate is 2.666667e+00.

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 1.83105e-09, poor iterative algorithm performance. Nonlinear convergence failure rate is 2.250000e+00.

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 2.44141e-09, poor iterative algorithm performance. Nonlinear convergence failure rate is 1.800000e+00.

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 2.74658e-09, poor iterative algorithm performance. Nonlinear convergence failure rate is 1.666667e+00.

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 3.35693e-09, poor iterative algorithm performance. Nonlinear convergence failure rate is 1.428571e+00.

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 3.66211e-09, poor iterative algorithm performance. Nonlinear convergence failure rate is 1.375000e+00.

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 4.27246e-09, poor iterative algorithm performance. Nonlinear convergence failure rate is 1.222222e+00.

[WARNING][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida_ls.c:1644][idaLsPerf] Warning: at t = 4.57764e-09, poor iterative algorithm performance. Nonlinear convergence failure rate is 1.200000e+00.

[ERROR][rank 0][C:\Users\mvs\source\repos\sundials711\src\ida\ida.c:1341][IDASolve] At t = 4.44031e-07, , mxstep steps taken before reaching tout.

 

Thank you again.

 

Siva

[Balos, Cody]

Hi Siva,

 

Can you provide some details about the problem you are trying to solve?

[Mettupalayam Sivaselvan]

Hi Cody,

 

Please find attached a short description of the problem I am trying to solve.

 

I have computed a solution for these equations using ode15i in MATLAB, both without and with providing a Jacobian.

 

I wonder if as a way of troubleshooting, I should do this with IDA with a direct linear solver first – first having IDA compute the Jacobian, then providing the Jacobian, then eventually using the custom linear solver, so that I can identify where I might be going wrong.

 

Thank you.