[Commits] [svn:einsteintoolkit] GRHydro_InitData/ (Rev. 163)
rhaas at tapir.caltech.edu
rhaas at tapir.caltech.edu
Mon Sep 17 11:27:23 CDT 2012
User: rhaas
Date: 2012/09/17 11:27 AM
Modified:
/trunk/
param.ccl
/trunk/src/
GRHydro_Bondi_new.F90
Log:
GRHydro_InitData: Move the metric fudge from Bondi routine to SmoothedPuncture thorn
From: Bruno Coutinho Mundim <bcmsma at astro.rit.edu>
Directory Changes:
Directory: /svn:executable/
===========================
+ *
File Changes:
Directory: /trunk/src/
======================
File [modified]: GRHydro_Bondi_new.F90
Delta lines: +11 -63
===================================================================
--- trunk/src/GRHydro_Bondi_new.F90 2012-09-17 16:27:20 UTC (rev 162)
+++ trunk/src/GRHydro_Bondi_new.F90 2012-09-17 16:27:22 UTC (rev 163)
@@ -67,7 +67,7 @@
rs = r_sonicpt_bondi
gam = gl_gamma
- !!write(*,*) 'Bondi_pars:',M,mdot_sonicpt_bondi,r_sonicpt_bondi,gl_gamma
+ write(*,*) 'Bondi_pars:',M,mdot_sonicpt_bondi,r_sonicpt_bondi,gl_gamma
rmin_bondi = M * bondi_rmin(1)
rmax_bondi = M * bondi_rmax(1)
@@ -91,7 +91,7 @@
logrmin = log10(rmin_bondi)
dlogr = (log10(rmax_bondi) - logrmin)/(1.*(N_points-1))
- !!write(*,*)'More pars:',cs,vs,rhos,hs,Kval,Qdot,logrmin,dlogr
+ write(*,*)'More pars:',cs,vs,rhos,hs,Kval,Qdot,logrmin,dlogr
rhotmp=1.0d30
imin=1
@@ -134,10 +134,10 @@
v_bondi(i) = vtmp
enddo
- !!write(*,*)"i=1:",r_bondi(1),rho_bondi(1),u_bondi(1),v_bondi(1)
- !!write(*,*)"i=100:",r_bondi(100),rho_bondi(100),u_bondi(100),v_bondi(100)
- !!write(*,*)"i=1000:",r_bondi(1000),rho_bondi(1000),u_bondi(1000),v_bondi(1000)
- !!write(*,*)"i=1500:",r_bondi(1500),rho_bondi(1500),u_bondi(1500),v_bondi(1500)
+ write(*,*)"i=1:",r_bondi(1),rho_bondi(1),u_bondi(1),v_bondi(1)
+ write(*,*)"i=100:",r_bondi(100),rho_bondi(100),u_bondi(100),v_bondi(100)
+ write(*,*)"i=1000:",r_bondi(1000),rho_bondi(1000),u_bondi(1000),v_bondi(1000)
+ write(*,*)"i=1500:",r_bondi(1500),rho_bondi(1500),u_bondi(1500),v_bondi(1500)
!!$ // find the derivative near r=M
rnew = 2.25 * M
@@ -152,7 +152,7 @@
call find_bondi_solution( rnew, rhocheck2, ucheck2, vcheck2, rs, rhos, M, Mdot, Kval, gam, Qdot )
drhodr = 100.0*(rhocheck2-rhocheck)/M
- !!write(6,*)'Rhocheck:',rhocheck,rhocheck2,drhodr
+ write(6,*)'Rhocheck:',rhocheck,rhocheck2,drhodr
nx = cctk_lsh(1)
ny = cctk_lsh(2)
@@ -172,58 +172,6 @@
zhat = zp/riso
roverm = riso/M
- if(bondi_use_smooth_puncture_data.ne.0 .and. roverm.lt.(bondi_smooth_puncture_radius*ONEmTINY)) then
-
- rsm=bondi_smooth_puncture_radius
- !!write(6,*)'smoothing at a radius',rsm
-
- if(bondi_smooth_gxx.ne.0) then
-
- f=(1.d0+0.5d0/rsm)**4
- df=-1.d0*(2.d0*rsm+1.d0)**3/4.d0/rsm**5
- ddf=(2.d0*rsm+1.d0)**2*(4.d0*rsm+5.d0)/4.d0/rsm**6
-
- a=f-0.625d0*rsm*df+0.125d0*rsm**2*ddf
- b=0.75d0/rsm*df-0.25d0*ddf
- c=-0.125d0/rsm**3*df+0.125d0/rsm**2*ddf
-
- !!write(6,*)'smooth gxx:',rsm,f,df,ddf,a,b,c
-
- gxx(i,j,k) = a + b*roverm**2 + c*roverm**4
- gyy(i,j,k) = a + b*roverm**2 + c*roverm**4
- gzz(i,j,k) = a + b*roverm**2 + c*roverm**4
-
- else
- f=1.d0+0.5d0/rsm
- df=-0.5d0/rsm**2
- ddf=1.d0/rsm**3
-
- a=f-0.625d0*rsm*df+0.125d0*rsm**2*ddf
- b=0.75d0/rsm*df-0.25d0*ddf
- c=-0.125d0/rsm**3*df+0.125d0/rsm**2*ddf
-
- !!write(6,*)'smooth psi:',rsm,f,df,ddf,a,b,c
-
- psinew=a + b*roverm**2 + c*roverm**4
- gxx(i,j,k) = psinew**4
- gyy(i,j,k) = psinew**4
- gzz(i,j,k) = psinew**4
- endif
- gxy(i,j,k) = 0.0
- gxz(i,j,k) = 0.0
- gyz(i,j,k) = 0.0
-
- else if(riso < 1.0e-7) then
- gxx(i,j,k) = 1.0e4
- gyy(i,j,k) = 1.0e4
- gzz(i,j,k) = 1.0e4
- gxy(i,j,k) = 0.0
- gxz(i,j,k) = 0.0
- gyz(i,j,k) = 0.0
- endif
-
- alp(i,j,k) = 1.0/sqrt(gxx(i,j,k))
-
if(roverm > ONEmTINY) then
rsch = 0.25 * ( 2.*riso + M)**2 / riso
jb = floor( 0.5 + (log10(rsch) - logrmin) / dlogr )
@@ -259,10 +207,10 @@
velx(i,j,k),vely(i,j,k),velz(i,j,k), &
eps(i,j,k),press(i,j,k),w_lorentz(i,j,k))
- !! if(riso.gt.1.014d0.and.riso.lt.1.015)write(6,*)'Point to check:', &
- !! x(i,j,k),y(i,j,k),z(i,j,k),riso,gxx(i,j,k),dens(i,j,k),tau(i,j,k),&
- !! sx(i,j,k),sy(i,j,k),sz(i,j,k),rho(i,j,k),eps(i,j,k),&
- !! velx(i,j,k),vely(i,j,k),velz(i,j,k)
+ if(riso.gt.1.014d0.and.riso.lt.1.015)write(6,*)'Point to check:', &
+ x(i,j,k),y(i,j,k),z(i,j,k),riso,gxx(i,j,k),dens(i,j,k),tau(i,j,k),&
+ sx(i,j,k),sy(i,j,k),sz(i,j,k),rho(i,j,k),eps(i,j,k),&
+ velx(i,j,k),vely(i,j,k),velz(i,j,k)
end do
end do
Property changes on: trunk/src/GRHydro_Bondi_new.F90
___________________________________________________________________
Directory: /trunk/
==================
File [modified]: param.ccl
Delta lines: +0 -14
===================================================================
--- trunk/param.ccl 2012-09-17 16:27:20 UTC (rev 162)
+++ trunk/param.ccl 2012-09-17 16:27:22 UTC (rev 163)
@@ -99,20 +99,6 @@
0.0:* :: "Anything positive"
} 1.0
-BOOLEAN bondi_use_smooth_puncture_data "Smooth puncture solution inside R=M using 4th order matching polynomial"
-{
-} "no"
-
-
-BOOLEAN bondi_smooth_gxx "Smooth puncture data by matching to gxx not psi"
-{
-} "yes"
-
-REAL bondi_smooth_puncture_radius "Radius inside which to mathc functions, in units of R/M"
-{
- 0:0.5 :: "Less than/equal to 0.5"
-} 0.5
-
BOOLEAN change_shock_direction "Change the shock direction"
{
} "no"
Property changes on: trunk/param.ccl
___________________________________________________________________
More information about the Commits
mailing list