Skip to content

Conversation

@maqsoodrajput
Copy link

I have designed a new solver that is faster than the default Rosenbrock solvers. This new solver is inspired by the IMEX solver. I developed it using the standard BE and Sdirk2a because an explicit method would not work for stiff problems. This new solver gave me a 5X speed-up when I used it with the ICON-ART model. The speed-up outweighed the minor loss in accuracy.

There are two files added:

  • implicit12.f90
  • implicit12.def

I have designed a new solver that is faster than the default Rosenbrock solvers.  
This new solver is inspired by the IMEX solver. I developed it using the standard BE and Sdirk2a because an explicit method would not work for stiff problems. 
This new solver gave me a 5X speed-up when I used it with the ICON-ART model. The speed-up outweighed the minor loss in accuracy.    

There are two files added:
- implicit12.f90
- implicit12.def
@yantosca yantosca self-assigned this Jun 10, 2024
@yantosca
Copy link
Contributor

Thanks @maqsoodrajput for this. I am at the International GEOS-Chem Conference this week so I may be slow to get to this PR. But this sounds like a great new feature for KPP. I can also try to set up the C-I tests for this integrator.

Also tagging @RolfSander @obin1 @msl3v @jimmielin

@yantosca yantosca requested review from RolfSander and yantosca June 10, 2024 15:58
@yantosca yantosca added the integrators Related to numerical integrators label Jun 10, 2024
@yantosca yantosca changed the base branch from main to dev June 10, 2024 15:59
@msl3v
Copy link
Contributor

msl3v commented Jun 10, 2024

@maqsoodrajput excellent!

@RolfSander
Copy link
Contributor

Thanks @maqsoodrajput for providing the new KPP integrator! I will run
it with my MECCA chemistry and compare it to Rosenbrock. Hopefully, I'll
find the time next week for this...

@yantosca
Copy link
Contributor

@maqsoodrajput @RolfSander: I have been running into an issue where the updated time isn't getting passed back via RSTATUS. I set up a test case with the small_strato mechanism.

   0.0%. T=0.432E+05  O1D= 0.9906E+02; O= 0.6624E+09; O3= 0.5326E+12; NO= 0.8725E+09; NO2= 0.2240E+09; O2= 0.1697E+17; 
 ### TIN, TOUT, in main    43200.000000000000        44100.000000000000 
 ### ierr:, in sdirk           0 
 ### tinitial, tfinal, in sdirk before calling sdirk_integrator:    43200.000000000000        44100.000000000000  
 ### RSTATUS in Sdirk:   0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000     
 ### RSTATE in Main    0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000        0.0000000000000000     
 ### T   0.0000000000000000 

I'm still digging but this means that the simulation is basically stuck in an infinite loop on the same timestep. Feedback greatly appreciated! My updates are in the feature/imex-solver branch.

@yantosca
Copy link
Contributor

I'm still digging but this means that the simulation is basically stuck in an infinite loop on the same timestep. Feedback greatly appreciated! My updates are in the feature/imex-solver branch.

I've added the following code at the end of SDirk_Integrator, which solves the issue in the previous comment:

!~~~>  End of simplified Newton iterations

      ! Successful return
      Ierr  = 1

      ! Return updated time
      RSTATUS(Ntexit) = Tfinal
      RSTATUS(Nhexit) = H

      END SUBROUTINE SDIRK_Integrator

Let me know if you all think this isn't the right thing to do.

@yantosca yantosca self-requested a review June 13, 2024 16:22
Copy link
Contributor

@yantosca yantosca left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I now approve as I was able to get the C-I test for implicit 12 to work properly.

@RolfSander
Copy link
Contributor

@maqsoodrajput: Would you like to write a few lines about your new
solver which we can add to our KPP documentation? Probably somewhere
here:
https://kpp.readthedocs.io/en/stable/tech_info/07_numerical_methods.html#sdirk

@yantosca yantosca added the feature New feature or request label Jan 27, 2025
@yantosca
Copy link
Contributor

@maqsoodrajput: I may try your new integrator within the GEOS-Chem mechanism. Would you have a paper reference that I can cite (and that we can add to the KPP documentation)?

@yantosca yantosca added this to the 3.2.0 milestone Jan 27, 2025
@maqsoodrajput
Copy link
Author

Hi @yantosca, apologies for not responding.

I could not secure funding for this project, therefore I could not finish the documentation/paper.
I can try to finish the documentation for this new solver for KPP in a few weeks.

I would be happy to participate in a paper. I have tested the solver in the box model and it also works well for ICON-ART. It would be great to run it with GEOS-Chem mechanism.

@maqsoodrajput
Copy link
Author

@maqsoodrajput: Would you like to write a few lines about your new solver which we can add to our KPP documentation? Probably somewhere here: https://kpp.readthedocs.io/en/stable/tech_info/07_numerical_methods.html#sdirk

Hi, apologies for super late response. The project did get extension so I had to stop on this one. I can find some time now and write documentation.

@RolfSander
Copy link
Contributor

@maqsoodrajput mentioned that he can find some time now to write the
documentation. Should we wait with the merge until the documentation is
available?

@yantosca
Copy link
Contributor

yantosca commented Feb 5, 2025

@RolfSander yes we can wait.

@yantosca yantosca mentioned this pull request Feb 27, 2025
@yantosca yantosca modified the milestones: 3.2.0, 3.3.0 Feb 27, 2025
@yantosca yantosca modified the milestones: 3.3.0, 3.4.0 Jul 17, 2025
@yantosca yantosca removed this from the 3.4.0 milestone Oct 9, 2025
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

feature New feature or request integrators Related to numerical integrators

Projects

None yet

Development

Successfully merging this pull request may close these issues.

4 participants