KineticPreProcessor / KPP

The KPP kinetic preprocessor is a software tool that assists the computer simulation of chemical kinetic systems
GNU General Public License v3.0
21 stars 11 forks source link

Add #RETURNRATES switch to allow user to return equation rates from Fun() for diagnostic purposes #35

Closed yantosca closed 2 years ago

yantosca commented 2 years ago

We have added a new KPP switch #RETURNRATES on|off that will allow the user to return the equation rates direclty for diagnostic purposes. When #RETURNRATES on is added to a mechanism, a new optional argument Aout will be added to the Fun routine in the KPP_ROOT_Function module:

! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
! 
! Fun - time derivatives of variables - Aggregate form
!   Arguments :
!      V         - Concentrations of variable species (local)
!      F         - Concentrations of fixed species (local)
!      RCT       - Rate constants (local)
!      Vdot      - Time derivative of variable species concentrations
!      Aout      - Optional argument to return equation rate constants
! 
! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

SUBROUTINE Fun ( V, F, RCT, Vdot, Aout )

! V - Concentrations of variable species (local)
  REAL(kind=dp) :: V(NVAR)
! F - Concentrations of fixed species (local)
  REAL(kind=dp) :: F(NFIX)
! RCT - Rate constants (local)
  REAL(kind=dp) :: RCT(NREACT)
! Vdot - Time derivative of variable species concentrations
  REAL(kind=dp) :: Vdot(NVAR)
! Aout - Optional argument to return equation rate constants
  REAL(kind=dp), OPTIONAL :: Aout(NREACT)

! Computation of equation rates
  A(1) = RCT(1)*F(2)
  A(2) = RCT(2)*V(2)*F(2)
  .. etc. ..

!### Use Aout to return equation rates
  IF ( PRESENT( Aout ) ) Aout = A

! Aggregate function
  Vdot(1) = A(5)-A(6)-A(7)
  Vdot(2) = 2*A(1)-A(2)+A(3)-A(4)+A(6)-A(9)+A(10)
  .. etc ..

If this switch is omitted, then no modification is made to the Fun() routine:

SUBROUTINE Fun ( V, F, RCT, Vdot )

! V - Concentrations of variable species (local)
  REAL(kind=dp) :: V(NVAR)
! F - Concentrations of fixed species (local)
  REAL(kind=dp) :: F(NFIX)
! RCT - Rate constants (local)
  REAL(kind=dp) :: RCT(NREACT)
! Vdot - Time derivative of variable species concentrations
  REAL(kind=dp) :: Vdot(NVAR)
... etc...

This switch is needed to provide support for the use case of running GEOS-Chem within the NASA GEOS Earth System Model.

msl3v commented 2 years ago

What are the consequences of just leaving the OPTIONAL statement w/o a switch? Is there a risk that it won't work with some compilers, etc?

RolfSander commented 2 years ago

@msl3v has a good point. I don't think that adding an OPTIONAL parameter would do any harm. I think that Aout can be added without the need for #RETURNRATES.

RolfSander commented 2 years ago

In the future, we may have variables that are both POINTERs and OPTIONAL. Maybe we could define ATTR_F90_OPT=4 (instead of 3) and then apply a similar logic as for ICNTRL(15) and Update_*. Then we can create an OPTIONAL POINTER with:

DefineVariable( name, VELM, bt, 0, 0, cmt, ATTR_F90_PTR+ATTR_F90_OPT )

yantosca commented 2 years ago

I can add that in too.

yantosca commented 2 years ago

What are the consequences of just leaving the OPTIONAL statement w/o a switch? Is there a risk that it won't work with some compilers, etc

No, it was just that I was thinking not to modify the original code. I could of course add Aout as OPTIONAL without the #RETURNRATES switch.

yantosca commented 2 years ago

I'll go ahead and modify it to add Aout as optional without #RETURNRATES and also add a macro for defining an optional pointer.

Then after that I'll start working on adding the ReadTheDocs documentation. Then I think we can release KPP 2.5.0, which will have added all of the "hacks" that we made for GEOS-Chem as optional switches.

KPP 3.0 will be the introduction of the autoreduction code, which is still undergoing work.

Thank you for the feedback!

RolfSander commented 2 years ago

it was just that I was thinking not to modify the original code.

Yes, that's very important, indeed. For backward compatibility, we need to make sure that those subroutines return the same results as in previous versions. However, when adding an optional parameter doesn't change any results (as long as you don't supply that optional parameter), I think we're on the safe side.

RolfSander commented 2 years ago

Then I think we can release KPP 2.5.0, which will have added all of the "hacks" that we made for GEOS-Chem as optional switches.

Yes, I think we have all the new features then. Before releasing 2.5.0, I'd like to do some code cleanup. Nothing major, but a couple of small things:

yantosca commented 2 years ago

Thanks @RolfSander for the extra points. I agree we should clean those up. I will start a new feature request so we don't lose sight of these. I can try to hammer some of these out starting now.

RolfSander commented 2 years ago

Thanks, @yantosca!

For the gcc warnings I certainly do need help. Some of the warnings are about FAMILIES and some about ATTRF90* comparisons between pointer and integer.

For the other items in the list, I'm happy to take care of them (unless you're faster :-)

yantosca commented 2 years ago

I can take care of these, @RolfSander.

Also, in dev, I should have fixed the ATTRF90 problem. I had declared var->attr as an int instead of as an int. Removing the * fixed the issue.

yantosca commented 2 years ago

also see #36