/*@@ @file CartGrid3D.F77 @date Thu Feb 18 @author Gabrielle Allen @desc Set up coordinates for a 3D Cartesian grid @enddesc @@*/ #define CARTGRID3D_DEBUG #include "cctk.h" #include "cctk_arguments.h" #include "cctk_parameters.h" subroutine CartGrid3D(CCTK_FARGUMENTS) implicit none DECLARE_CCTK_FARGUMENTS DECLARE_CCTK_PARAMETERS integer CCTK_Equals integer iconv, i, j, k CCTK_REAL x_origin,y_origin,z_origin CCTK_REAL this_dx,this_dy,this_dz CCTK_REAL xmax1,ymax1,zmax1 CCTK_REAL lower,upper,olower,oupper iconv = 2**(cctk_convlevel-1) c -------------------------------------------------------------- c c BYRANGE c c User gives: minimum and maximum values of coordinates and c the number of gridpoints on the coarse grid c c -------------------------------------------------------------- if (CCTK_Equals(type,'byrange') .eq. 1) then if (xyzmax .gt. 0) then xmax1 = xyzmax ymax1 = xyzmax zmax1 = xyzmax else xmax1 = xmax ymax1 = ymax zmax1 = zmax end if if (CCTK_EQUALS(symmetry,'octant')) then c Grid spacing on coarsest grid if (no_origin.eq.1) then coarse_dx = xmax/(DBLE(cctk_gsh(1))-1.5d0) coarse_dy = ymax/(DBLE(cctk_gsh(2))-1.5d0) coarse_dz = zmax/(DBLE(cctk_gsh(3))-1.5d0) else coarse_dx = xmax/(DBLE(cctk_gsh(1))-1.0d0) coarse_dy = ymax/(DBLE(cctk_gsh(2))-1.0d0) coarse_dz = zmax/(DBLE(cctk_gsh(3))-1.0d0) end if c Grid spacing on this grid this_dx = coarse_dx/cctk_levfac(1) this_dy = coarse_dy/cctk_levfac(2) this_dz = coarse_dz/cctk_levfac(3) c Minimum coordinate values on this grid if (no_origin.eq.1) then x_origin = (0.5-dble(cctk_nghostzones(1)))*this_dx y_origin = (0.5-dble(cctk_nghostzones(2)))*this_dy z_origin = (0.5-dble(cctk_nghostzones(3)))*this_dz else x_origin = (-dble(cctk_nghostzones(1)))*this_dx y_origin = (-dble(cctk_nghostzones(2)))*this_dy z_origin = (-dble(cctk_nghostzones(3)))*this_dz end if else if (CCTK_EQUALS(symmetry,'quadrant')) then print *,"FIXME" else if (CCTK_Equals(symmetry,'bitant') .eq. 1) then print *,"FIXME" else if (CCTK_Equals(symmetry,'full') .eq. 1) then c Set minimum values of coordinates x_origin = xmin y_origin = ymin z_origin = zmin c dx,dy,dz on the coarsest grid of each GH coarse_dx = (xmax-xmin)/max(cctk_gsh(1)-1,1) coarse_dy = (ymax-ymin)/max(cctk_gsh(2)-1,1) coarse_dz = (zmax-zmin)/max(cctk_gsh(3)-1,1) end if c dx,dy,dz on the grid we are on this_dx = coarse_dx*iconv/cctk_levfac(1) this_dy = coarse_dy*iconv/cctk_levfac(2) this_dz = coarse_dz*iconv/cctk_levfac(3) write(*,'(1X,A,1X,3(A,G12.7,2X))') & 'Grid by range: Setting ','dx=>',this_dx,'dy=>', & this_dy,'dz=>',this_dz c ----------------------------------------------------------- c c BYSPACING c c User gives: grid spacing on the coarsest GH and c the number of gridpoints on the coarsest GH c c ----------------------------------------------------------- else if (CCTK_Equals(type,'byspacing') .eq. 1) then c Dx, Dy, Dx on the coarsest grid if (dxyz .gt. 0) then coarse_dx = dxyz coarse_dy = dxyz coarse_dz = dxyz else coarse_dx = dx coarse_dy = dy coarse_dz = dz end if c dx, dy, dz on the grid we are on this_dx = coarse_dx*iconv/cctk_levfac(1) this_dy = coarse_dy*iconv/cctk_levfac(2) this_dz = coarse_dz*iconv/cctk_levfac(3) if (CCTK_Equals(symmetry,'bitant') .eq. 1) then if (no_origin.eq.1) then x_origin = (0.5 - cctk_gsh(1)/2)*this_dx y_origin = (0.5 - cctk_gsh(2)/2)*this_dy z_origin = (-dble(cctk_nghostzones(3))+0.5d0)*this_dz else x_origin = (- cctk_gsh(1)/2)*this_dx y_origin = (- cctk_gsh(2)/2)*this_dy z_origin = (-dble(cctk_nghostzones(3)))*this_dz end if else if (CCTK_Equals(symmetry,'quadrant') .eq. 1) then if (no_origin.eq.1) then x_origin = (-dble(cctk_nghostzones(1))+0.5d0)*this_dx y_origin = (-dble(cctk_nghostzones(2))+0.5d0)*this_dy z_origin = (0.5D0 - cctk_gsh(3)/2)*this_dz else x_origin = (-dble(cctk_nghostzones(1)))*this_dx y_origin = (-dble(cctk_nghostzones(2)))*this_dy z_origin = (- cctk_gsh(3)/2)*this_dz end if else if (CCTK_Equals(symmetry,'octant') .eq. 1) then if (no_origin.eq.1) then x_origin = (-dble(cctk_nghostzones(1))+0.5d0)*this_dx y_origin = (-dble(cctk_nghostzones(2))+0.5d0)*this_dy z_origin = (-dble(cctk_nghostzones(3))+0.5d0)*this_dz else x_origin = (-dble(cctk_nghostzones(1)))*this_dx y_origin = (-dble(cctk_nghostzones(2)))*this_dy z_origin = (-dble(cctk_nghostzones(3)))*this_dz end if else if (CCTK_Equals(symmetry,'full').eq.1) then if (no_origin.eq.1) then x_origin = (0.5 - cctk_gsh(1)/2)*this_dx y_origin = (0.5 - cctk_gsh(2)/2)*this_dy z_origin = (0.5 - cctk_gsh(3)/2)*this_dz else x_origin = (- cctk_gsh(1)/2)*this_dx y_origin = (- cctk_gsh(2)/2)*this_dy z_origin = (- cctk_gsh(3)/2)*this_dz end if end if c -------------------------------------------------------------- c c BOX (-0.5 to 0.5) c c User gives: number of gridpoints on the coarse grid c c -------------------------------------------------------------- elseif (CCTK_Equals(type,'box') .eq. 1) then c Coordinates are all -0.5 to 0.5 x_origin = -0.5 y_origin = -0.5 z_origin = -0.5 c dx,dy,dz on the coarsest grid of each GH coarse_dx = 1.d0/max(cctk_gsh(1)-1,1) coarse_dy = 1.d0/max(cctk_gsh(1)-1,1) coarse_dz = 1.d0/max(cctk_gsh(1)-1,1) c dx,dy,dz on the grid we are on this_dx = coarse_dx*iconv/cctk_levfac(1) this_dy = coarse_dy*iconv/cctk_levfac(2) this_dz = coarse_dz*iconv/cctk_levfac(3) c Special cases if (cctk_gsh(1) .eq. 1) x_origin = 0.0D0 if (cctk_gsh(2) .eq. 1) y_origin = 0.0D0 if (cctk_gsh(3) .eq. 1) z_origin = 0.0D0 write(*,'(1X,A,1X,3(A,G12.7,2X))') & 'Box grid: Setting ','dx=>',this_dx,'dy=>',this_dy, & 'dz=>',this_dz end if c Set spatial coordinates c ----------------------- do i=1,cctk_lsh(1) do j=1,cctk_lsh(2) do k=1,cctk_lsh(3) x(i,j,k) = this_dx*(i+cctk_lbnd(1)-1) + x_origin y(i,j,k) = this_dy*(j+cctk_lbnd(2)-1) + y_origin z(i,j,k) = this_dz*(k+cctk_lbnd(3)-1) + z_origin r(i,j,k)=sqrt(x(i,j,k)**2+y(i,j,k)**2+z(i,j,k)**2) end do end do end do cctk_delta_space(1) = this_dx cctk_delta_space(2) = this_dy cctk_delta_space(3) = this_dz cctk_origin_space(1) = x_origin cctk_origin_space(2) = y_origin cctk_origin_space(3) = z_origin lower = x_origin upper = x_origin+this_dx*(cctk_gsh(1)-1) call CCTK_RegisterCoordRange(cctkGH,lower,upper,"x") lower = y_origin upper = y_origin+this_dy*(cctk_gsh(2)-1) call CCTK_RegisterCoordRange(cctkGH,lower,upper,"y") lower = z_origin upper = z_origin+this_dz*(cctk_gsh(3)-1) call CCTK_RegisterCoordRange(cctkGH,lower,upper,"z") #ifdef CARTGRID3D_DEBUG write(*,*) write(*,*) " CartGrid3D" write(*,*) " ----------" write(*,100) coarse_dx,coarse_dy,coarse_dz 100 format(2X,"Dx,Dy,Dz on coarse grid ",f6.3,f6.3,f6.3) write(*,101) cctk_delta_space(1),cctk_delta_space(2),cctk_delta_space(3) 101 format(2X,"Dx,Dy,Dz on this grid ",f6.3,f6.3,f6.3) write(*,*) " Convergence level = ",cctk_convlevel write(*,*) " Grid level = ",cctk_levfac(1),cctk_levfac(2),cctk_levfac(3) write(*,102) x_origin,y_origin,z_origin 102 format(2X,"Minimum global coords ",f6.3,f6.3,f6.3) write(*,103) x_origin+this_dx*(cctk_gsh(1)-1), & y_origin+this_dy*(cctk_gsh(2)-1),z_origin+this_dz*(cctk_gsh(3)-1) 103 format(2X,"Maximum global coords ",f6.3,f6.3,f6.3) write(*,104) x(1,1,1),y(1,1,1),z(1,1,1) 104 format(2X,"Minimum local coords ",f6.3,f6.3,f6.3) write(*,105) x(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3)), & y(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3)),z(cctk_lsh(1),cctk_lsh(2),cctk_lsh(3)) 105 format(2X,"Maximum local coords ",f6.3,f6.3,f6.3) #endif end