[svn:einsteintoolkit] incoming/EOS_Omni/ (Rev. 12) - Commits

3 Nov 2010

User: cott
Date: 2010/11/03 11:49 AM
Modified:
 /EOS_Omni/
  configuration.ccl
 /EOS_Omni/src/
  EOS_Omni_Handles.c, EOS_Omni_Module.F90, EOS_Omni_SingleVarCalls.F90
 /EOS_Omni/src/nuc_eos/
  readtable.F90
Log:
 * fix a number of small issues
File Changes:
Directory: /EOS_Omni/
=====================
File [modified]: configuration.ccl
Delta lines: +5 -0
===================================================================

--- EOS_Omni/configuration.ccl	2010-10-24 14:06:52 UTC (rev 11)
+++ EOS_Omni/configuration.ccl	2010-11-03 16:49:03 UTC (rev 12)
@@ -1,2 +1,7 @@
 REQUIRES HDF5
+PROVIDES EOS_Omni
+{
+  SCRIPT
+  LANG
+}
Directory: /EOS_Omni/src/nuc_eos/
=================================
File [modified]: readtable.F90
Delta lines: +5 -5
===================================================================
--- EOS_Omni/src/nuc_eos/readtable.F90	2010-10-24 14:06:52 UTC (rev 11)
+++ EOS_Omni/src/nuc_eos/readtable.F90	2010-11-03 16:49:03 UTC (rev 12)
@@ -21,13 +21,13 @@
   real*8 buffer1,buffer2,buffer3,buffer4
   accerr=0
-  write(*,*) "Reading Nuclear EOS Table"
+!  write(*,*) "Reading Nuclear EOS Table"
call h5open_f(error)
call h5fopen_f (trim(adjustl(eos_filename)), H5F_ACC_RDONLY_F, file_id, error)
-  write(6,*) trim(adjustl(eos_filename))
+!  write(6,*) trim(adjustl(eos_filename))
! read scalars
   dims1(1)=1
@@ -57,8 +57,8 @@
      stop "Could not read EOS table file"
   endif
-  write(message,"(a25,1P3i5)") "We have nrho ntemp nye: ", nrho,ntemp,nye
-  write(*,*) message
+!  write(message,"(a25,1P3i5)") "We have nrho ntemp nye: ", nrho,ntemp,nye
+!  write(*,*) message
allocate(alltables(nrho,ntemp,nye,nvars))
@@ -226,7 +226,7 @@
   eos_tempmin = 10.0d0**logtemp(1)
   eos_tempmax = 10.0d0**logtemp(ntemp)
-  write(6,*) "Done reading eos tables"
+!  write(6,*) "Done reading eos tables"
end subroutine nuc_eos_readtable
Directory: /EOS_Omni/src/
=========================
File [modified]: EOS_Omni_Handles.c
Delta lines: +2 -0
===================================================================
--- EOS_Omni/src/EOS_Omni_Handles.c	2010-10-24 14:06:52 UTC (rev 11)
+++ EOS_Omni/src/EOS_Omni_Handles.c	2010-11-03 16:49:03 UTC (rev 12)
@@ -10,6 +10,8 @@
         return 2;
     if (CCTK_EQUALS(name, "Hybrid"))
         return 3;
+    if (CCTK_EQUALS(name, "nuc_eos"))
+        return 4;
     return 0;
 }
File [modified]: EOS_Omni_Module.F90
Delta lines: +3 -4
===================================================================
--- EOS_Omni/src/EOS_Omni_Module.F90	2010-10-24 14:06:52 UTC (rev 11)
+++ EOS_Omni/src/EOS_Omni_Module.F90	2010-11-03 16:49:03 UTC (rev 12)
@@ -3,23 +3,22 @@
   implicit none
real*8,parameter :: rho_gf = 1.61620075314614d-18
-!  real*8,parameter :: rho_gf = 1.61930347d-18
   real*8,parameter :: press_gf = 1.7982953469278d-39
-!  real*8,parameter :: press_gf = 1.80171810d-39
   real*8,parameter :: eps_gf = 1.11265006d-21
   real*8,parameter :: time_gf = 2.03001708d+05
   real*8,parameter :: mass_gf = 5.02765209d-34
   real*8,parameter :: length_gf = 6.77140812d-06
