User: rhaas Date: 2012/08/09 01:26 AM

Modified: /trunk/src/ GRHydro_AdvectedLoopM.F90

Log: Revised parameters for advected loop test

From: Josh Faber jfaberrit@jafsmamacbook.main.ad.rit.edu

File Changes:

Directory: /trunk/src/ ======================

File [modified]: GRHydro_AdvectedLoopM.F90 Delta lines: +60 -12 =================================================================== --- trunk/src/GRHydro_AdvectedLoopM.F90 2012-08-09 06:26:42 UTC (rev 141) +++ trunk/src/GRHydro_AdvectedLoopM.F90 2012-08-09 06:26:44 UTC (rev 142) @@ -67,7 +67,7 @@ CCTK_REAL :: range_x,range_y,range_z,range_d CCTK_REAL :: cos_theta, sin_theta, tan_theta CCTK_REAL :: Bvecx_d, Bvecy_d, Bvecz_d - CCTK_REAL :: x_d, y_d, z_d + CCTK_REAL :: x_d, y_d, z_d,diaglength CCTK_REAL :: dx_d, dy_d, dz_d, dx_x, dz_x

!!$Adiabatic index for test: @@ -83,19 +83,47 @@ rhoval = 1.0d0 pressval = 3.0d0

+ if (CCTK_EQUALS(advectedloop_type,"2D")) then + + !!$ Vx, Vy and Vz values: - if (CCTK_EQUALS(advectedloop_case,"V^z/=0")) then - vxval=0.2d0/sqrt(6.0d0) - vyval=0.5d0*vxval - vzval=vyval - else if (CCTK_EQUALS(advectedloop_case,"V^z=0")) then - vxval=0.2d0/sqrt(6.0d0) - vyval=0.5d0*vxval - vzval=0.0d0 - else - call CCTK_WARN(0,"V^z component case not recognized!") - end if + if (CCTK_EQUALS(advectedloop_case,"V^z/=0")) then + +!!$vxval=0.2d0/sqrt(6.0d0) +!!$ This new choice yields a crossing time of exactly t=24 +!!$ assuming -1<x<1; -0.5<y<0.5 + vxval=1.d0/12.d0 + vyval=0.5d0*vxval + vzval=vyval + + else if (CCTK_EQUALS(advectedloop_case,"V^z=0")) then +!!$ vxval=0.2d0/sqrt(6.0d0) + vxval=1.d0/12.d0 + vyval=0.5d0*vxval + vzval=0.0d0 + else + call CCTK_WARN(0,"V^z component case not recognized!") + end if

+ else if (CCTK_EQUALS(advectedloop_type,"3D")) then + + vxval=0.2d0*sqrt(2.0d0) + vyval=0.2d0 + + if (CCTK_EQUALS(advectedloop_case,"V^z/=0")) then + +!!$vxval=0.2d0/sqrt(6.0d0) +!!$ This new choice yields a crossing time of exactly t=5 +!!$ assuming -1<x<1; -0.5<y<0.5 + vzval=0.1 + + else if (CCTK_EQUALS(advectedloop_case,"V^z=0")) then + vzval=0.0d0 + else + call CCTK_WARN(0,"V^z component case not recognized!") + end if + + endif nx = cctk_lsh(1) ny = cctk_lsh(2) nz = cctk_lsh(3) @@ -110,7 +138,10 @@ range_z = (cctk_gsh(3)-1)*dz

range_d = sqrt(range_z**2+range_x**2) + diaglength = range_x*range_z/range_d

+!!$ For 3-d case, the grid is going to be assumed to be a cube + cos_theta = range_z/range_d sin_theta = range_x/range_d tan_theta = sin_theta/cos_theta @@ -162,16 +193,33 @@

else if (CCTK_EQUALS(advectedloop_type,"3D")) then

+ +!!$ tangential velocity should be parallel to (1,1,1) plus a normal component (-1,0,1) + velx(i,j,k)=cos_theta*vxval-sin_theta*vzval vely(i,j,k)=vyval velz(i,j,k)=cos_theta*vzval+sin_theta*vxval

Bvecz_d=0.0d0

+!!$ x_d = (x+z)/sqrt(2) => x_d=0 is equivalent to x+z=0 + x_d = cos_theta*x(i,j,k)+sin_theta*z(i,j,k) y_d = y(i,j,k) z_d = cos_theta*z(i,j,k)-sin_theta*x(i,j,k)

+!!$ need to make x_d periodic! + + if(x_d.gt.1.5*diaglength) then + x_d=x_d-2.0*diaglength + else if (x_d.gt.0.5*diaglength .and. x_d.lt.1.5*diaglength) then + x_d=x_d-diaglength + else if(x_d.lt.-1.5*diaglength) then + x_d=x_d+2.0*diaglength + else if (x_d.lt.(-0.5*diaglength) .and. x_d.gt.(-1.5*diaglength)) then + x_d=x_d+diaglength + endif + radius = sqrt(x_d**2+y_d**2)

if (CCTK_EQUALS(advectedloop_delA,"Exact")) then