/*@@ @file ADM_PreloopF @date Aug 1998 @author Gabrielle Allen @desc Report on parameters used for ADM_Evolution @enddesc @@*/ #include "cctk.h" #include "cctk_Parameters.h" #include "cctk_Arguments.h" /*@@ @routine ADM_evolveLapseFlag @date June 2000 @author Gerd Lanfermann @desc Routine which checks the current slicing (announced in the grid scalar "active_slicing_handle" and set a flag which tells the evolution code whether it needs to calc the lapse or not. @enddesc @@*/ subroutine ADM_evolveLapseFlag(active_slicing_handle) USE ADM_Scalars implicit none CCTK_INT active_slicing_handle INTEGER Einstein_GetSlicingHandle if ((active_slicing_handle .eq. Einstein_GetSlicingHandle("maximal")) .OR. $ (active_slicing_handle .eq. Einstein_GetSlicingHandle("static")).OR. $ (active_slicing_handle .eq. Einstein_GetSlicingHandle("geodesic"))) then ADM_evolveLapse = .FALSE. else ADM_evolveLapse = .TRUE. endif return end /*@@ @routine ADM_Preloop @date Aug 1998 @author Gabrielle Allen @desc Report on parameters used for ADM_Evolution Basically just a copy of ADM_Lapse.F with output @enddesc @calls @history @hauthor Gerd Lanfermann @hdate July 99 @hdesc For the slicing, the routine have been changed to use the slicing handles. Read ./Einstein/Slicing.c for the details. There is a grid scalar called "active_slicing_handle" which defines what slicing to use in the next iteration(it's a number). Every slicing checks against this number by calling GetSlicingHandle, if it matches, it does the slicing. This basically substituted the prev. CCTK_Equals(slicing," ") comparison. To support mixed slicings (changing the slicing during an evoltion) PreLoop needs to evaluate the current slicing by calling ADM_evolveLapseFlag. @endhistory @@*/ subroutine ADM_Setup(CCTK_ARGUMENTS) USE ADM_Scalars implicit none DECLARE_CCTK_ARGUMENTS DECLARE_CCTK_PARAMETERS c Perhaps this and others should go into cctk.h INTEGER CCTK_Equals INTEGER Einstein_GetSlicingHandle c Diffusion coefficients for lapse evolution CCTK_REAL :: alp_Geps,alp_nablaeps CCTK_REAL :: dx,dy,dz,dt CHARACTER*200 :: infoline CHARACTER*200 :: fortslicing CCTK_INT :: nslicing integer,save:: prev_slicing_handle = -1 integer new_slice c Slicing management c ---------------------------------------------------- c Set the flag, if a slicing needs to calculated or not call ADM_evolveLapseFlag(active_slicing_handle) c Check if slicing has changed from prev. calls if (prev_slicing_handle.ne.active_slicing_handle) then new_slice = 1 prev_slicing_handle = active_slicing_handle else new_slice = 0 endif c Get a fortran string for the slicing call CCTK_FortranString(nslicing,slicing,fortslicing) c Set up the limits of loops for the evolution methods c ---------------------------------------------------- c (stencil width = 1) i1 = 2 i2 = cctk_lsh(1)-1 j1 = 2 j2 = cctk_lsh(2)-1 k1 = 2 k2 = cctk_lsh(3)-1 c We can switch this off to avoid evolving the lapse in the ADM routines c ---------------------------------------------------------------------- dt = cctk_delta_time dx = cctk_delta_space(1) dy = cctk_delta_space(2) dz = cctk_delta_space(3) c GEODESIC SLICING c ---------------- if (active_slicing_handle .eq. Einstein_GetSlicingHandle("geodesic")) then if (new_slice.eq.1) then call CCTK_INFO("Geodesic slicing (lapse set to unity)") call CCTK_INFO(" geodesic lapse will not be evolved in ADM routines") endif c STATIC SLICING c -------------- else if (active_slicing_handle .eq. Einstein_GetSlicingHandle("static")) then if (new_slice.eq.1) then call CCTK_INFO("Static slicing (lapse fixed at initial values)") call CCTK_INFO(" lapse will not be evolved in ADM routines") endif c 1+LOG SLICING c ------------- else if (active_slicing_handle .eq. Einstein_GetSlicingHandle("1+log")) then if (new_slice.eq.1) then call CCTK_INFO("1+log slicing") endif alp_Geps = flat_dissip alp_nablaeps = dissip c In the NCSA style if (CCTK_Equals(slicing_flavour,'kleing') == 1) then if (new_slice.eq.1) then call CCTK_INFO(" in the kleing style") if (store_initial_lapse == 1) then call CCTK_INFO(" using the initial lapse coeff") else call CCTK_INFO(" but not using the initial lapse coeff") endif endif ADM_flatNablaAlpCoeff = alp_Geps/Dt ADM_nablaAlpCoeff = 0D0 if (new_slice.eq.1) then if (ADM_flatNablaAlpCoeff .NE. 0D0) THEN write(infoline,'(A21,G12.7)') & 'flat space diffusion ',real(ADM_flatNablaAlpCoeff) call CCTK_INFO(infoline) else call CCTK_INFO(' with no flat space diffusion') end if call CCTK_INFO( ' and no curved space diffusion') endif else ADM_nablaAlpCoeff = alp_nablaeps*dx**2/dt*0.5D0 ADM_flatNablaAlpCoeff = alp_Geps*dx**2/dt*0.5D0 if (new_slice.eq.1) then if (ADM_flatNablaAlpCoeff .NE. 0D0) then write(infoline,'(A21,G12.7)') & 'flat space diffusion ',real(ADM_flatNablaAlpCoeff) call CCTK_INFO(infoline) end if if (ADM_NablaAlpCoeff .NE. 0D0) then write(infoline,'(A25,G12.7)') & ' curved space diffusion ',real(ADM_NablaAlpCoeff) call CCTK_INFO(infoline) endif if (alp_nablaeps .gt. 0.5D0*dx**2/dt) then call CCTK_INFO(" you are probably too diffusive") else if (alp_nablaeps .lt. 1e-10) then call CCTK_INFO(" curved space diffusion = 0 for efficiency") end if endif endif c BONA-MASSO TYPE HARMONIC SLICINGS c --------------------------------- else if (active_slicing_handle .eq. Einstein_GetSlicingHandle("harmonic")) then if (new_slice.eq.1) then call CCTK_INFO("Bona-Masso family slicing") write(infoline,'(A10,G12.7)') & ' with f = ',real(harmonic_f) call CCTK_INFO(infoline) endif c NCSA TYPE SLICING WITH g(gamma)=sqrt(det g) c ------------------------------------------- else if (active_slicing_handle .eq. Einstein_GetSlicingHandle("detg")) then if (new_slice.eq.1) then call CCTK_INFO("detg ( = plain harmonic) slicing") endif c EXTERNAL SLICING c ---------------- else if ((active_slicing_handle .eq. Einstein_GetSlicingHandle("external"))) then if (new_slice.eq.1) then call CCTK_INFO("External slicing") call CCTK_INFO(" Lapse will not be evolved in ADM routines") endif c EXTERNAL SLICING c ---------------- else if ((active_slicing_handle .eq. Einstein_GetSlicingHandle("maximal"))) then if (new_slice.eq.1) then call CCTK_INFO("Maximal slicing") call CCTK_INFO(" Lapse will not be evolved in ADM routines") endif c SLICING CONDITION NOT SUPPLIED OR MISSPELT c ------------------------------------------ else if (Einstein_GetSlicingHandle(fortslicing(1:nslicing)).lt.0) then call CCTK_WARN(1,"WARNING in thorn ADM: slicing in parfile not registered by ADM") else call CCTK_WARN(1,"WARNING slicing registered but not coded!") end if return end