real*8,parameter :: inv_rho_gf = 6.18735016707159d17
-!  real*8,parameter :: inv_rho_gf = 6.17549470205236d17
   real*8,parameter :: inv_press_gf = 5.56082181777535d38
-!  real*8,parameter :: inv_press_gf = 5.55025783445257d38
   real*8,parameter :: inv_eps_gf = 8.98755175549085d20
   real*8,parameter :: inv_time_gf = 4.92606692747629d-6
   real*8,parameter :: inv_mass_gf = 1.98899999860571d33
   real*8,parameter :: inv_length_gf = 147679.77092481d0
+  real*8,parameter :: clite = 2.99792458d10
+  real*8,parameter :: cliteinv2 = 1.11265005605362d-21
+
   real*8 :: poly_k_cgs = 0.0d0
   real*8 :: gl_k_cgs   = 0.0d0
File [modified]: EOS_Omni_SingleVarCalls.F90
Delta lines: +130 -5
===================================================================
--- EOS_Omni/src/EOS_Omni_SingleVarCalls.F90	2010-10-24 14:06:52 UTC (rev 11)
+++ EOS_Omni/src/EOS_Omni_SingleVarCalls.F90	2010-11-03 16:49:03 UTC (rev 12)
@@ -30,6 +30,10 @@
   real*8           :: hybrid_local_gamma, hybrid_local_k_cgs, &
                       hybrid_p_poly, hybrid_p_th
   real*8,parameter :: zero = 0.0d0
+  ! temporary vars for nuc_eos
+  real*8           :: xrho,xye,xtemp,xenr,xent
+  real*8           :: xprs,xmunu,xcs2
+  real*8           :: xdedt,xdpderho,xdpdrhoe
anyerr    = 0
   keyerr(:) = 0
@@ -84,13 +88,32 @@
            press(i) = hybrid_p_poly + hybrid_p_th
         enddo
-!!!     case (4)
-!!!
-!!!        do i=1,npoints
-!!!
-!!!        enddo
+     case (4)
+        do i=1,npoints
+           xrho = rho(i) * inv_rho_gf
+           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,&
+                keytemp,keyerr(i),rf_precision)
+           
+           if(keyerr(i).ne.0) then
+              anyerr = 1
+           endif
+
+           if(keytemp.eq.1) then
+              eps(i) = xenr * eps_gf
+           else
+              temp(i) = xtemp
+           endif
+
+           press(i) = xprs * press_gf
+
+        enddo
+
      case DEFAULT
         write(warnstring,*) "eoskey ",eoskey," not implemented!"
         call CCTK_WARN(0,warnstring)
@@ -121,6 +144,10 @@
   real*8           :: hybrid_local_gamma, hybrid_local_k_cgs, &
                       hybrid_dp_poly, hybrid_dp_th1, hybrid_dp_th2
   real*8,parameter :: zero = 0.0d0
+  ! temporary vars for nuc_eos
+  real*8           :: xrho,xye,xtemp,xenr,xent
+  real*8           :: xprs,xmunu,xcs2
+  real*8           :: xdedt,xdpderho,xdpdrhoe
anyerr    = 0
   keyerr(:) = 0
@@ -181,6 +208,28 @@
enddo
+     case (4)
+        do i=1,npoints
+           xrho = rho(i) * inv_rho_gf
+           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,&
+                keytemp,keyerr(i),rf_precision)
+           
+           if(keyerr(i).ne.0) then
+              anyerr = 1
+           endif
+
+           if(keytemp.eq.1) then
+              eps(i) = xenr * eps_gf
+           else
+              temp(i) = xtemp
+           endif
+
+           dpdrhoe(i) = xdpdrhoe * press_gf * inv_rho_gf
+        enddo
      case DEFAULT
         write(warnstring,*) "eoskey ",eoskey," not implemented!"
         call CCTK_WARN(0,warnstring)
@@ -206,6 +255,11 @@
   ! local vars
   integer          :: i
   character(256)   :: warnstring
