diff options
author | cott <cott@c83d129a-5a75-4d5a-9c4d-ed3a5855bf45> | 2011-08-02 21:32:00 +0000 |
---|---|---|
committer | cott <cott@c83d129a-5a75-4d5a-9c4d-ed3a5855bf45> | 2011-08-02 21:32:00 +0000 |
commit | 415ab7729def9fad475a4f40e7906ebeb49a29ee (patch) | |
tree | cdf3ba087941dd4909fc382ec8610f4f748c68f2 /src/GRHydro_TmunuM.F90 | |
parent | 4b2e68805d07a963e081b8d4284ee1d6e621efb9 (diff) |
* Optimize: remove support for shift_state = 0 (except for shock tubes and
Cowling calculations of spherically symmetric objects, there is no reason
not to have storage for the shift).
git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinEvolve/GRHydro/trunk@259 c83d129a-5a75-4d5a-9c4d-ed3a5855bf45
Diffstat (limited to 'src/GRHydro_TmunuM.F90')
-rw-r--r-- | src/GRHydro_TmunuM.F90 | 39 |
1 files changed, 12 insertions, 27 deletions
diff --git a/src/GRHydro_TmunuM.F90 b/src/GRHydro_TmunuM.F90 index dfb7ab5..e4df1e0 100644 --- a/src/GRHydro_TmunuM.F90 +++ b/src/GRHydro_TmunuM.F90 @@ -100,35 +100,20 @@ !!$ Calculate lower components and square of shift vector. - if (shift_state .ne. 0) then - - betaxlow = gxx(i,j,k)*betax(i,j,k) + gxy(i,j,k)*betay(i,j,k) + gxz(i,j,k)*betaz(i,j,k) - betaylow = gxy(i,j,k)*betax(i,j,k) + gyy(i,j,k)*betay(i,j,k) + gyz(i,j,k)*betaz(i,j,k) - betazlow = gxz(i,j,k)*betax(i,j,k) + gyz(i,j,k)*betay(i,j,k) + gzz(i,j,k)*betaz(i,j,k) - beta2 = betax(i,j,k)*betaxlow + betay(i,j,k)*betaylow + betaz(i,j,k)*betazlow - - bdotbeta = betaxlow*Bvecx(i,j,k)+betaylow*Bvecy(i,j,k)+betazlow*Bvecz(i,j,k) - vdotbeta = betaxlow*velx(i,j,k)+betaylow*vely(i,j,k)+betazlow*velz(i,j,k) - -!!$ u0 low is missing the w_lorentz factor (see below)!! - utlow = -1.d0*alp(i,j,k) + vdotbeta - - btlow = -1.0d0*w_lorentz(i,j,k)*alp(i,j,k)*bdotv + & - bdotbeta/w_lorentz(i,j,k) + w_lorentz(i,j,k)*bdotv*vdotbeta - - - else - - betaxlow = 0.0D0 - betaylow = 0.0D0 - betazlow = 0.0D0 - beta2 = 0.0D0 + + betaxlow = gxx(i,j,k)*betax(i,j,k) + gxy(i,j,k)*betay(i,j,k) + gxz(i,j,k)*betaz(i,j,k) + betaylow = gxy(i,j,k)*betax(i,j,k) + gyy(i,j,k)*betay(i,j,k) + gyz(i,j,k)*betaz(i,j,k) + betazlow = gxz(i,j,k)*betax(i,j,k) + gyz(i,j,k)*betay(i,j,k) + gzz(i,j,k)*betaz(i,j,k) + beta2 = betax(i,j,k)*betaxlow + betay(i,j,k)*betaylow + betaz(i,j,k)*betazlow + + bdotbeta = betaxlow*Bvecx(i,j,k)+betaylow*Bvecy(i,j,k)+betazlow*Bvecz(i,j,k) + vdotbeta = betaxlow*velx(i,j,k)+betaylow*vely(i,j,k)+betazlow*velz(i,j,k) !!$ u0 low is missing the w_lorentz factor (see below)!! - utlow = -1.0*alp(i,j,k) - btlow = utlow*w_lorentz(i,j,k)*bdotv - - end if + utlow = -1.d0*alp(i,j,k) + vdotbeta + + btlow = -1.0d0*w_lorentz(i,j,k)*alp(i,j,k)*bdotv + & + bdotbeta/w_lorentz(i,j,k) + w_lorentz(i,j,k)*bdotv*vdotbeta !!$ Calculate the specific relativistic enthalpy times rho + the mag. field contribution times the !!$ square of the lorentz factor. |