c -------------------------------------------------------------- #include "cctk.h" #include "cctk_Parameters.h" #include "cctk_Arguments.h" #include "CactusEinstein/Einstein/src/Einstein.h" c -------------------------------------------------------------- ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc ccc ADM_predcorr ccc Perform a 2nd order predictor-corrector step ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc subroutine ADM_predcorr(CCTK_ARGUMENTS) USE ADM_Scalars implicit none DECLARE_CCTK_ARGUMENTS DECLARE_CCTK_PARAMETERS DECLARE_CCTK_FUNCTIONS integer :: i,j,k integer :: ierror CCTK_REAL :: & adm_rhsK_metric_xx,adm_rhsK_metric_xy,adm_rhsK_metric_xz, & adm_rhsK_metric_yy,adm_rhsK_metric_yz,adm_rhsK_metric_zz, & adm_rhsK_matter_xx,adm_rhsK_matter_xy,adm_rhsK_matter_xz, & adm_rhsK_matter_yy,adm_rhsK_matter_yz,adm_rhsK_matter_zz c Temporaries for corrector step CCTK_REAL :: & kxx_n,kxy_n,kxz_n,kyy_n,kyz_n,kzz_n, & gxx_n,gxy_n,gxz_n,gyy_n,gyz_n,gzz_n, & alp_n c Declarations for macros from ADM Evolution #include "CactusEinstein/ADM/src/macro/ALPSOURCES_declare.h" #include "CactusEinstein/ADM/src/macro/GSOURCES_declare.h" #include "CactusEinstein/ADM/src/macro/KSOURCES_declare.h" CCTK_REAL, PARAMETER :: one=1D0 CCTK_REAL :: pi CCTK_REAL :: dx,dy,dz,dt c Set up grid spscings dt=cctk_delta_time dx=cctk_delta_space(1) dy=cctk_delta_space(2) dz=cctk_delta_space(3) pi = ACOS(-1D0) c On iteration zero blank out the sources so we do not get FPEs c when we use (but clobber) them at boundaries later. adms_kxx = 0D0; adms_kxy = 0D0; adms_kxz = 0D0 adms_kyy = 0D0; adms_kyz = 0D0; adms_kzz = 0D0 adms_gxx = 0D0; adms_gxy = 0D0; adms_gxz = 0D0 adms_gyy = 0D0; adms_gyz = 0D0; adms_gzz = 0D0 adms_alp = 0D0 c ================================================================== c c PREDICTOR STEP c c ================================================================== c **CCTK** HACK c CALL ADM_Lapse(CCTK_ARGUMENTS) do k = k1, k2 do j = j1, j2 do i = i1, i2 #include "CactusEinstein/ADM/src/macro/KSOURCES_guts.h" #include "CactusEinstein/ADM/src/macro/GSOURCES_guts.h" if (ADM_evolveLapse) then #include "CactusEinstein/ADM/src/macro/ALPSOURCES_guts.h" end if end do end do end do #include "CactusEinstein/ADM/src/macro/KSOURCES_undefine.h" #include "CactusEinstein/ADM/src/macro/GSOURCES_undefine.h" #include "CactusEinstein/ADM/src/macro/ALPSOURCES_undefine.h" #include "adm_prepred_ADM.inc" #include "adm_preloop_ADM.inc" do k = 1, cctk_lsh(3) do j = 1, cctk_lsh(2) do i = 1, cctk_lsh(1) #include "CactusEinstein/ADM/src/macro/FTCSSTEP_ALP_guts.h" #include "CactusEinstein/ADM/src/macro/FTCSSTEP_G_guts.h" #include "CactusEinstein/ADM/src/macro/FTCSSTEP_K_guts.h" end do end do end do #include "CactusEinstein/ADM/src/macro/FTCSSTEP_ALP_undefine.h" #include "CactusEinstein/ADM/src/macro/FTCSSTEP_G_undefine.h" #include "CactusEinstein/ADM/src/macro/FTCSSTEP_K_undefine.h" CALL CCTK_SyncGroup(ierror,cctkGH,"einstein::lapse") CALL CCTK_SyncGroup(ierror,cctkGH,"einstein::metric") CALL CCTK_SyncGroup(ierror,cctkGH,"einstein::curv") call ADM_Boundary(cctkGH,"einstein::lapse") call ADM_Boundary(cctkGH,"einstein::metric") call ADM_Boundary(cctkGH,"einstein::curv") c ================================================================== c c CORRECTOR STEP STEP c c ================================================================== c **CCTK** HACK c CALL ADM_Lapse(CCTK_ARGUMENTS) do k = k1, k2 do j = j1, j2 do i = i1, i2 #include "CactusEinstein/ADM/src/macro/KSOURCES_guts.h" #include "CactusEinstein/ADM/src/macro/GSOURCES_guts.h" if (ADM_evolveLapse) then #include "CactusEinstein/ADM/src/macro/ALPSOURCES_guts.h" end if end do end do end do #include "CactusEinstein/ADM/src/macro/ALPSOURCES_undefine.h" #include "CactusEinstein/ADM/src/macro/KSOURCES_undefine.h" #include "CactusEinstein/ADM/src/macro/GSOURCES_undefine.h" #include "adm_precorr_ADM.inc" #include "adm_preloop_ADM.inc" do k=1,cctk_lsh(3) do j=1,cctk_lsh(2) do i=1,cctk_lsh(1) gxx_n = gxx(i,j,k)+adms_gxx(i,j,k)*dt gxx(i,j,k) = 0.5D0*(gxx_n+ADM_gxx_p(i,j,k)) gxy_n = gxy(i,j,k)+adms_gxy(i,j,k)*dt gxy(i,j,k) = 0.5D0*(gxy_n+ADM_gxy_p(i,j,k)) gxz_n = gxz(i,j,k)+adms_gxz(i,j,k)*dt gxz(i,j,k) = 0.5D0*(gxz_n+ADM_gxz_p(i,j,k)) gyy_n = gyy(i,j,k)+adms_gyy(i,j,k)*dt gyy(i,j,k) = 0.5D0*(gyy_n+ADM_gyy_p(i,j,k)) gyz_n = gyz(i,j,k)+adms_gyz(i,j,k)*dt gyz(i,j,k) = 0.5D0*(gyz_n+ADM_gyz_p(i,j,k)) gzz_n = gzz(i,j,k)+adms_gzz(i,j,k)*dt gzz(i,j,k) = 0.5D0*(gzz_n+ADM_gzz_p(i,j,k)) kxx_n = kxx(i,j,k)+adms_kxx(i,j,k)*dt kxx(i,j,k) = 0.5D0*(kxx_n+ADM_kxx_p(i,j,k)) kxy_n = kxy(i,j,k)+adms_kxy(i,j,k)*dt kxy(i,j,k) = 0.5D0*(kxy_n+ADM_kxy_p(i,j,k)) kxz_n = kxz(i,j,k)+adms_kxz(i,j,k)*dt kxz(i,j,k) = 0.5D0*(kxz_n+ADM_kxz_p(i,j,k)) kyy_n = kyy(i,j,k)+adms_kyy(i,j,k)*dt kyy(i,j,k) = 0.5D0*(kyy_n+ADM_kyy_p(i,j,k)) kyz_n = kyz(i,j,k)+adms_kyz(i,j,k)*dt kyz(i,j,k) = 0.5D0*(kyz_n+ADM_kyz_p(i,j,k)) kzz_n = kzz(i,j,k)+adms_kzz(i,j,k)*dt kzz(i,j,k) = 0.5D0*(kzz_n+ADM_kzz_p(i,j,k)) alp_n = alp(i,j,k)+dt*adms_alp(i,j,k) alp(i,j,k) = 0.5D0*(alp_n+ADM_alp_p(i,j,k)) end do end do end do CALL CCTK_SyncGroup(ierror,cctkGH,"einstein::lapse") CALL CCTK_SyncGroup(ierror,cctkGH,"einstein::metric") CALL CCTK_SyncGroup(ierror,cctkGH,"einstein::curv") call ADM_Boundary(cctkGH,"einstein::lapse") call ADM_Boundary(cctkGH,"einstein::metric") call ADM_Boundary(cctkGH,"einstein::curv") end subroutine ADM_predcorr