+  ! temporary vars for nuc_eos
+  real*8           :: xrho,xye,xtemp,xenr,xent
+  real*8           :: xprs,xmunu,xcs2
+  real*8           :: xdedt,xdpderho,xdpdrhoe
+
   anyerr    = 0
   keyerr(:) = 0
@@ -241,6 +295,32 @@
            dpdepsrho(i) = (hybrid_gamma_th - 1.0d0) * rho(i)
         enddo
+     case (4)
+        ! nuc_eos
+        do i=1,npoints
+
+           xrho = rho(i) * inv_rho_gf
+           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,&
+                keytemp,keyerr(i),rf_precision)
+           
+           if(keyerr(i).ne.0) then
+              anyerr = 1
+           endif
+
+           if(keytemp.eq.1) then
+              eps(i) = xenr * eps_gf
+           else
+              temp(i) = xtemp
+           endif
+
+           dpdepsrho(i) = xdpderho * press_gf * inv_eps_gf 
+
+        enddo
+
      case DEFAULT
         write(warnstring,*) "eoskey ",eoskey," not implemented!"
         call CCTK_WARN(0,warnstring)
@@ -270,6 +350,10 @@
   real*8           :: hybrid_local_gamma, hybrid_local_k_cgs, &
                       hybrid_p_poly, hybrid_p_th
   real*8,parameter :: zero = 0.0d0
+  ! temporary vars for nuc_eos
+  real*8           :: xrho,xye,xtemp,xenr,xent
+  real*8           :: xprs,xmunu,xcs2
+  real*8           :: xdedt
anyerr    = 0
   keyerr(:) = 0
@@ -331,6 +415,35 @@
            cs2(i) = (hybrid_local_gamma * hybrid_p_poly + hybrid_gamma_th * hybrid_p_th) / &
                 rho(i) / (1.0d0 + eps(i) + xpress/rho(i))
         enddo
+     case(4)
+        ! nuc_eos
+
+        do i=1,npoints
+           
+           xrho = rho(i) * inv_rho_gf
+           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,&
+                keytemp,keyerr(i),rf_precision)
+           
+           if(keyerr(i).ne.0) then
+              anyerr = 1
+           endif
+
+           if(keytemp.eq.1) then
+              eps(i) = xenr * eps_gf
+           else
+              temp(i) = xtemp
+           endif
+
+           cs2(i) = xcs2 * cliteinv2 / &
+                (1.0d0 + eps(i) + xprs * press_gf / rho(i))
+
+        enddo
+
+
      case DEFAULT
         write(warnstring,*) "eoskey ",eoskey," not implemented!"
         call CCTK_WARN(0,warnstring)
@@ -357,6 +470,10 @@
   ! local vars
   integer        :: i
   character(256) :: warnstring
+  ! temporary vars for nuc_eos
+  real*8           :: xrho,xye,xtemp,xenr,xent
+  real*8           :: xprs,xmunu,xcs2
+  real*8           :: xdedt,xdpderho,xdpdrhoe
if(keytemp.eq.1) then
      anyerr = 1
@@ -381,6 +498,10 @@
         ! hybrid EOS
         write(warnstring,*) "EOS_Omni_EpsFromPress call not supported for hybrid EOS"
         call CCTK_WARN(0,warnstring)
+     case (4)
+        ! nuc EOS
+        write(warnstring,*) "EOS_Omni_EpsFromPress call not supported for nuc_eos yet"
+        call CCTK_WARN(0,warnstring)
      case DEFAULT
         write(warnstring,*) "eoskey ",eoskey," not implemented!"
         call CCTK_WARN(0,warnstring)
@@ -431,6 +552,10 @@
         ! hybrid EOS
         write(warnstring,*) "EOS_Omni_RestMassDensityFromEpsPress not supported for hybrid EOS"
         call CCTK_WARN(0,warnstring)
+     case (4)
+        ! nuc EOS
+        write(warnstring,*) "EOS_Omni_RestMassDensityFromEpsPress call not supported for nuc_eos yet"
+        call CCTK_WARN(0,warnstring)
      case DEFAULT
         write(warnstring,*) "eoskey ",eoskey," not implemented!"
         call CCTK_WARN(0,warnstring)

    