[Commits] [svn:einsteintoolkit] incoming/EOS_Omni/src/ (Rev. 17)

cott at tapir.caltech.edu cott at tapir.caltech.edu
Fri Nov 12 22:13:28 CST 2010


User: cott
Date: 2010/11/12 10:13 PM

Modified:
 /EOS_Omni/src/
  EOS_Omni_SingleVarCalls.F90
 /EOS_Omni/src/nuc_eos/
  nuc_eos.F90

Log:
 * small optimization: introduce call to nuc_eos that just returns the pressure (and eps)

File Changes:

Directory: /EOS_Omni/src/nuc_eos/
=================================

File [modified]: nuc_eos.F90
Delta lines: +101 -0
===================================================================
--- EOS_Omni/src/nuc_eos/nuc_eos.F90	2010-11-12 00:35:41 UTC (rev 16)
+++ EOS_Omni/src/nuc_eos/nuc_eos.F90	2010-11-13 04:13:28 UTC (rev 17)
@@ -264,6 +264,87 @@
 
 end subroutine nuc_eos_short
 
+
+subroutine nuc_eos_press_eps(xrho,xtemp,xye,xenr,xprs,&
+     keytemp,keyerr,rfeps)
+
+  use eosmodule
+  implicit none
+
+  real*8, intent(in)    :: xrho,xye
+  real*8, intent(inout) :: xtemp,xenr
+  real*8, intent(in)    :: rfeps
+  real*8, intent(out)   :: xprs
+  integer, intent(in)   :: keytemp
+  integer, intent(out)  :: keyerr
+
+  ! local variables
+  real*8 :: xcs2,xdpderho,xdpdrhoe
+  real*8 :: lr,lt,y,xx,xeps,leps,xs
+  real*8 :: d1,d2,d3,ff(8)
+  integer :: keyerrt = 0
+
+  if(xrho.gt.eos_rhomax) then
+     stop "nuc_eos: rho > rhomax"
+  endif
+
+  if(xrho.lt.eos_rhomin*1.2d0) then
+     call nuc_low_eos(xrho,xenr,xprs,xcs2,xdpderho,xdpdrhoe,keytemp)
+     return
+  endif
+
+  if(xye.gt.eos_yemax) then
+     stop "nuc_eos: ye > yemax"
+  endif
+
+  if(xye.lt.eos_yemin) then
+     stop "nuc_eos: ye < yemin"
+  endif
+
+  if(keytemp.eq.1) then
+     if(xtemp.gt.eos_tempmax) then
+        stop "nuc_eos: temp > tempmax"
+     endif
+     
+     if(xtemp.lt.eos_tempmin) then
+        call nuc_low_eos(xrho,xenr,xprs,xcs2,xdpderho,xdpdrhoe,keytemp)
+        return
+     endif
+  endif
+
+  lr = log10(xrho)
+  lt = log10(xtemp)
+  y = xye
+  xeps = xenr + energy_shift
+  leps = log10(max(xeps,1.0d0))
+
+  keyerr = 0
+
+  if(keytemp.eq.0) then
+     !need to find temperature based on xeps
+     call findtemp(lr,lt,y,leps,keyerrt,rfeps)
+     if(keyerrt.ne.0) then
+        keyerr = keyerrt
+        return
+     endif
+     xtemp = 10.0d0**lt
+
+  elseif(keytemp.eq.2) then
+     stop "eos_nuc_press does not support keytemp.eq.2"
+  endif
+
+  ! have temperature, proceed:
+  call findall_press_eps(lr,lt,y,ff)
+  xprs = 10.0d0**ff(1)
+
+  xeps = 10.0d0**ff(2) - energy_shift
+  if((keytemp.eq.1).or.(keytemp.eq.2)) then
+     xenr = xeps
+  endif
+
+end subroutine nuc_eos_press_eps
+
+
 subroutine findthis(lr,lt,y,value,array,d1,d2,d3)
 
   use eosmodule
@@ -321,3 +402,23 @@
   ff(:) = ffx(:,1)
 
 end subroutine findall_short
+
+
+subroutine findall_press_eps(lr,lt,y,ff)
+
+  use eosmodule
+  implicit none
+
+  real*8 ffx(8,1)
+  real*8 ff(8)
+  real*8 lr,lt,y
+  integer i
+  integer :: nvarsx = 2
+
+
+! Ewald's interpolator           
+  call intp3d_many(lr,lt,y,ffx,1,alltables(:,:,:,1:2), &
+       nrho,ntemp,nye,nvarsx,logrho,logtemp,ye)
+  ff(:) = ffx(:,1)
+
+end subroutine findall_press_eps

Directory: /EOS_Omni/src/
=========================

File [modified]: EOS_Omni_SingleVarCalls.F90
Delta lines: +1 -2
===================================================================
--- EOS_Omni/src/EOS_Omni_SingleVarCalls.F90	2010-11-12 00:35:41 UTC (rev 16)
+++ EOS_Omni/src/EOS_Omni_SingleVarCalls.F90	2010-11-13 04:13:28 UTC (rev 17)
@@ -96,8 +96,7 @@
            xtemp = temp(i)
            xye = ye(i)
            xenr = eps(i) * inv_eps_gf
-           call nuc_eos_short(xrho,xtemp,xye,xenr,xprs,& 
-                xent,xcs2,xdedt,xdpderho,xdpdrhoe,xmunu,&
+           call nuc_eos_press_eps(xrho,xtemp,xye,xenr,xprs,& 
                 keytemp,keyerr(i),rf_precision)
            
            if(keyerr(i).ne.0) then



More information about the Commits mailing list