Skip to content

Commit

Permalink
change rate inteperpolation to logRate instead of Rate, for log-log i…
Browse files Browse the repository at this point in the history
…nterpolation
  • Loading branch information
Debraheem committed Sep 6, 2024
1 parent 474f301 commit a7d6526
Show file tree
Hide file tree
Showing 2 changed files with 37 additions and 26 deletions.
30 changes: 15 additions & 15 deletions net/test/test_output
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@

eps_nuc 0.839654151E+13

d_lneps_dlnT 9.328757835
d_lneps_dlnT 9.328757799
d_lneps_dlnRho 1.009292766


Expand All @@ -34,7 +34,7 @@

eps_nuc 0.302124143E+11

d_lneps_dlnT 11.356884345
d_lneps_dlnT 11.356884263
d_lneps_dlnRho 1.015667220


Expand All @@ -59,7 +59,7 @@

eps_nuc 0.835361769E+13

d_lneps_dlnT 9.326560302
d_lneps_dlnT 9.326560266
d_lneps_dlnRho 1.009292884


Expand All @@ -84,7 +84,7 @@

eps_nuc 0.536017372E+09

d_lneps_dlnT 11.559205087
d_lneps_dlnT 11.559204995
d_lneps_dlnRho 1.013842798


Expand All @@ -109,7 +109,7 @@

eps_nuc 0.150660462E+25

d_lneps_dlnT -0.253424018
d_lneps_dlnT -0.253423995
d_lneps_dlnRho 1.008418067


Expand Down Expand Up @@ -151,7 +151,7 @@

eps_nuc 0.150660462E+25

d_lneps_dlnT -0.253424018
d_lneps_dlnT -0.253423995
d_lneps_dlnRho 1.008418067


Expand Down Expand Up @@ -194,22 +194,22 @@
test_one_zone_burn_small_net
number of species 21
large final abundances 1.0000000000000000D-02
c12 1 9.8714342254068388D-01
he4 2 1.2802186836379657D-02
c12 1 9.8714342283429213D-01
he4 2 1.2802186549072116D-02

xsum 1.0000000000006330D+00
xsum 1.0000000000005884D+00


test_one_zone_burn_const_P
number of species 21
large final abundances 1.0000000000000000D-02
o16 1 7.1798810044269901D-01
si28 2 1.6142410416068431D-01
s32 3 6.0439378324441878D-02
ar36 4 2.9692711889885499D-02
ca40 5 1.6560528697491497D-02
o16 1 7.1798810044342576D-01
si28 2 1.6142410415591074D-01
s32 3 6.0439378327695317D-02
ar36 4 2.9692711893868861D-02
ca40 5 1.6560528702226556D-02

xsum 1.0000000000000000D+00
xsum 1.0000000000000002D+00



33 changes: 22 additions & 11 deletions rates/private/rates_support.f90
Original file line number Diff line number Diff line change
Expand Up @@ -167,17 +167,14 @@ subroutine get_rates_from_table(r1, r2)

do i = r1,r2

rate_raw(i) = &
(rattab_f(1,k,i) + dt*(rattab_f(2,k,i) + &
dt*(rattab_f(3,k,i) + dt*rattab_f(4,k,i))) &
) * dtab(i)
rate_raw(i) = exp10(rattab_f(1,k,i) + dt*(rattab_f(2,k,i) + &
dt*(rattab_f(3,k,i) + dt*rattab_f(4,k,i))) ) * dtab(i)

rate_raw_dRho(i) = rate_raw(i) * ddtab(i) / dtab(i)
! Derivative with respect to density (dRho)
rate_raw_dRho(i) = rate_raw(i) * (ddtab(i) / dtab(i))

rate_raw_dT(i) = &
(rattab_f(2,k,i) + 2*dt*(rattab_f(3,k,i) + &
1.5d0*dt*rattab_f(4,k,i)) &
) * dtab(i) / (btemp * ln10)
! Derivative with respect to temperature (dT)
rate_raw_dT(i) = rate_raw(i) * (rattab_f(2,k,i) + 2*dt*rattab_f(3,k,i) + 3*dt**2*rattab_f(4,k,i)) / (btemp)

end do

Expand All @@ -201,7 +198,7 @@ subroutine do_make_rate_tables( &
integer, intent(out) :: ierr

integer :: i, j, operr, num_to_add_to_cache,thread_num
real(dp) :: logT, btemp
real(dp) :: logT, btemp, rate_logR
real(dp), pointer :: work1(:)=>null(), f1(:)=>null(), rattab_f(:,:,:)=>null()
integer, pointer :: reaction_id(:) =>null()
real(dp), allocatable, target :: work(:,:)
Expand All @@ -211,7 +208,7 @@ subroutine do_make_rate_tables( &
include 'formats'

ierr = 0

rate_logR = 0
rattab_f(1:4,1:nrattab,1:num_reactions) => &
rattab_f1(1:4*nrattab*num_reactions)

Expand Down Expand Up @@ -282,6 +279,20 @@ subroutine do_make_rate_tables( &
write(*, '(a,i4,2x,a)') 'missing raw rate for ', &
j, trim(reaction_Name(reaction_id(j)))
a_okay = .false.
else
! **New Addition: Store the logarithm of the rate**
if (rattab(j, i) > 0.0_dp) then ! Only take log of positive values
rate_logR = log10(rattab(j, i))
else if (rattab(j, i) == 0.0_dp) then
rate_logR = -999
else
! warning for rates set to -1
! like rni56ec_to_co56,rco56ec_to_fe56
!write(*, '(a,i4,2x,a)') 'Warning: non-positive value for rate in ', &
! j, trim(reaction_Name(reaction_id(j)))
rate_logR = rattab(j, i) ! don't convert to log space.
end if
rattab(j, i) = rate_logR ! Store value in rattab
end if
end do
if (.not. a_okay) all_okay = .false.
Expand Down

0 comments on commit a7d6526

Please sign in to comment.