| 1 | !!****if* source/Grid/GridMain/UG/Grid_init
|
|---|
| 2 | !!
|
|---|
| 3 | !! NAME
|
|---|
| 4 | !!
|
|---|
| 5 | !! Grid_init
|
|---|
| 6 | !!
|
|---|
| 7 | !!
|
|---|
| 8 | !! SYNOPSIS
|
|---|
| 9 | !!
|
|---|
| 10 | !! Grid_init()
|
|---|
| 11 | !!
|
|---|
| 12 | !!
|
|---|
| 13 | !! DESCRIPTION
|
|---|
| 14 | !!
|
|---|
| 15 | !! Initialize the runtime parameters of the Grid unit.
|
|---|
| 16 | !!
|
|---|
| 17 | !!
|
|---|
| 18 | !! ARGUMENTS
|
|---|
| 19 | !!
|
|---|
| 20 | !!
|
|---|
| 21 | !! PARAMETERS
|
|---|
| 22 | !!
|
|---|
| 23 | !! iguard [INTEGER]
|
|---|
| 24 | !! number of guard cells along IAXIS
|
|---|
| 25 | !! jguard [INTEGER]
|
|---|
| 26 | !! number of guard cells along JAXIS
|
|---|
| 27 | !! kguard [INTEGER]
|
|---|
| 28 | !! number of guard cells along KAXIS
|
|---|
| 29 | !! iGridSize [INTEGER]
|
|---|
| 30 | !! the global number of interior grid cells along IAXIS
|
|---|
| 31 | !! only needed when operating in NONFIXEDBLOSIZE mode
|
|---|
| 32 | !! jGridSize [INTEGER]
|
|---|
| 33 | !! the global number of interior grid cells along JAXIS
|
|---|
| 34 | !! only needed when operating in NONFIXEDBLOSIZE mode
|
|---|
| 35 | !! kGridSize [INTEGER]
|
|---|
| 36 | !! the global number of interior grid cells along KAXIS
|
|---|
| 37 | !! only needed when operating in NONFIXEDBLOSIZE mode
|
|---|
| 38 | !!
|
|---|
| 39 | !! flux_correct [BOOLEAN]
|
|---|
| 40 | !! turns flux correction on or off
|
|---|
| 41 | !! should always be false in UG since all blocks are on same level
|
|---|
| 42 | !! compute_grid_size [BOOLEAN]
|
|---|
| 43 | !! compute grid size in the case of non-fixed-block size,
|
|---|
| 44 | !! non fixed block size mode means block dims are not specified at compile time
|
|---|
| 45 | !! smalle [REAL]
|
|---|
| 46 | !! Cutoff value for energy (used in some boundary condition handling)
|
|---|
| 47 | !! smallx [REAL]
|
|---|
| 48 | !! Cutoff value for abundances
|
|---|
| 49 | !!
|
|---|
| 50 | !!***
|
|---|
| 51 |
|
|---|
| 52 | #ifdef DEBUG_ALL
|
|---|
| 53 | #define DEBUG_GRID
|
|---|
| 54 | #endif
|
|---|
| 55 |
|
|---|
| 56 | subroutine Grid_init()
|
|---|
| 57 |
|
|---|
| 58 | use physicalData
|
|---|
| 59 | use Driver_interface, ONLY : Driver_abortFlash, Driver_getMype, &
|
|---|
| 60 | Driver_getNumProcs, Driver_getComm
|
|---|
| 61 | use RuntimeParameters_interface, ONLY : RuntimeParameters_get, &
|
|---|
| 62 | RuntimeParameters_mapStrToInt
|
|---|
| 63 |
|
|---|
| 64 | use Driver_interface, ONLY : Driver_getMype, Driver_getNumProcs, Driver_getComm
|
|---|
| 65 | use Grid_data
|
|---|
| 66 | use gr_sbInterface, ONLY : gr_sbInit
|
|---|
| 67 |
|
|---|
| 68 | implicit none
|
|---|
| 69 |
|
|---|
| 70 | #include "Flash.h"
|
|---|
| 71 | #include "constants.h"
|
|---|
| 72 |
|
|---|
| 73 |
|
|---|
| 74 | character(len=MAX_STRING_LENGTH) :: xl_bcString,xr_bcString
|
|---|
| 75 | character(len=MAX_STRING_LENGTH) :: yl_bcString,yr_bcString
|
|---|
| 76 | character(len=MAX_STRING_LENGTH) :: zl_bcString,zr_bcString
|
|---|
| 77 | character(len=MAX_STRING_LENGTH) :: eosModeString, grav_boundary_type
|
|---|
| 78 |
|
|---|
| 79 | logical :: useProtonImaging
|
|---|
| 80 |
|
|---|
| 81 | integer :: i, localnxb, localnyb, localnzb
|
|---|
| 82 |
|
|---|
| 83 | !! Get all of the parallel environment information from the driver
|
|---|
| 84 |
|
|---|
| 85 | call Driver_getMype(GLOBAL_COMM, gr_globalMe)
|
|---|
| 86 | call Driver_getNumProcs(GLOBAL_COMM, gr_globalNumProcs)
|
|---|
| 87 | call Driver_getComm(GLOBAL_COMM, gr_globalComm)
|
|---|
| 88 |
|
|---|
| 89 | call Driver_getMype(MESH_COMM, gr_meshMe)
|
|---|
| 90 | call Driver_getNumProcs(MESH_COMM, gr_meshNumProcs)
|
|---|
| 91 | call Driver_getComm(MESH_COMM, gr_meshComm)
|
|---|
| 92 |
|
|---|
| 93 | call Driver_getMype(MESH_ACROSS_COMM, gr_meshAcrossMe)
|
|---|
| 94 | call Driver_getNumProcs(MESH_ACROSS_COMM, gr_meshAcrossNumProcs)
|
|---|
| 95 | call Driver_getComm(MESH_ACROSS_COMM, gr_meshAcrossComm)
|
|---|
| 96 |
|
|---|
| 97 | do i = 1, MDIM
|
|---|
| 98 | call Driver_getMype(AXIS_COMM, gr_axisMe(i),i)
|
|---|
| 99 | call Driver_getNumProcs(AXIS_COMM, gr_axisNumProcs(i),i)
|
|---|
| 100 | call Driver_getComm(AXIS_COMM, gr_axisComm(i),i)
|
|---|
| 101 | end do
|
|---|
| 102 |
|
|---|
| 103 | call RuntimeParameters_get("geometry", gr_str_geometry)
|
|---|
| 104 | call RuntimeParameters_mapStrToInt(gr_str_geometry, gr_geometry)
|
|---|
| 105 | call RuntimeParameters_get("geometryOverride",gr_geometryOverride)
|
|---|
| 106 |
|
|---|
| 107 | call RuntimeParameters_get("bndPriorityOne",gr_bndOrder(1))
|
|---|
| 108 | call RuntimeParameters_get("bndPriorityTwo",gr_bndOrder(2))
|
|---|
| 109 | call RuntimeParameters_get("bndPriorityThree",gr_bndOrder(3))
|
|---|
| 110 |
|
|---|
| 111 | call RuntimeParameters_get("iGridSize", gr_gIndexSize(IAXIS))
|
|---|
| 112 | call RuntimeParameters_get("jGridSize", gr_gIndexSize(JAXIS))
|
|---|
| 113 | call RuntimeParameters_get("kGridSize", gr_gIndexSize(KAXIS))
|
|---|
| 114 |
|
|---|
| 115 | call RuntimeParameters_get("compute_grid_size",gr_compute_grid_size)
|
|---|
| 116 |
|
|---|
| 117 |
|
|---|
| 118 | !get the boundary conditions stored as strings in the flash.par file
|
|---|
| 119 | call RuntimeParameters_get("xl_boundary_type", xl_bcString)
|
|---|
| 120 | call RuntimeParameters_get("xr_boundary_type", xr_bcString)
|
|---|
| 121 | call RuntimeParameters_get("yl_boundary_type", yl_bcString)
|
|---|
| 122 | call RuntimeParameters_get("yr_boundary_type", yr_bcString)
|
|---|
| 123 | call RuntimeParameters_get("zl_boundary_type", zl_bcString)
|
|---|
| 124 | call RuntimeParameters_get("zr_boundary_type", zr_bcString)
|
|---|
| 125 |
|
|---|
| 126 | !map the string boundary conditions to integer constants defined in constants.h
|
|---|
| 127 | call RuntimeParameters_mapStrToInt(xl_bcString,gr_domainBC(LOW,IAXIS))
|
|---|
| 128 | call RuntimeParameters_mapStrToInt(xr_bcString,gr_domainBC(HIGH,IAXIS))
|
|---|
| 129 | call RuntimeParameters_mapStrToInt(yl_bcString,gr_domainBC(LOW,JAXIS))
|
|---|
| 130 | call RuntimeParameters_mapStrToInt(yr_bcString,gr_domainBC(HIGH,JAXIS))
|
|---|
| 131 | call RuntimeParameters_mapStrToInt(zl_bcString,gr_domainBC(LOW,KAXIS))
|
|---|
| 132 | call RuntimeParameters_mapStrToInt(zr_bcString,gr_domainBC(HIGH,KAXIS))
|
|---|
| 133 |
|
|---|
| 134 | call RuntimeParameters_get('xmin', gr_imin)
|
|---|
| 135 | call RuntimeParameters_get('xmax', gr_imax)
|
|---|
| 136 | call RuntimeParameters_get('ymin', gr_jmin)
|
|---|
| 137 | call RuntimeParameters_get('ymax', gr_jmax)
|
|---|
| 138 | call RuntimeParameters_get('zmin', gr_kmin)
|
|---|
| 139 | call RuntimeParameters_get('zmax', gr_kmax)
|
|---|
| 140 |
|
|---|
| 141 | gr_globalDomain(LOW,IAXIS) = gr_imin
|
|---|
| 142 | gr_globalDomain(LOW,JAXIS) = gr_jmin
|
|---|
| 143 | gr_globalDomain(LOW,KAXIS) = gr_kmin
|
|---|
| 144 | gr_globalDomain(HIGH,IAXIS) = gr_imax
|
|---|
| 145 | gr_globalDomain(HIGH,JAXIS) = gr_jmax
|
|---|
| 146 | gr_globalDomain(HIGH,KAXIS) = gr_kmax
|
|---|
| 147 |
|
|---|
| 148 | call RuntimeParameters_get('smalle', gr_smalle)
|
|---|
| 149 | call RuntimeParameters_get('smallx', gr_smallx)
|
|---|
| 150 | call RuntimeParameters_get('smlrho', gr_smallrho)
|
|---|
| 151 |
|
|---|
| 152 | ! Determine the geometries of the individual dimensions, and scale
|
|---|
| 153 | ! angle value parameters that are expressed in degrees to radians.
|
|---|
| 154 | ! This call must be made after gr_geometry, gr_domainBC, and gr_{j,k}{min,max}
|
|---|
| 155 | ! have been set based on the corresponding runtime parameters.
|
|---|
| 156 | call gr_initGeometry()
|
|---|
| 157 |
|
|---|
| 158 |
|
|---|
| 159 | #ifdef FLASH_EOS
|
|---|
| 160 | call RuntimeParameters_get("eosMode", eosModeString)
|
|---|
| 161 | call RuntimeParameters_mapStrToInt(eosModeString,gr_eosMode)
|
|---|
| 162 | call RuntimeParameters_get("eosModeInit", eosModeString)
|
|---|
| 163 | call RuntimeParameters_mapStrToInt(eosModeString,gr_eosModeInit)
|
|---|
| 164 | #else
|
|---|
| 165 | gr_eosMode = 1
|
|---|
| 166 | gr_eosModeInit = 1
|
|---|
| 167 | #endif
|
|---|
| 168 | gr_eosModeNow = gr_eosModeInit ! may change after initialization is done
|
|---|
| 169 |
|
|---|
| 170 | #ifdef FLASH_PARTICLES
|
|---|
| 171 | call RuntimeParameters_get('useParticles',gr_useParticles)
|
|---|
| 172 | #else
|
|---|
| 173 | gr_useParticles=.false.
|
|---|
| 174 | #endif
|
|---|
| 175 |
|
|---|
| 176 | #ifdef FLASH_EDEP
|
|---|
| 177 | call RuntimeParameters_get('useEnergyDeposition', gr_useEnergyDeposition)
|
|---|
| 178 | gr_useParticles = gr_useEnergyDeposition
|
|---|
| 179 | #else
|
|---|
| 180 | gr_useEnergyDeposition = .false.
|
|---|
| 181 | #endif
|
|---|
| 182 |
|
|---|
| 183 | #ifdef FLASH_GRID_PARTICLES
|
|---|
| 184 | call RuntimeParameters_get('useProtonImaging',useProtonImaging)
|
|---|
| 185 | if (useProtonImaging) then
|
|---|
| 186 | gr_useParticles=.true.
|
|---|
| 187 | end if
|
|---|
| 188 | #endif
|
|---|
| 189 |
|
|---|
| 190 | gr_globalNumBlocks = gr_meshNumProcs
|
|---|
| 191 | gr_globalOffset = gr_meshMe
|
|---|
| 192 |
|
|---|
| 193 |
|
|---|
| 194 | gr_justExchangedGC = .false.
|
|---|
| 195 | !! do i = UNK_VARS_BEGIN,UNK_VARS_END
|
|---|
| 196 | !! gr_vars(i)=i
|
|---|
| 197 | !! end do
|
|---|
| 198 |
|
|---|
| 199 | #ifdef FIXEDBLOCKSIZE
|
|---|
| 200 | gr_gIndexSize(IAXIS)=gr_axisNumProcs(IAXIS)*NXB
|
|---|
| 201 | gr_gIndexSize(JAXIS)=gr_axisNumProcs(JAXIS)*NYB
|
|---|
| 202 | gr_gIndexSize(KAXIS)=gr_axisNumProcs(KAXIS)*NZB
|
|---|
| 203 | #else
|
|---|
| 204 |
|
|---|
| 205 | if (any((gr_gIndexSize / gr_axisNumProcs) <= 0)) then
|
|---|
| 206 | if(gr_meshMe == MASTER_PE) print*,'the local block size is:', gr_gIndexSize / gr_axisNumProcs
|
|---|
| 207 | call Driver_abortFlash("[Grid_init] Must set runtime parameters iGridSize, jGridSize, kGridSize so that" // &
|
|---|
| 208 | " a block contains at least one cell.")
|
|---|
| 209 | end if
|
|---|
| 210 |
|
|---|
| 211 | #endif
|
|---|
| 212 |
|
|---|
| 213 | call gr_create_surr_blks()
|
|---|
| 214 |
|
|---|
| 215 |
|
|---|
| 216 | gr_guard = NGUARD
|
|---|
| 217 | gr_guard(JAXIS) = gr_guard(JAXIS)*K2D
|
|---|
| 218 | gr_guard(KAXIS) = gr_guard(KAXIS)*K3D
|
|---|
| 219 |
|
|---|
| 220 | gr_allPeriodic = .true.
|
|---|
| 221 | do i = 1,NDIM
|
|---|
| 222 | if(gr_domainBC(LOW,i)/=PERIODIC)gr_allPeriodic=.false.
|
|---|
| 223 | if(gr_domainBC(HIGH,i)/=PERIODIC)gr_allPeriodic=.false.
|
|---|
| 224 | end do
|
|---|
| 225 |
|
|---|
| 226 | !Check if there are gravitational isolated boundary conditions
|
|---|
| 227 | !in order to determine which solvers to intialize.
|
|---|
| 228 | call RuntimeParameters_get("grav_boundary_type", grav_boundary_type)
|
|---|
| 229 | gr_isolatedBoundaries = (grav_boundary_type=="isolated")
|
|---|
| 230 |
|
|---|
| 231 |
|
|---|
| 232 | !! Now find the minimum cell size (gr_minCellSize)
|
|---|
| 233 | !! - only considering non-angle coordinates
|
|---|
| 234 |
|
|---|
| 235 | gr_minCellSizes(IAXIS) = (gr_imax - gr_imin) / gr_gIndexSize(IAXIS)
|
|---|
| 236 | gr_minCellSize = gr_minCellSizes(IAXIS)
|
|---|
| 237 |
|
|---|
| 238 | if (NDIM >= 2) then
|
|---|
| 239 | gr_minCellSizes(JAXIS) = (gr_jmax - gr_jmin) / gr_gIndexSize(JAXIS)
|
|---|
| 240 | if (.not.gr_dirIsAngular(JAXIS)) then
|
|---|
| 241 | gr_minCellSize = min(gr_minCellSize,gr_minCellSizes(JAXIS))
|
|---|
| 242 | end if
|
|---|
| 243 | end if
|
|---|
| 244 |
|
|---|
| 245 | if (NDIM == 3) then
|
|---|
| 246 | gr_minCellSizes(KAXIS) = (gr_kmax - gr_kmin) / gr_gIndexSize(KAXIS)
|
|---|
| 247 | if (.not. gr_dirIsAngular(KAXIS)) then
|
|---|
| 248 | gr_minCellSize = min(gr_minCellSize,gr_minCellSizes(KAXIS))
|
|---|
| 249 | end if
|
|---|
| 250 | end if
|
|---|
| 251 |
|
|---|
| 252 |
|
|---|
| 253 |
|
|---|
| 254 | call gr_setDataStructInfo()
|
|---|
| 255 | gr_offset(:,:)=0
|
|---|
| 256 | gr_offset(FACEX_DATATYPE,IAXIS)=1
|
|---|
| 257 | gr_offset(FACEY_DATATYPE,JAXIS)=1
|
|---|
| 258 | gr_offset(FACEZ_DATATYPE,KAXIS)=1
|
|---|
| 259 | call gr_createDomain()
|
|---|
| 260 | call gr_ptInit()
|
|---|
| 261 | call gr_bcInit()
|
|---|
| 262 | call gr_ptMapInit()
|
|---|
| 263 |
|
|---|
| 264 | call gr_sbInit()
|
|---|
| 265 |
|
|---|
| 266 | gr_region=0.0
|
|---|
| 267 | end subroutine Grid_init
|
|---|
| 268 |
|
|---|
| 269 |
|
|---|
| 270 | subroutine gr_create_surr_blks
|
|---|
| 271 | use Grid_data, ONLY : surr_blks, gr_axisNumProcs, gr_meshMe, gr_domainBC, &
|
|---|
| 272 | gr_meshNumProcs, gr_axisMe
|
|---|
| 273 | implicit none
|
|---|
| 274 | integer :: surr_dist(-1:1,-K2D:K2D,-K3D:K3D)
|
|---|
| 275 | integer :: iMe, jMe, kMe, iNeigh, jNeigh, kNeigh
|
|---|
| 276 | integer :: i, j, k, n , lx, ly, lz
|
|---|
| 277 | logical :: outsideI, outsideJ, outsideK
|
|---|
| 278 |
|
|---|
| 279 | lx = gr_axisNumProcs(1)
|
|---|
| 280 | if (NDIM < 2) then
|
|---|
| 281 | ly = 1
|
|---|
| 282 | else
|
|---|
| 283 | ly = gr_axisNumProcs(2)
|
|---|
| 284 | end if
|
|---|
| 285 | if (NDIM < 3) then
|
|---|
| 286 | lz = 1
|
|---|
| 287 | else
|
|---|
| 288 | lz = gr_axisNumProcs(3)
|
|---|
| 289 | end if
|
|---|
| 290 |
|
|---|
| 291 | do k = -K3D, K3D
|
|---|
| 292 | kMe = gr_axisMe(3)
|
|---|
| 293 | kNeigh = kMe + k
|
|---|
| 294 | outsideK = (kNeigh < 0 .OR. kNeigh .GE. lz)
|
|---|
| 295 | if (outsideK .AND. k==(-1)) then
|
|---|
| 296 | if (gr_domainBC(LOW,KAXIS) == PERIODIC) then
|
|---|
| 297 | kNeigh = kNeigh + lz
|
|---|
| 298 | else
|
|---|
| 299 | kNeigh = gr_domainBC(LOW,KAXIS)
|
|---|
| 300 | end if
|
|---|
| 301 | else if (outsideK .AND. k==1) then
|
|---|
| 302 | if (gr_domainBC(LOW,KAXIS) == PERIODIC) then
|
|---|
| 303 | kNeigh = kNeigh - lz
|
|---|
| 304 | else
|
|---|
| 305 | kNeigh = gr_domainBC(HIGH,KAXIS)
|
|---|
| 306 | end if
|
|---|
| 307 | end if
|
|---|
| 308 | do j = -K2D, K2D
|
|---|
| 309 | jMe = gr_axisMe(2)
|
|---|
| 310 | jNeigh = jMe + j
|
|---|
| 311 | outsideJ = (jNeigh < 0 .OR. jNeigh .GE. ly)
|
|---|
| 312 | if (outsideJ .AND. j==(-1)) then
|
|---|
| 313 | if (gr_domainBC(LOW,JAXIS) == PERIODIC) then
|
|---|
| 314 | jNeigh = jNeigh + ly
|
|---|
| 315 | else
|
|---|
| 316 | jNeigh = gr_domainBC(LOW,JAXIS)
|
|---|
| 317 | end if
|
|---|
| 318 | else if (outsideJ .AND. j==1) then
|
|---|
| 319 | if (gr_domainBC(LOW,JAXIS) == PERIODIC) then
|
|---|
| 320 | jNeigh = jNeigh - ly
|
|---|
| 321 | else
|
|---|
| 322 | jNeigh = gr_domainBC(HIGH,JAXIS)
|
|---|
| 323 | end if
|
|---|
| 324 | end if
|
|---|
| 325 | do i = -1, 1
|
|---|
| 326 | iMe = gr_axisMe(1)
|
|---|
| 327 | iNeigh = iMe + i
|
|---|
| 328 | outsideI = (iNeigh < 0 .OR. iNeigh .GE. ly)
|
|---|
| 329 | if (outsideI .AND. i==(-1)) then
|
|---|
| 330 | if (gr_domainBC(LOW,IAXIS) == PERIODIC) then
|
|---|
| 331 | iNeigh = iNeigh + lx
|
|---|
| 332 | else
|
|---|
| 333 | iNeigh = gr_domainBC(LOW,IAXIS)
|
|---|
| 334 | end if
|
|---|
| 335 | else if (outsideI .AND. i==1) then
|
|---|
| 336 | if (gr_domainBC(LOW,IAXIS) == PERIODIC) then
|
|---|
| 337 | iNeigh = iNeigh - lx
|
|---|
| 338 | else
|
|---|
| 339 | iNeigh = gr_domainBC(HIGH,IAXIS)
|
|---|
| 340 | end if
|
|---|
| 341 | end if
|
|---|
| 342 | if (iNeigh < 0) then
|
|---|
| 343 | n = iNeigh
|
|---|
| 344 | else if (jNeigh < 0) then
|
|---|
| 345 | n = jNeigh
|
|---|
| 346 | else if (kNeigh < 0) then
|
|---|
| 347 | n = kNeigh
|
|---|
| 348 | else
|
|---|
| 349 | n = iNeigh + (jNeigh + kNeigh * ly) * lx
|
|---|
| 350 | end if
|
|---|
| 351 | surr_dist(i,j,k) = n
|
|---|
| 352 | end do
|
|---|
| 353 | end do
|
|---|
| 354 | end do
|
|---|
| 355 |
|
|---|
| 356 | surr_blks(1,:,:,:,1) = 1
|
|---|
| 357 | surr_blks(2,:,:,:,1) = surr_dist(:,:,:)
|
|---|
| 358 | surr_blks(3,:,:,:,1) = 1
|
|---|
| 359 |
|
|---|
| 360 | #if NDIM == 1
|
|---|
| 361 | #define THREE_POW 3
|
|---|
| 362 | #elif NDIM == 2
|
|---|
| 363 | #define THREE_POW 9
|
|---|
| 364 | #else
|
|---|
| 365 | #define THREE_POW 27
|
|---|
| 366 | #endif
|
|---|
| 367 | 999 format("Neighbors on PE",I5,", lb=",I2,",surr_blks=", &
|
|---|
| 368 | THREE_POW("(",I4,",",I4,",",I2,")"))
|
|---|
| 369 |
|
|---|
| 370 | !print 999,gr_meshMe,1,surr_blks(1:3,:,1:1+2*K2D,1:1+2*K3D,1)
|
|---|
| 371 | end subroutine gr_create_surr_blks
|
|---|