search

AMReX-Astro / Microphysics

common astrophysical microphysics routines with interfaces for the different AMReX codes
https://amrex-astro.github.io/Microphysics
Other
35 stars 35 forks source link

centered numerical Jac has wrong scaling? #179

Closed zingale closed 5 years ago

zingale commented 5 years ago

In numerical_jacobian.F90, for the centered_diff part, it looks like we scale the composition terms with aion, but we don't do this in the forward difference section. My informal tests show the centered diff version crashy, so this is likely wrong.

dwillcox commented 5 years ago

Yeah, it looks like we should be doing this in the forward difference section as well.

dwillcox commented 5 years ago

Recapping the Slack discussion --

For comparison, the analytic jacobian routine returns derivatives in terms of Y, which are converted to X.

Species wrt species: Jij = d(dYi/dt)/dYj Species wrt temperature: Jij = d(dYi/dt)/dT

Then the conversion is done by the calling routine as:

    do n = 1, nspec_evolve                                                                                                                                                                                         
       state % jac(n,:) = state % jac(n,:) * aion(n)                                                                                                                                                               
       state % jac(:,n) = state % jac(:,n) * aion_inv(n)                                                                                                                                                           
    enddo

By looking at the calls to numerical_jac in VBDF and BS, those integrators expect numerical_jac to return the following, corresponding to the integration variables. There is no conversion from Y to X by the calling routines.

Species wrt species: Jij = d(dXi/dt)/dXj Species wrt temperature: Jij = d(dXi/dt)/dT

The centered difference numerical jacobian appears to be returning the following:

Species wrt species: Jij = d(dXi/dt)/dXj Species wrt temperature: Jij = d(dYi/dt)/dT

The forward difference numerical jacobian appears to be returning the following:

Species wrt species: Jij = d(dYi/dt)/dXj Species wrt temperature: Jij = d(dYi/dt)/dT

So both versions of the numerical jacobian are incorrect. Both return the species derivatives wrt temperature in terms of Y, and the foward difference also returns species derivatives wrt species in terms of mixed Y and X.

dwillcox commented 5 years ago

Closed in #180