diff options
author | knarf <knarf@e296648e-0e4f-0410-bd07-d597d9acff87> | 2012-12-19 15:12:36 +0000 |
---|---|---|
committer | knarf <knarf@e296648e-0e4f-0410-bd07-d597d9acff87> | 2012-12-19 15:12:36 +0000 |
commit | 1c980c2cf1278260feb6bb9b613f8af0b22382ce (patch) | |
tree | 2ede115336a741780133ccbceeb823223f939553 /src/gauge.F | |
parent | 076e916c60d9a50dbd84932ae4d891977d21989a (diff) |
Fix compiler warnings.
Most of them could be fixed by renaming .F77 files to .F
Some had to be fixed by explicitly declaring some variables using
CCTK_DECLARE() (which also only works for .F, not for .F77)
git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinInitialData/Exact/trunk@287 e296648e-0e4f-0410-bd07-d597d9acff87
Diffstat (limited to 'src/gauge.F')
-rw-r--r-- | src/gauge.F | 187 |
1 files changed, 187 insertions, 0 deletions
diff --git a/src/gauge.F b/src/gauge.F new file mode 100644 index 0000000..93ceb6a --- /dev/null +++ b/src/gauge.F @@ -0,0 +1,187 @@ +C This routine sets the lapse and/or shift by calling a routine +C that does it pointwise. Note that it could be easily modified +C to set the Bona-Masso variables B_xx etc. +C $Header$ + +#include "cctk.h" +#include "cctk_Parameters.h" +#include "cctk_Arguments.h" +#include "cctk_Functions.h" + + subroutine Exact__gauge(CCTK_ARGUMENTS) + + implicit none + + DECLARE_CCTK_ARGUMENTS + DECLARE_CCTK_PARAMETERS + DECLARE_CCTK_FUNCTIONS + + integer i,j,k + integer nx,ny,nz + logical set_lapse, set_dtlapse, set_shift, set_dtshift + + CCTK_REAL tt, xx, yy, zz + CCTK_REAL gxxtmp, gyytmp, gzztmp, + $ gxytmp, gyztmp, gxztmp, + $ hxxtmp, hyytmp, hzztmp, + $ hxytmp, hyztmp, hxztmp, + $ dxgxxtmp, dxgyytmp, dxgzztmp, + $ dxgxytmp, dxgyztmp, dxgxztmp, + $ dygxxtmp, dygyytmp, dygzztmp, + $ dygxytmp, dygyztmp, dygxztmp, + $ dzgxxtmp, dzgyytmp, dzgzztmp, + $ dzgxytmp, dzgyztmp, dzgxztmp, + $ alptmp, dtalptmp, axtmp, aytmp, aztmp, + $ betaxtmp, betaytmp, betaztmp, + $ dtbetaxtmp, dtbetaytmp, dtbetaztmp, + $ bxxtmp, bxytmp, bxztmp, + $ byxtmp, byytmp, byztmp, + $ bzxtmp, bzytmp, bzztmp + CCTK_REAL + $ exact_psi, + $ exact_psix, exact_psiy, exact_psiz, + $ exact_psixx, exact_psiyy, exact_psizz, + $ exact_psixy, exact_psiyz, exact_psixz + LOGICAL is_initial_slice, is_later_slice + +C are we on the initial slice or some later slice? +C n.b. the logical expressions later in this function involving +C these flags below would be *so* much nicer if Fortran +C grokked C-style conditional expressions... :) :) + is_initial_slice = cctk_iteration .eq. 0 + is_later_slice = cctk_iteration .ne. 0 + +C Grid parameters. + nx = cctk_lsh(1) + ny = cctk_lsh(2) + nz = cctk_lsh(3) + +C This code used to set t = time + dt/2 to get 2nd order accuracy, +C but this leads to the initial data being set at the wrong time. :( +C In the context of MoL, we want to set variables at the standard Cactus +C time (cctk_time), because MoL takes care of calling us at each MoL +C iteration, and updating the field variables appropriately. +C +C Alas, setting at cctk_time probably gives O(dt) errors for non-MoL +C evoutions where Exact is used to set stuff at each time step. +C Fixing this [unless we just declare all non-MoL stuff obselete :) ] +C probably requires cleaning up our (++messy) schedule.ccl , which +C is why this remains a bug for now... :( :( + +CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC + +C +C Set lapse and/or shift? +C + if ( is_initial_slice ) then + set_lapse = CCTK_Equals(initial_lapse, "exact").ne.0 + set_shift = CCTK_Equals(initial_shift, "exact").ne.0 + set_dtlapse = CCTK_Equals(initial_dtlapse, "exact").ne.0 + set_dtshift = CCTK_Equals(initial_dtshift, "exact").ne.0 + end if + if ( is_later_slice ) then + set_lapse = CCTK_Equals(lapse_evolution_method, "exact").ne.0 + set_shift = CCTK_Equals(shift_evolution_method, "exact").ne.0 + set_dtlapse = CCTK_Equals(dtlapse_evolution_method, "exact").ne.0 + set_dtshift = CCTK_Equals(dtshift_evolution_method, "exact").ne.0 + end if + + if ( set_lapse .or. set_shift .or. set_dtlapse .or. set_dtshift) then + +C$omp parallel do private ( +C$omp$ i, j, k, +C$omp$ tt, xx, yy, zz, +C$omp$ alptmp, dtalptmp, axtmp, aytmp, aztmp, +C$omp$ betaxtmp, betaytmp, betaztmp, +C$omp$ dtbetaxtmp, dtbetaytmp, dtbetaztmp, +C$omp$ bxxtmp, bxytmp, bxztmp, +C$omp$ byxtmp, byytmp, byztmp, +C$omp$ bzxtmp, bzytmp, bzztmp, +C$omp$ dxgxxtmp, dxgyytmp, dxgzztmp, +C$omp$ dxgxytmp, dxgyztmp, dxgxztmp, +C$omp$ dygxxtmp, dygyytmp, dygzztmp, +C$omp$ dygxytmp, dygyztmp, dygxztmp, +C$omp$ dzgxxtmp, dzgyytmp, dzgzztmp, +C$omp$ dzgxytmp, dzgyztmp, dzgxztmp, +C$omp$ exact_psi, +C$omp$ exact_psix, exact_psiy, exact_psiz, +C$omp$ exact_psixx, exact_psiyy, exact_psizz, +C$omp$ exact_psixy, exact_psiyz, exact_psixz) + do k=1,nz + do j=1,ny + do i=1,nx + + tt = cctk_time + xx = x(i,j,k) - cctk_time * shift_add_x + yy = y(i,j,k) - cctk_time * shift_add_y + zz = z(i,j,k) - cctk_time * shift_add_z + +C Initialize the psi of exact +C (also to tell the models about the conformal_state) + if (conformal_state .ne. 0) then + exact_psi = 1.0D0 + else + exact_psi = 0.0D0 + end if + exact_psix = 0.0D0 + exact_psiy = 0.0D0 + exact_psiz = 0.0D0 + exact_psixx = 0.0D0 + exact_psiyy = 0.0D0 + exact_psizz = 0.0D0 + exact_psixy = 0.0D0 + exact_psiyz = 0.0D0 + exact_psixz = 0.0D0 + + call Exact__Bona_Masso_data( + $ decoded_exact_model, + $ xx, yy, zz, tt, + $ gxxtmp, gyytmp, gzztmp, + $ gxytmp, gyztmp, gxztmp, + $ hxxtmp, hyytmp, hzztmp, + $ hxytmp, hyztmp, hxztmp, + $ exact_psi, + $ exact_psix, exact_psiy, exact_psiz, + $ exact_psixx, exact_psiyy, exact_psizz, + $ exact_psixy, exact_psiyz, exact_psixz, + $ dxgxxtmp, dxgyytmp, dxgzztmp, + $ dxgxytmp, dxgyztmp, dxgxztmp, + $ dygxxtmp, dygyytmp, dygzztmp, + $ dygxytmp, dygyztmp, dygxztmp, + $ dzgxxtmp, dzgyytmp, dzgzztmp, + $ dzgxytmp, dzgyztmp, dzgxztmp, + $ alptmp, dtalptmp, axtmp, aytmp, aztmp, + $ betaxtmp, betaytmp, betaztmp, + $ dtbetaxtmp, dtbetaytmp, dtbetaztmp, + $ bxxtmp, bxytmp, bxztmp, + $ byxtmp, byytmp, byztmp, + $ bzxtmp, bzytmp, bzztmp) + + if ( set_lapse ) then + alp(i,j,k) = alptmp + end if + if ( set_shift ) then + betax(i,j,k) = betaxtmp + shift_add_x + betay(i,j,k) = betaytmp + shift_add_y + betaz(i,j,k) = betaztmp + shift_add_z + end if + if ( set_dtlapse ) then + dtalp(i,j,k) = dtalptmp + end if + if ( set_dtshift ) then + dtbetax(i,j,k) = dtbetaxtmp + dtbetay(i,j,k) = dtbetaytmp + dtbetaz(i,j,k) = dtbetaztmp + end if + end do + end do + end do + +CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC + + else + call CCTK_WARN(1,'Exact__gauge has been called without doing anything') + end if + + return + end |