{"owner":"eschnett","repo":"cactusamrex","channel":"main","channels":["main","master"],"can_delete":true,"inode":{"File":{"path":[{"basename":"CartGrid3D","pos":"ANFQ3BZCQVIJK.BUAAA","meta":1005},{"basename":"src","pos":"ANFQ3BZCQVIJK.YOYAC","meta":1005},{"basename":"CartGrid3D.c","pos":"ANFQ3BZCQVIJK.2OYAC","meta":420}],"file":"/*@@\n @file CartGrid3D.c\n @date Thu Oct 7 13:20:06 1999\n @author Tom Goodale\n @desc\n Set up coordinates for a 3D Cartesian grid.\n C Conversion of Fortran routine written by Gab.\n @version $Id$\n @enddesc\n@@*/\n\n#include \n#include \n\n#include \"cctk.h\"\n#include \"cctk_Arguments.h\"\n#include \"cctk_Parameters.h\"\n#include \"util_Table.h\"\n\n#include \"Symmetry.h\"\n#include \"CoordBase.h\"\n#include \"loopcontrol.h\"\n\nstatic const char *rcsid = \"$Header$\";\n\nCCTK_FILEVERSION(CactusBase_CartGrid3D_CartGrid3D_c);\n\n/********************************************************************\n ********************* Macro Definitions **************************\n ********************************************************************/\n#define max(a, b) ((a) > (b) ? (a) : (b))\n#define SQR(a) ((a) * (a))\n\n/********************************************************************\n ********************* Scheduled Routine Prototypes ***************\n ********************************************************************/\nvoid CartGrid3D_SetRanges(CCTK_ARGUMENTS);\nvoid CartGrid3D_SetCoordinates(CCTK_ARGUMENTS);\n\n/********************************************************************\n ********************* External Routine Prototypes ****************\n ********************************************************************/\nvoid DecodeSymParameters3D(int sym[6]);\n\n/********************************************************************\n ********************* Local Routine Prototypes *******************\n ********************************************************************/\n\n/********************************************************************\n ********************* Scheduled Routines *************************\n ********************************************************************/\n/*@@\n @routine CartGrid3D_SetRanges\n @date Oct 1999?\n @author Tom Goodale? Gabrielle Allen? Thomas Radke\n @date Mon 3 Jan 2005\n @author Thomas Radke\n @desc\n Sets up ranges for Cartesian coordinates.\n @enddesc\n @calls DecodeSymParameters3D, CCTK_Equals, CCTK_WARN,\n CCTK_CoordRegisterRange, CCTK_CoordRegisterRangePhysIndex,\n CCTK_INFO, Coord_CoordHandle, Util_TableSetReal,\n Util_TableSetString, Util_TableSetInt\n\n @var CCTK_ARGUMENTS\n @vdesc Cactus argument list\n @vtype\n @vio in/out\n @endvar\n@@*/\nvoid CartGrid3D_SetRanges(CCTK_ARGUMENTS) {\n int i;\n int coord_handle, ierr;\n CCTK_REAL this_delta[3], origin[3], min1[3], max1[3];\n CCTK_REAL *coarse_delta[3];\n double lower[3], upper[3];\n int domainsym[6], cntstag[3], loweri[3], upperi[3], do_periodic[3];\n char coord_name[16];\n DECLARE_CCTK_ARGUMENTS;\n DECLARE_CCTK_PARAMETERS;\n\n if (CCTK_EQUALS(set_coordinate_ranges_on, \"first level\")) {\n /* Ranges must be set up only once, and this must happen on the\n coarse grid. However, the coarse grid itself may not actually\n exist; in this case, use the coarsest existing grid. We assume\n that this is the first grid for which this routine is\n called. */\n static int is_coarsest_refinement_level = 1;\n if (!is_coarsest_refinement_level) {\n return;\n }\n is_coarsest_refinement_level = 0;\n } else {\n /* Ranges need to be set up only once (or once per map), on the\n coarsest refinement level */\n int const is_coarsest_refinement_level =\n cctk_levfac[0] == 1 && cctk_levfac[1] == 1 && cctk_levfac[2] == 1;\n if (!is_coarsest_refinement_level) {\n return;\n }\n }\n\n coarse_delta[0] = coarse_dx;\n coarse_delta[1] = coarse_dy;\n coarse_delta[2] = coarse_dz;\n\n /* Calculate the coordinate ranges only for the coarsest level */\n /* Avoid origin? Default is yes */\n cntstag[0] = no_origin && no_originx && avoid_origin && avoid_originx;\n cntstag[1] = no_origin && no_originy && avoid_origin && avoid_originy;\n cntstag[2] = no_origin && no_originz && avoid_origin && avoid_originz;\n\n /* Determine symmetries of domain */\n DecodeSymParameters3D(domainsym);\n\n do_periodic[0] = periodic && periodic_x;\n do_periodic[1] = periodic && periodic_y;\n do_periodic[2] = periodic && periodic_z;\n\n /* Calculate physical indices, using symmetries and periodicity */\n for (i = 0; i < 3; i++) {\n loweri[i] = 0;\n upperi[i] = cctk_gsh[i] - 1;\n\n if (domainsym[2 * i + 0] || do_periodic[i]) {\n loweri[i] += cctk_nghostzones[i];\n }\n if (domainsym[2 * i + 1] || do_periodic[i]) {\n upperi[i] -= cctk_nghostzones[i];\n }\n }\n\n /****************************************************************\n *\n * BYRANGE\n *\n * User gives: minimum and maximum values of coordinates and\n * the number of gridpoints on the coarse grid\n *\n ***************************************************************/\n /**************************************************************\n *\n * BOX (-0.5 to 0.5)\n *\n * User gives: number of gridpoints on the coarse grid\n *\n **************************************************************/\n\n if (CCTK_Equals(type, \"byrange\") || CCTK_Equals(type, \"box\")) {\n if (CCTK_Equals(type, \"box\")) {\n /* Coordinates are all -0.5 to 0.5 */\n min1[0] = min1[1] = min1[2] = -0.5;\n max1[0] = max1[1] = max1[2] = 0.5;\n } else {\n if (xyzmin != -424242) {\n min1[0] = min1[1] = min1[2] = xyzmin;\n } else {\n min1[0] = xmin;\n min1[1] = ymin;\n min1[2] = zmin;\n }\n\n if (xyzmax != -424242) {\n max1[0] = max1[1] = max1[2] = xyzmax;\n } else {\n max1[0] = xmax;\n max1[1] = ymax;\n max1[2] = zmax;\n }\n }\n\n /* Grid spacing on coarsest grid */\n for (i = 0; i < 3; i++) {\n if (domainsym[2 * i + 0]) {\n if (cntstag[i]) {\n *coarse_delta[i] =\n max1[i] / (cctk_gsh[i] - cctk_nghostzones[i] - 0.5);\n origin[i] = -(cctk_nghostzones[i] - 0.5) * *coarse_delta[i];\n } else {\n *coarse_delta[i] = max1[i] / (cctk_gsh[i] - cctk_nghostzones[i] - 1);\n origin[i] = -cctk_nghostzones[i] * *coarse_delta[i];\n }\n } else if (domainsym[2 * i + 1]) {\n if (cntstag[i]) {\n *coarse_delta[i] =\n fabs(min1[i]) / (cctk_gsh[i] - cctk_nghostzones[i] - 0.5);\n } else {\n *coarse_delta[i] =\n fabs(min1[i]) / (cctk_gsh[i] - cctk_nghostzones[i] - 1);\n }\n origin[i] = min1[i];\n } else {\n if (cntstag[i]) {\n CCTK_VWarn(4, __LINE__, __FILE__, CCTK_THORNSTRING,\n \"Ignoring request to avoid origin in %c-direction, \"\n \"it is not relevant for this grid type\",\n 'x' + i);\n }\n *coarse_delta[i] = (max1[i] - min1[i]) / max(cctk_gsh[i] - 1, 1);\n origin[i] = min1[i];\n }\n\n this_delta[i] = *coarse_delta[i];\n }\n }\n\n /**************************************************************\n * BYSPACING\n *\n * User gives: grid spacing on the coarsest GH and\n * the number of gridpoints on the coarsest GH\n *\n **************************************************************/\n else if (CCTK_Equals(type, \"byspacing\")) {\n /* Dx, Dy, Dx on the coarsest grid */\n if (dxyz > 0) {\n *coarse_delta[0] = *coarse_delta[1] = *coarse_delta[2] = dxyz;\n } else {\n *coarse_delta[0] = dx;\n *coarse_delta[1] = dy;\n *coarse_delta[2] = dz;\n }\n\n for (i = 0; i < 3; i++) {\n this_delta[i] = *coarse_delta[i];\n\n /* Set minimum values of coordinates */\n if (domainsym[2 * i + 0]) {\n origin[i] = -(cctk_nghostzones[i] - cntstag[i] * 0.5);\n } else if (domainsym[2 * i + 1]) {\n origin[i] = -(cctk_gsh[i] - 1 - cctk_nghostzones[i] + cntstag[i] * 0.5);\n } else {\n origin[i] = -0.5 * (cctk_gsh[i] - 1 - cntstag[i] * cctk_gsh[i] % 2);\n }\n origin[i] *= this_delta[i];\n }\n }\n\n /**************************************************************\n * COORDBASE\n *\n * CoordBase gives: grid spacing on the coarsest GH and\n * minimum and maximum values of coordinates and\n * the number of gridpoints on the coarsest GH\n *\n **************************************************************/\n else if (CCTK_Equals(type, \"coordbase\")) {\n CCTK_REAL physical_min[3];\n CCTK_REAL physical_max[3];\n CCTK_REAL interior_min[3];\n CCTK_REAL interior_max[3];\n CCTK_REAL exterior_min[3];\n CCTK_REAL exterior_max[3];\n CCTK_REAL spacing[3];\n int d;\n\n ierr = GetDomainSpecification(3, physical_min, physical_max, interior_min,\n interior_max, exterior_min, exterior_max,\n spacing);\n if (ierr)\n CCTK_WARN(0, \"error returned from function GetDomainSpecification\");\n\n /* Adapt to convergence level */\n for (d = 0; d < 3; ++d) {\n spacing[d] *= pow(cctkGH->cctk_convfac, cctkGH->cctk_convlevel);\n }\n\n ierr = ConvertFromPhysicalBoundary(3, physical_min, physical_max,\n interior_min, interior_max, exterior_min,\n exterior_max, spacing);\n if (ierr)\n CCTK_WARN(0, \"error returned from function ConvertFromPhysicalBoundary\");\n\n for (d = 0; d < 3; ++d) {\n origin[d] = exterior_min[d];\n this_delta[d] = spacing[d];\n *coarse_delta[d] = this_delta[d];\n }\n }\n /**************************************************************\n * MULTIPATCH\n *\n * MultiPatch gives: grid spacing on the coarsest GH and\n * minimum and maximum values of coordinates and\n * the number of gridpoints on the coarsest GH\n *\n **************************************************************/\n else if (CCTK_Equals(type, \"multipatch\")) {\n CCTK_REAL physical_min[3];\n CCTK_REAL physical_max[3];\n CCTK_REAL interior_min[3];\n CCTK_REAL interior_max[3];\n CCTK_REAL exterior_min[3];\n CCTK_REAL exterior_max[3];\n CCTK_REAL spacing[3];\n CCTK_INT map;\n int d;\n\n map = MultiPatch_GetMap(cctkGH);\n if (map < 0)\n CCTK_WARN(0, \"error returned from function MultiPatch_GetMap\");\n\n ierr = MultiPatch_GetDomainSpecification(\n map, 3, physical_min, physical_max, interior_min, interior_max,\n exterior_min, exterior_max, spacing);\n if (ierr)\n CCTK_WARN(\n 0, \"error returned from function MultiPatch_GetDomainSpecification\");\n\n /* Adapt to convergence level */\n for (d = 0; d < 3; ++d) {\n spacing[d] *= pow(cctkGH->cctk_convfac, cctkGH->cctk_convlevel);\n }\n\n if (CCTK_IsFunctionAliased(\"MultiPatch_ConvertFromPhysicalBoundary\")) {\n ierr = MultiPatch_ConvertFromPhysicalBoundary(\n map, 3, physical_min, physical_max, interior_min, interior_max,\n exterior_min, exterior_max, spacing);\n if (ierr)\n CCTK_WARN(0, \"error returned from function \"\n \"MultiPatch_ConvertFromPhysicalBoundary\");\n } else {\n ierr = ConvertFromPhysicalBoundary(3, physical_min, physical_max,\n interior_min, interior_max,\n exterior_min, exterior_max, spacing);\n if (ierr)\n CCTK_WARN(0,\n \"error returned from function ConvertFromPhysicalBoundary\");\n }\n\n for (d = 0; d < 3; ++d) {\n origin[d] = exterior_min[d];\n this_delta[d] = spacing[d];\n *coarse_delta[d] = this_delta[d];\n }\n }\n\n else {\n if (0)\n CCTK_WARN(0, \"type is out of bounds\");\n }\n\n /* Register the coordinate ranges */\n for (i = 0; i < 3; i++) {\n cctkGH->cctk_delta_space[i] = this_delta[i];\n cctkGH->cctk_origin_space[i] = origin[i];\n\n lower[i] = origin[i];\n upper[i] = origin[i] + this_delta[i] * (cctk_gsh[i] - 1);\n if (CCTK_CoordRegisterRange(cctkGH, lower[i], upper[i], i + 1, NULL,\n \"cart3d\") < 0) {\n CCTK_VWarn(0, __LINE__, __FILE__, CCTK_THORNSTRING,\n \"Failed to register %c-coordinate computational range\",\n 'x' + i);\n }\n if (CCTK_CoordRegisterRangePhysIndex(cctkGH, loweri[i], upperi[i], i + 1,\n NULL, \"cart3d\") < 0) {\n CCTK_VWarn(0, __LINE__, __FILE__, CCTK_THORNSTRING,\n \"Failed to register %c-coordinate physical range\", 'x' + i);\n }\n }\n\n CCTK_INFO(\"Grid Spacings:\");\n CCTK_VInfo(CCTK_THORNSTRING, \"dx=>%12.7e dy=>%12.7e dz=>%12.7e\",\n (double)cctk_delta_space[0], (double)cctk_delta_space[1],\n (double)cctk_delta_space[2]);\n CCTK_INFO(\"Computational Coordinates:\");\n CCTK_VInfo(CCTK_THORNSTRING,\n \"x=>[%6.3f,%6.3f] y=>[%6.3f,%6.3f] z=>[%6.3f,%6.3f]\", lower[0],\n upper[0], lower[1], upper[1], lower[2], upper[2]);\n CCTK_INFO(\"Indices of Physical Coordinates:\");\n CCTK_VInfo(CCTK_THORNSTRING, \"x=>[%d,%d] y=>[%d,%d] z=>[%d,%d]\", loweri[0],\n upperi[0], loweri[1], upperi[1], loweri[2], upperi[2]);\n\n if ((domainsym[0] == GFSYM_ROTATION_Y || domainsym[2] == GFSYM_ROTATION_X) &&\n (lower[2] + upper[2] > 1e-12)) {\n CCTK_WARN(0, \"minimum z must equal maximum z for rotation symmetry\");\n }\n\n if ((domainsym[0] == GFSYM_ROTATION_Z || domainsym[4] == GFSYM_ROTATION_X) &&\n (lower[1] + upper[1] > 1e-12)) {\n CCTK_WARN(0, \"minimum y must equal maximum y for rotation symmetry\");\n }\n\n if ((domainsym[2] == GFSYM_ROTATION_Z || domainsym[4] == GFSYM_ROTATION_Y) &&\n (lower[0] - upper[0] > 1e-12)) {\n CCTK_WARN(0, \"minimum x must equal maximum x for rotation symmetry\");\n }\n\n /* cart3d */\n for (i = 0; i < 3; i++) {\n sprintf(coord_name, \"%c\", 'x' + i);\n\n coord_handle = Coord_CoordHandle(cctkGH, coord_name, \"cart3d\");\n if (coord_handle < 0) {\n CCTK_VWarn(0, __LINE__, __FILE__, CCTK_THORNSTRING,\n \"Error retreiving coordinate handle for '%s' of cart3d\",\n coord_name);\n }\n sprintf(coord_name, \"grid::%c\", 'x' + i);\n ierr = Util_TableSetInt(coord_handle, loweri[i], \"PHYSICALMIN\");\n ierr += Util_TableSetReal(coord_handle, lower[i], \"COMPMIN\");\n ierr += Util_TableSetInt(coord_handle, upperi[i], \"PHYSICALMAX\");\n ierr += Util_TableSetReal(coord_handle, upper[i], \"COMPMAX\");\n ierr += Util_TableSetString(coord_handle, \"uniform\", \"TYPE\");\n ierr += Util_TableSetString(coord_handle, \"no\", \"TIMEDEPENDENT\");\n ierr += Util_TableSetString(coord_handle, \"CCTK_REAL\", \"DATATYPE\");\n ierr +=\n Util_TableSetInt(coord_handle, CCTK_VarIndex(coord_name), \"GAINDEX\");\n ierr += Util_TableSetReal(coord_handle, cctk_delta_space[i], \"DELTA\");\n }\n\n /* cart2d */\n for (i = 0; i < 2; i++) {\n sprintf(coord_name, \"%c\", 'x' + i);\n coord_handle = Coord_CoordHandle(cctkGH, coord_name, \"cart2d\");\n if (coord_handle < 0) {\n CCTK_VWarn(0, __LINE__, __FILE__, CCTK_THORNSTRING,\n \"Error retreiving coordinate handle for '%s' of cart2d\",\n coord_name);\n }\n sprintf(coord_name, \"grid::%c\", 'x' + i);\n ierr = Util_TableSetReal(coord_handle, lower[i], \"PHYSICALMIN\"); /*??*/\n ierr += Util_TableSetReal(coord_handle, lower[i], \"COMPMIN\");\n ierr += Util_TableSetReal(coord_handle, upper[i], \"PHYSICALMAX\"); /*??*/\n ierr += Util_TableSetReal(coord_handle, upper[i], \"COMPMAX\");\n ierr += Util_TableSetString(coord_handle, \"uniform\", \"TYPE\");\n ierr += Util_TableSetString(coord_handle, \"no\", \"TIMEDEPENDENT\");\n ierr += Util_TableSetString(coord_handle, \"CCTK_REAL\", \"DATATYPE\");\n ierr +=\n Util_TableSetInt(coord_handle, CCTK_VarIndex(coord_name), \"GAINDEX\");\n ierr += Util_TableSetReal(coord_handle, cctk_delta_space[i], \"DELTA\");\n }\n\n /* cart1d */\n coord_handle = Coord_CoordHandle(cctkGH, \"x\", \"cart1d\");\n if (coord_handle < 0) {\n CCTK_WARN(0, \"Error retreiving coordinate handle for x of cart1d\");\n }\n ierr = Util_TableSetReal(coord_handle, lower[0], \"PHYSICALMIN\"); /*??*/\n ierr += Util_TableSetReal(coord_handle, lower[0], \"COMPMIN\");\n ierr += Util_TableSetReal(coord_handle, upper[0], \"PHYSICALMAX\"); /*??*/\n ierr += Util_TableSetReal(coord_handle, upper[0], \"COMPMAX\");\n ierr += Util_TableSetString(coord_handle, \"uniform\", \"TYPE\");\n ierr += Util_TableSetString(coord_handle, \"no\", \"TIMEDEPENDENT\");\n ierr += Util_TableSetString(coord_handle, \"CCTK_REAL\", \"DATATYPE\");\n ierr += Util_TableSetInt(coord_handle, CCTK_VarIndex(\"grid::x\"), \"GAINDEX\");\n ierr += Util_TableSetReal(coord_handle, cctk_delta_space[0], \"DELTA\");\n\n /* Set up coordinate tables */\n /* Should this be done in a function?\n WriteCoordinateTable(cctkGH, \"cart3d\"); */\n}\n\n/*@@\n @routine CartGrid3D_SetCoordinates\n @date 2004-06-17\n @author Christian Ott\n @desc\n Sets up Cartesian coordinates.\n @enddesc\n @var CCTK_ARGUMENTS\n @vdesc Cactus argument list\n @vtype\n @vio in/out\n @endvar\n@@*/\nvoid CartGrid3D_SetCoordinates(CCTK_ARGUMENTS) {\n DECLARE_CCTK_ARGUMENTS;\n\n /* CCTK_VInfo(CCTK_THORNSTRING,\"Resetting coordinates after regridding.\"); */\n\n const int di = 1;\n const int dj = di * (cctk_ash[0] + 1);\n const int dk = dj * (cctk_ash[1] + 1);\n\n CCTK_LOOP3_ALL(CartGrid3D, cctkGH, i,j,k) {\n int idx = i*di + j*dj + k*dk;\n x[idx] = CCTK_DELTA_SPACE(0) * (i + cctk_lbnd[0]) + CCTK_ORIGIN_SPACE(0);\n y[idx] = CCTK_DELTA_SPACE(1) * (j + cctk_lbnd[1]) + CCTK_ORIGIN_SPACE(1);\n z[idx] = CCTK_DELTA_SPACE(2) * (k + cctk_lbnd[2]) + CCTK_ORIGIN_SPACE(2);\n r[idx] = sqrt(SQR(x[idx]) + SQR(y[idx]) + SQR(z[idx]));\n } CCTK_ENDLOOP3_ALL(CartGrid3D);\n}\n","blame":{"changes":["ANFQ3BZCQVIJLNWXETGDK26W3QLGIXH6V5FHSB5DEXWJY35YTUKAC","Q5LTJSV55Y3VKO7T7G325DWOAP2R4G25IC7ENJBYOGL2Q3MDHYTQC"],"lines":[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,1,0]}}},"login":null}