Scrambler
1
|
00001 !######################################################################### 00002 ! 00003 ! Copyright (C) 2003-2012 Department of Physics and Astronomy, 00004 ! University of Rochester, 00005 ! Rochester, NY 00006 ! 00007 ! amr_control.f90 is part of AstroBEAR. 00008 ! 00009 ! AstroBEAR is free software: you can redistribute it and/or modify 00010 ! it under the terms of the GNU General Public License as published by 00011 ! the Free Software Foundation, either version 3 of the License, or 00012 ! (at your option) any later version. 00013 ! 00014 ! AstroBEAR is distributed in the hope that it will be useful, 00015 ! but WITHOUT ANY WARRANTY; without even the implied warranty of 00016 ! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00017 ! GNU General Public License for more details. 00018 ! 00019 ! You should have received a copy of the GNU General Public License 00020 ! along with AstroBEAR. If not, see <http://www.gnu.org/licenses/>. 00021 ! 00022 !######################################################################### 00025 00028 00031 00035 00038 MODULE AmrControl 00039 USE TreeLevelOps, ONLY : InheritNeighborsChildren, InheritOldNodeOverlapsChildren, InheritNewNodeOverlapsChildren, & 00040 InheritOverlapsOldChildren, InheritOverlapsNewChildren, CreateChildrens, UpdateOverlaps, NullifyNeighbors, & 00041 AgeNodesChildren, AgeNodes, DestroyNodes, DestroyOldNodes 00042 00043 USE DataLevelOps, ONLY : ProlongateParentsData, ApplyOverlaps, ApplyChildrenData, ApplyInitialChildrenData, SyncFluxes, & 00044 InitInfos, InitialInitInfos, AfterOverlaps, UpdateChildMasks, SetErrFlags, RestrictionFixups, & 00045 AccumulateFluxes, CoarsenDataForParents, CoarsenInitialDataForParents, InitGrids, AdvanceGrids, & 00046 ApplyPhysicalBCs, AfterFixups, ChildMaskOverlaps, ScheduledAdvanceGrids, CompleteAdvanceGrids, WaitingAdvances, InitRestartGrids, UpdateTimeDerivs, ClearFixupFluxes, UpdateMeanDensity, ClearParentFixups, ClearChildFixups, TransferTimeDerivs 00047 00048 USE DataDeclarations 00049 USE DistributionControl 00050 USE ModuleControl, ONLY : BeforeGlobalStep 00051 USE TimeStep 00052 USE TreeDeclarations 00053 USE TreeNodeOps 00054 USE CommunicationControl 00055 USE ExplicitControl 00056 USE IOControl 00057 USE Timing 00058 USE ParticleControl 00059 USE ProcessingControl 00060 USE ProcessingDeclarations 00061 # if defined HYPRE 00062 USE EllipticControl, ONLY : Elliptic, InitialElliptic, EllipticInit, ApplyEllipticBC 00063 # endif 00064 # if defined PTHREADS 00065 USE ThreadControl 00066 # endif 00067 # if defined PTH 00068 USE PthControl 00069 # endif 00070 00071 IMPLICIT NONE 00072 PRIVATE 00073 PUBLIC AMRInit, AMRAdvance, AMRStart, DomainInit, ClearAllNodeLists, DestroyAllnodes, LevelsInit, AMR, BackupData, RestoreData, BoundaryZoneInit, PrintAllocations 00074 INTEGER :: framecounter=0 00075 REAL(KIND=qPREC) :: temp 00076 SAVE 00077 PUBLIC MpiTest 00078 00079 00080 00081 CONTAINS 00084 00086 SUBROUTINE AmrInit 00087 INTEGER :: i 00088 CALL ReadInGlobalData() 00089 ! ALLOCATE(leveldX(-2:MaxLevel), leveldt(0:MaxLevel), leveltnow(0:MaxLevel), levelgmbc(-2:MaxLevel,2), levelMX(-2:MaxLevel,3)) 00090 ALLOCATE(Nodes(-2:MaxLevel), OldNodes(-2:MaxLevel), ExternalNodes(-2:MaxLevel), & 00091 OldExternalNodes(-2:MaxLevel), LastLocalNode(-2:MaxLevel), LastExternalNode(-2:MaxLevel), & 00092 LastOldLocalNode(-2:MaxLevel), LastOldExternalNode(-2:MaxLevel), BackupNodes(-2:MaxLevel), BackupExternalNodes(-2:MaxLevel)) 00093 DO i=-2,MaxLevel 00094 NULLIFY(Nodes(i)%p) 00095 NULLIFY(OldNodes(i)%p) 00096 NULLIFY(ExternalNodes(i)%p) 00097 NULLIFY(OldExternalNodes(i)%p) 00098 NULLIFY(LastLocalNode(i)%p) 00099 NULLIFY(LastExternalNode(i)%p) 00100 NULLIFY(LastOldLocalNode(i)%p) 00101 NULLIFY(LastOldExternalNode(i)%p) 00102 NULLIFY(BackupNodes(i)%p) 00103 NULLIFY(BackupExternalNodes(i)%p) 00104 END DO 00105 GmGlobal(:,1)=1 00106 GmGlobal(:,2)=GmX 00107 IF (ANY(ANY(Gmthbc(1:nDim,:)==2,2) .AND. .NOT. ALL(Gmthbc(1:nDim,:)==2,2))) THEN 00108 PRINT*, 'Bad Periodic Boundary conditions. If one side is periodic, the other side must be as well. Stopping...' 00109 STOP 00110 END IF 00111 lHydroPeriodic(nDim+1:)=.false. 00112 lHydroPeriodic(1:nDim)=(Gmthbc(1:nDim,1)==2) 00113 lAnyPeriodic=lAnyPeriodic .OR. lHydroPeriodic 00114 # if defined PTH 00115 CALL PthInit() 00116 # else 00117 IF (iThreaded == THREADED) THEN 00118 IF (MPI_ID == 0) PRINT*, 'Error - You must compile with PTHREADFLAG=2 to use iThreaded == ', THREADED 00119 STOP 00120 END IF 00121 # endif 00122 lRestart=lRestart .OR. lPostProcess 00123 END SUBROUTINE AmrInit 00124 00126 SUBROUTINE LevelsInit 00127 USE HyperbolicDeclarations 00128 INTEGER :: i, step 00129 levels(-2:MaxLevel)%id=(/(i,i=-2,MaxLevel)/) 00130 levels(-2)%dx=(GxBounds(1,2)-GxBounds(1,1))/GmX(1) 00131 levels(-2)%mx=Gmx 00132 levels(-2)%steps=1 00133 levels(-2)%step=1 00134 DO i=-1,MaxLevel 00135 levels(i)%dx = levels(i-1)%dx / levels(i-1)%CoarsenRatio 00136 levels(i)%mx(1:nDim)=levels(i-1)%mX(1:nDim)*levels(i-1)%CoarsenRatio 00137 levels(i)%mx(nDim+1:3)=1 00138 levels(i)%steps=levels(i-1)%CoarsenRatio 00139 levels(i)%step=1 00140 END DO 00141 CALL SinkParticleInit 00142 00143 ! Now calculate mbc's staring with the finest level and working backwards. 00144 ! (lost) (needed) 00145 00146 ! levels(MaxLevel)%ombc(2)=max(levels(MaxLevel-1)%steps*hyperbolic_mbc+source_mbc,particle_mbc) !ghost cells in overlap 2 00147 00148 # if defined PTHREADS 00149 CALL ThreadsInit() 00150 # endif 00151 # if defined HYPRE 00152 CALL EllipticInit() 00153 # else 00154 IF(lElliptic)THEN 00155 IF(MPI_ID==0) PRINT*, "Elliptic solver requested, but AstroBEAR is not compiled with hypre. Stopping... " 00156 STOP 00157 END IF 00158 # endif 00159 CALL ExplicitInit() 00160 END SUBROUTINE LevelsInit 00161 00162 00164 SUBROUTINE BoundaryZoneInit() 00165 INTEGER :: i, step 00166 00167 IF (lPostProcess) THEN 00168 hyperbolic_mbc=0d0 00169 particle_mbc=0d0 00170 elliptic_mbc=0d0 00171 afterstep_mbc=0d0 00172 END IF 00173 levels(MaxLevel)%ambc(1)=(levels(MaxLevel)%steps-1)*(hyperbolic_mbc+afterstep_mbc)+afterstep_mbc !Size of ghost zone region to update in first advance 00174 levels(MaxLevel)%ambc(2)=afterstep_mbc 00175 00176 levels(MaxLevel)%egmbc(1)=max(levels(MaxLevel)%ambc(1)+hyperbolic_mbc+elliptic_mbc, particle_mbc, Processing_mbc) ! Initial elliptic prolongate and ghost 00177 levels(MaxLevel)%gmbc(1)=max(levels(MaxLevel)%ambc(1)+hyperbolic_mbc,particle_mbc, processing_mbc) ! Initial q prolongate and ghost 00178 00179 levels(MaxLevel)%egmbc(2)=max(levels(MaxLevel)%ambc(2)+hyperbolic_mbc+elliptic_mbc, particle_mbc, processing_mbc) ! Second elliptic ghost 00180 levels(MaxLevel)%gmbc(2)=max(levels(MaxLevel)%ambc(2)+hyperbolic_mbc, particle_mbc, processing_mbc) ! Second q ghost 00181 00182 00183 levels(MaxLevel)%ombc(1)=levels(MaxLevel)%egmbc(1) !Initial Overlap range 00184 levels(MaxLevel)%nmbc=levels(MaxLevel)%egmbc(2) 00185 levels(MaxLevel)%ombc(2)=levels(MaxLevel)%nmbc !Overlaps become neighbors on second step 00186 00187 levels(MaxLevel)%pmbc=levels(MaxLevel)%ombc(1)+modulo(levels(MaxLevel)%ombc(1),levels(MaxLevel-1)%CoarsenRatio) !ensures it is divisible by the coarsen ratio 00188 00189 DO i=MaxLevel-1,-2,-1 00190 step=levels(i)%steps 00191 00192 IF (i > -2) levels(i)%pmbc=ceiling(real(levels(i+1)%ombc(1))/real(levels(i)%CoarsenRatio)) 00193 00194 levels(i)%ambc(step)=afterstep_mbc 00195 levels(i)%egmbc(step)=max(levels(i)%ambc(step)+hyperbolic_mbc+elliptic_mbc, ceiling(real(particle_mbc)/real(product(Levels(i:MaxLevel-1)%CoarsenRatio))), levels(i)%pmbc, Processing_mbc) ! Second elliptic ghost 00196 levels(i)%gmbc(step)=max(levels(i)%ambc(step)+hyperbolic_mbc, ceiling(real(particle_mbc)/real(product(Levels(i:MaxLevel-1)%CoarsenRatio))), levels(i)%pmbc, Processing_mbc) ! Second q ghost 00197 levels(i)%ombc(step)=max(levels(i)%egmbc(step), levels(i)%pmbc) !Initial Neighbor range. Ensures that no missed nieghbors children could be within ombc(1) on the next level 00198 levels(i)%nmbc=levels(i)%ombc(step) !Overlaps become neighbors on second step 00199 00200 IF (levels(i)%steps == 1) THEN 00201 levels(i)%ambc(2)=0 00202 levels(i)%egmbc(2)=0 00203 levels(i)%gmbc(2)=0 00204 levels(i)%ombc(2)=0 00205 ELSE 00206 DO step=levels(i)%steps-1, 1, -1 00207 levels(i)%ambc(step)=levels(i)%ombc(step+1)+afterstep_mbc 00208 levels(i)%egmbc(step)=max(levels(i)%ambc(step)+hyperbolic_mbc+elliptic_mbc, ceiling(real(particle_mbc)/real(product(Levels(i:MaxLevel-1)%CoarsenRatio))), Processing_mbc) ! Initial elliptic prolongate and ghost 00209 levels(i)%gmbc(step)=max(levels(i)%ambc(step)+hyperbolic_mbc,ceiling(real(particle_mbc)/real(product(Levels(i:MaxLevel-1)%CoarsenRatio))), Processing_mbc) ! Initial q prolongate and ghost 00210 levels(i)%ombc(step)=max(levels(i)%egmbc(step), ceiling(real(levels(i+1)%ombc(1))/real(levels(i)%CoarsenRatio))) !Initial Overlap range 00211 END DO 00212 END IF 00213 END DO 00214 00215 ! IF (MPI_ID == 0) THEN 00216 ! DO i=-2,MaxLevel 00217 ! write(*,'(A8,I2,A8,2I2,A8,2I2,A8,2I2,A8,2I2,A8,1I2,A8,1I2)') "level", i, "ambc=", levels(i)%ambc, "gmbc=",levels(i)%gmbc(1:2),"egmbc=",levels(i)%egmbc,"ombc=",levels(i)%ombc, "nmbc=", levels(i)%nmbc, "pmbc=",levels(i)%pmbc 00218 ! END DO 00219 ! END IF 00220 00221 nPeriodic_overlaps(1:nDim)=levels(0)%ombc(1)/Gmx(1:nDim) + 1 00222 ! write(*,*) 'nPeriodic_overlaps = ', nPeriodic_overlaps, levels(0)%ombc(1) 00223 END SUBROUTINE BoundaryZoneInit 00224 00226 SUBROUTINE DomainInit 00227 INTEGER :: i, iErr 00228 TYPE(NodeBox), POINTER :: root_box 00229 TYPE(NodeDef), POINTER :: root_node 00230 TYPE(DomainDef) :: Domain 00231 NAMELIST /DomainData/ Domain 00232 00233 00234 NULLIFY(root_box) 00235 NULLIFY(root_node) 00236 00237 IF (MPI_ID == 0) THEN 00238 00239 CALL CreateNodeBox(GmGlobal, root_box, 0) 00240 00241 CALL AddNode(-2, root_box, root_node) 00242 !CALL InitialInitInfo(root_node%Info,-2, root_box%mGlobal) 00243 END IF 00244 00245 ALLOCATE(Domains(nDomains)) 00246 Domains(1)%mthbc=Gmthbc 00247 Domains(1)%mGlobal=1 00248 Domains(1)%mGlobal(1:nDim,2)=Gmx(1:nDim) 00249 00250 ! DO i=1,nDomains 00251 ! ! Read in level data values. 00252 ! READ(GLOBAL_DATA_HANDLE,NML=DomainData,IOStat=iErr) 00253 ! IF(iErr/=0) THEN 00254 ! PRINT *, "AmrInit() error: unable to read DomainData namelist from ", GLOBAL_DATA_FILE, "." 00255 ! STOP 00256 ! END IF 00257 ! Domains(i)=Domain 00258 ! END DO 00259 ! ALLOCATE(root_node%proclist(MPI_NP), root_node%proctime(MPI_NP)) 00260 ! root_node%proclist=(/(i-1,i=1,MPI_NP)/) 00261 ! root_node%proctime=1d0 00262 00263 END SUBROUTINE DomainInit 00264 00266 SUBROUTINE ReadInGlobalData() 00267 INTEGER :: iErr,i 00268 INTEGER, DIMENSION(0:MaxDepth) :: InterpOpt 00269 INTEGER, DIMENSION(0:MaxDepth) :: CoarsenRatio=2 00270 REAL(KIND=qPREC), DIMENSION(0:MaxDepth) :: qTolerance 00271 REAL(KIND=qPrec), DIMENSION(0:MaxDepth) :: DesiredFillRatios 00272 LOGICAL :: temp 00273 NAMELIST /LevelData/ qTolerance, DesiredFillRatios 00274 00275 MaintainAuxArrays=.true. 00276 ! Open data file (quit if not successful). 00277 OPEN(UNIT=GLOBAL_DATA_HANDLE,FILE=GLOBAL_DATA_FILE) 00278 00279 ! Read in global data values. 00280 !READ(GLOBAL_DATA_HANDLE,NML=GlobalData,IOStat=iErr) 00281 READ(GLOBAL_DATA_HANDLE,NML=GlobalData) 00282 00283 ! If regridding then also restarting 00284 lRestart=lRestart .OR. lRegrid 00285 00286 IF (MaintainAuxArrays) THEN 00287 MaintainAuxArrays = .false. 00288 CLOSE(UNIT=GLOBAL_DATA_HANDLE) 00289 00290 OPEN(UNIT=GLOBAL_DATA_HANDLE,FILE=GLOBAL_DATA_FILE) 00291 00292 ! Read in global data values. 00293 !READ(GLOBAL_DATA_HANDLE,NML=GlobalData,IOStat=iErr) 00294 READ(GLOBAL_DATA_HANDLE,NML=GlobalData) 00295 00296 IF (MaintainAuxArrays) THEN 00297 IF (MPI_ID == 0) THEN 00298 write(*,*) 'MaintainAuxArrays is no longer set in global.data and should be removed' 00299 write(*,*) 'It is instead automatically turned on if lMHD is true in physics.data' 00300 END IF 00301 END IF 00302 ELSE 00303 IF (MPI_ID == 0) THEN 00304 write(*,*) 'MaintainAuxArrays is no longer set in global.data and should be removed' 00305 write(*,*) 'It is instead automatically turned on if lMHD is true in physics.data' 00306 END IF 00307 END IF 00308 00309 ! Read in level data values. 00310 ALLOCATE(levels(-2:MaxLevel)) 00311 InterpOpt=-1 00312 qTolerance=-1 00313 DesiredFillRatios=-1 00314 READ(GLOBAL_DATA_HANDLE,NML=LevelData,IOStat=iErr) 00315 00316 IF(iErr/=0) THEN 00317 PRINT *, "ReadInGlobalData() error: unable to read LevelData namelist from ", GLOBAL_DATA_FILE, "." 00318 STOP 00319 END IF 00320 IF (ANY(InterpOpt > -1)) THEN 00321 IF (MPI_ID == 0) THEN 00322 write(*,*) 'InterpOpt in global.data is now depracated' 00323 write(*,*) 'Use InterpOpts in physics.data instead' 00324 END IF 00325 END IF 00326 DO i=1,MaxLevel 00327 IF (qTolerance(i) == -1) qTolerance(i)=qTolerance(i-1) 00328 IF (DesiredFillRatios(i) == -1) DesiredFillRatios(i)=DesiredFillRatios(i-1) 00329 END DO 00330 00331 00332 IF (iThreaded == NON_THREADED) LevelBalance=(/0d0,0d0/) 00333 DO i=-2,MIN(ROOT_LEVEL-1, MaxLevel) 00334 levels(i)%CoarsenRatio=1 00335 levels(i)%qTolerance=0 00336 levels(i)%DesiredFillRatios=1 00337 levels(i)%steps=1 00338 levels(i)%step=1 00339 END DO 00340 DO i=ROOT_LEVEL,MaxLevel 00341 levels(i)%CoarsenRatio=CoarsenRatio(i) 00342 levels(i)%qTolerance=qTolerance(i) 00343 levels(i)%DesiredFillRatios=DesiredFillRatios(i) 00344 levels(i)%steps=levels(i-1)%CoarsenRatio 00345 levels(i)%step=1 00346 END DO 00347 END SUBROUTINE ReadInGlobalData 00348 00350 00353 00357 SUBROUTINE AMRAdvance(tnext) 00358 REAL(KIND=qPREC) :: tnext,tback 00359 LOGICAL :: lasttimestep, lSuccess 00360 lasttimestep=.false. 00361 DO WHILE (.NOT. lasttimestep) !levels(ROOT_LEVEL)%tnow < tnext) 00362 CALL BackupData() 00363 tback=levels(ROOT_LEVEL)%tnow 00364 lSuccess=.false. 00365 DO WHILE (.not. lSuccess) 00366 CALL GetNextTimeStep(tnext-levels(ROOT_LEVEL)%tnow, lasttimestep) 00367 CALL AMR(BaseLevel) 00368 CALL PrintAllocations 00369 IF (lPrintDebugFrame) THEN 00370 CALL ProcessData() 00371 CALL WriteDataFrame(1000+current_frame) 00372 lPrintDebugFrame=.false. 00373 END IF 00374 CALL TestBadCFL(lSuccess) 00375 IF (.NOT. lSuccess) THEN 00376 levels(ROOT_LEVEL:MaxLevel)%tnow=tback 00377 CALL RestoreData() 00378 lasttimestep=.false. 00379 levels(ROOT_LEVEL)%CurrentLevelStep=levels(ROOT_LEVEL)%CurrentLevelStep-1 00380 END IF 00381 END DO 00382 END DO 00383 END SUBROUTINE AMRAdvance 00384 00386 00387 00388 00391 00394 RECURSIVE SUBROUTINE AMR(n,ladvance_opt) 00395 USE TreeLevelComms 00396 USE DataLevelComms 00397 INTEGER :: n, nSteps, step 00398 INTEGER :: iErr 00399 LOGICAL, OPTIONAL :: ladvance_opt 00400 LOGICAL :: lHasData, lHasChildren, lHasParents, lHasParentswData, lHasChildrenwData, lOverlapsAreOld 00401 LOGICAL :: lAdvance 00402 CALL StartTimer(iAMR, n) 00403 lHasData = (n >= 0) 00404 lHasParentswData = (n >= 1) 00405 lHasChildren = (n < MaxLevel) 00406 lHasParents = lRegridLevel(n) !(n > BaseLevel) 00407 nsteps=levels(n)%steps 00408 lAdvance=.true. 00409 IF (PRESENT(lAdvance_opt)) lAdvance=lAdvance_opt 00410 IF (.NOT. lAdvance) nSteps=1 00411 00412 ! ========================================================! 00413 ! Stage 1 ! 00414 ! Complete communication with parents and get level load ! 00415 ! ========================================================! 00416 IF (lRegridLevel(n)) THEN 00417 CALL CompRecvGridsFromParents(n) 00418 CALL CompSendGridsToChildren(n-1) 00419 CALL PostRecvOverlapsNeighbors(n) 00420 IF (lHasParentswData) CALL PostRecvParentsData(n) 00421 CALL SortNodes(n) 00422 00423 IF (lHasData) CALL InitInfos(n) 00424 00425 CALL CompRecvOverlapsNeighbors(n) 00426 CALL CompSendOverlapsNeighbors(n-1) 00427 00428 IF (lHasParentswData) THEN 00429 CALL CompRecvParentsData(n) 00430 CALL CompSendChildrenData(n-1) 00431 CALL ProlongateParentsData(n) 00432 END IF 00433 00434 IF (lHasData) CALL ChildMaskOverlaps(n) !need to know about neighbors for elliptic calls 00435 ! CALL CompRecvOldNodeOverlaps(n) 00436 ELSE 00437 IF (lHasParentswData) THEN 00438 CALL PostRecvParentsData(n) 00439 CALL CompRecvParentsData(n) 00440 CALL CompSendChildrenData(n-1) 00441 CALL ProlongateParentsData(n) 00442 END IF 00443 END IF 00444 IF (lHasData) CALL GetLevelLoad(n) 00445 IF (lHasParentswData) CALL ClearParentFixups(n) 00446 00447 ! ========================================================! 00448 ! Stage 2 ! 00449 ! Begin AMR steps ! 00450 ! ========================================================! 00451 00452 DO step=1,nsteps 00453 IF (lHasData) CALL ClearFixupFluxes(n) 00454 CALL UpdateLevelStep(n, step, lAdvance) 00455 lOverlapsAreOld = (step == 1 .AND. lRegridLevel(n)) ! On first step of non-persistent level, overlaps come from old grids 00456 00457 ! ========================================================! 00458 ! Stage 2a ! 00459 ! Fill Ghost zones with overlap data ! 00460 ! ========================================================! 00461 00462 IF (.NOT. lRegridLevel(n) .OR. step == 2) CALL UpdateOverlaps(n) ! On second and subsequent steps - ghost zones are filled from neighbors - so overlaps=>neighbors 00463 00464 IF (lHasData) THEN 00465 CALL PostRecvOverlaps(n) 00466 CALL BeforeGlobalStep(n) 00467 CALL PostSendOverlaps(n) 00468 CALL ClearCompletedWorkLoads(n) 00469 CALL ApplyOverlaps(n,step) 00470 CALL CompSendOverlaps(n) 00471 CALL CompRecvOverlaps(n) 00472 IF (lHasParentswData .AND. step == 1) CALL AfterOverlaps(n) 00473 CALL ParticlePreUpdate(n) 00474 CALL ApplyPhysicalBCs(n) 00475 # if defined HYPRE 00476 IF (lElliptic) CALL ApplyEllipticBC(n) 00477 IF (lElliptic .AND. (.NOT. lAdvance .OR. lRegridLevel(n))) THEN !Don't need to solve again if we have not regridded. 00478 IF (lNeedMeanDensity) CALL UpdateMeanDensity(n) 00479 IF (n == 0 .AND. lRegridLevel(0)) CALL TransferTimeDerivs(0) 00480 CALL InitialElliptic(n) 00481 END IF 00482 # endif 00483 IF(lExplicit) CALL InitialExplicit(n) 00484 END IF 00485 00486 ! ========================================================! 00487 ! Stage 2b ! 00488 ! Age children and create new level n+1 nodes ! 00489 ! ========================================================! 00490 00491 IF (lHasChildren) THEN 00492 IF (lRegridLevel(n+1)) THEN 00493 IF (lHasData) CALL SetErrFlags(n) 00494 CALL AgeNodesChildren(n) 00495 CALL AgeNodes(n+1) 00496 CALL CreateChildrens(n) 00497 #if defined PTH 00498 IF (iThreaded == THREADED) CALL pth_MPI_BARRIER2(MPI_COMM_WORLD, ierr) 00499 #endif 00500 CALL DistributeChildrens(n) 00501 00502 00503 ! ========================================================! 00504 ! Stage 2c ! 00505 ! Communicate child information with neighbors ! 00506 ! If iThreaded == NON_THREADED then we can do level -1 ! 00507 ! advances to give everyone a chance to catch up ! 00508 ! Also post send/recv of child information to child procs ! 00509 ! ========================================================! 00510 00511 CALL PostRecvNeighboringChildren(n) 00512 CALL PostSendNeighboringChildren(n) 00513 IF (lOverlapsAreOld) THEN 00514 CALL PostRecvOverlappingChildrenFromOldNodes(n) 00515 CALL PostRecvOverlappingChildrenFromNewNodes(n) 00516 CALL PostSendOverlappingChildrenToOldNodes(n) 00517 CALL PostSendOverlappingChildrenToNewNodes(n) 00518 CALL InheritOldNodeOverlapsChildren(n) 00519 CALL InheritNewNodeOverlapsChildren(n) 00520 CALL InheritNeighborsChildren(n) 00521 ELSE 00522 CALL InheritOverlapsOldChildren(n) 00523 CALL InheritNeighborsChildren(n) 00524 END IF 00525 END IF 00526 IF (lHasData) CALL ClearChildFixups(n) 00527 ! If we need to pre-advance coarse grids then do so now. 00528 IF (iThreaded == NON_THREADED .AND. lHasData .AND. lAdvance) THEN 00529 !advance and call elliptic and update qchild(iPhiDot) before sending children data 00530 CALL AdvanceGrids(n) 00531 #if defined HYPRE 00532 IF (lElliptic) THEN 00533 IF (lNeedMeanDensity) CALL UpdateMeanDensity(n) 00534 CALL Elliptic(n) 00535 END IF 00536 #endif 00537 IF(lExplicit) THEN 00538 CALL Explicit(n) 00539 END IF 00540 CALL ParticlePostUpdate(n) 00541 CALL PrintAdvance(n) 00542 CALL UpdateTimeDerivs(n) 00543 END IF 00544 00545 IF (lRegridLevel(n+1)) THEN 00546 00547 IF (lOverlapsAreOld) THEN 00548 CALL CompRecvOverlappingChildrenFromOldNodes(n) 00549 CALL CompRecvOverlappingChildrenFromNewNodes(n) 00550 CALL PostSendOverlapsToOldNodesChildren(n) 00551 CALL PostRecvOldNodeOverlaps(n+1) 00552 CALL CompRecvNeighboringChildren(n) 00553 CALL CompSendOverlappingChildrenToNewNodes(n) 00554 CALL CompSendOverlappingChildrenToOldNodes(n) 00555 CALL CompSendOverlapsToOldNodesChildren(n) 00556 CALL CompRecvOldNodeOverlaps(n+1) 00557 CALL CompSendNeighboringChildren(n) 00558 ELSE 00559 CALL CompRecvNeighboringChildren(n) 00560 CALL InheritOverlapsNewChildren(n) 00561 CALL PostRecvOldNodeOverlaps(n+1) 00562 CALL PostSendOverlapsToNodesOldChildren(n) 00563 CALL CompSendNeighboringChildren(n) 00564 CALL CompSendOverlapsToNodesOldChildren(n) 00565 CALL CompRecvOldNodeOverlaps(n+1) 00566 END IF 00567 00568 CALL PostSendGridsToChildren(n) 00569 CALL PostRecvGridsFromParents(n+1) 00570 00571 CALL PostSendOverlapsNeighbors(n) 00572 END IF 00573 00574 IF (lHasData) CALL PostSendChildrenData(n) 00575 IF (lAdvance .AND. lHasData) CALL PostRecvChildrenData(n) 00576 00577 IF (.NOT. lAdvance) THEN 00578 CALL ParticlePostUpdate(n) 00579 CALL AMR(n+1, lAdvance_opt=.false.) 00580 END IF 00581 END IF 00582 00583 00584 IF (lAdvance) THEN 00585 00586 ! ========================================================! 00587 ! Stage 2d ! 00588 ! CALL AMR and/or Advance Grids depending on threading ! 00589 ! choice ! 00590 ! ========================================================! 00591 00592 IF (iThreaded == NON_THREADED) THEN 00593 IF (lHasChildren) THEN 00594 CALL AMR(n+1) 00595 ELSE 00596 CALL AdvanceGrids(n) 00597 #if defined HYPRE 00598 IF (lNeedMeanDensity) CALL UpdateMeanDensity(n) 00599 IF (lElliptic) CALL Elliptic(n) 00600 #endif 00601 IF(lExplicit) CALL Explicit(n) 00602 CALL ParticlePostUpdate(n) 00603 CALL PrintAdvance(n) 00604 END IF 00605 ELSE IF (iThreaded == PSEUDO_THREADED) THEN 00606 IF (lHasData) CALL ScheduledAdvanceGrids(n) 00607 IF (lHasChildren) CALL AMR(n+1) 00608 IF (lHasData) CALL CompleteAdvanceGrids(n) 00609 #if defined PTHREADS || defined PTH 00610 ELSEIF (iThreaded == THREADED) THEN 00611 IF (lHasData) CALL LaunchAdvanceThread(n) 00612 IF (lHasChildren) CALL AMR(n+1) 00613 IF (lHasData) CALL JoinAdvanceThread(n) 00614 # endif 00615 END IF 00616 00617 00618 ! ========================================================! 00619 ! Stage 2e ! 00620 ! Synchronize Data with children ! 00621 ! ========================================================! 00622 IF (lHasChildren) THEN 00623 IF (lHasData) THEN 00624 CALL ApplyChildrenData(n) 00625 CALL CompSendParentsData(n+1) 00626 CALL CompRecvChildrenData(n) 00627 END IF 00628 END IF 00629 00630 00631 ! ========================================================! 00632 ! Stage 2e ! 00633 ! Synchronize Data with neighbors ! 00634 ! and if threaded do some advances to kill time ! 00635 ! ========================================================! 00636 00637 IF (lHasData) THEN 00638 CALL RestrictionFixups(n) 00639 CALL AfterFixups(n) 00640 CALL PostRecvFluxes(n) 00641 CALL PostSendFluxes(n) 00642 00643 IF (iThreaded == PSEUDO_THREADED .OR. iThreaded == THREADED) THEN 00644 IF (iThreaded == PSEUDO_THREADED .AND. n > 0) CALL WaitingAdvances(n) 00645 #if defined PTH 00646 IF (iThreaded == THREADED .AND. (lElliptic .OR. lParticles)) CALL pth_MPI_BARRIER(MPI_COMM_WORLD, ierr) 00647 #endif 00648 #if defined HYPRE 00649 IF (lElliptic) CALL Elliptic(n) 00650 #endif 00651 IF(lExplicit) CALL Explicit(n) 00652 IF (lParticles) CALL ParticlePostUpdate(n) 00653 CALL PrintAdvance(n) 00654 END IF 00655 IF (lHasChildren .AND. lRegridLevel(n+1)) CALL UpdateChildMasks(n) !In order to synchronize fluxes we need to know about neighbors children 00656 CALL SyncFluxes(n) 00657 CALL CompSendFluxes(n) 00658 CALL CompRecvFluxes(n) 00659 IF (lHasParentswData) CALL AccumulateFluxes(n) 00660 END IF 00661 ! IF (step == 2) CALL NullifyNeighbors(n) !Since neighbors became overlaps in UpdateOverlaps - we want to avoid deallocating the list twice 00662 IF (RestartStep) EXIT 00663 END IF 00664 END DO 00665 00666 00667 ! ========================================================! 00668 ! Stage 2f ! 00669 ! Coarsen Data for parents ! 00670 ! ========================================================! 00671 IF (lAdvance) THEN 00672 IF (lHasParentswData) THEN 00673 CALL CoarsenDataForParents(n) 00674 CALL PostSendParentsData(n) 00675 END IF 00676 END IF 00677 CALL StopTimer(iAMR, n) 00678 END SUBROUTINE AMR 00679 00680 00681 00682 00687 RECURSIVE SUBROUTINE AMRStart(n) 00688 USE TreeLevelComms 00689 USE DataLevelComms 00690 INTEGER :: n, nSteps, step, ierr 00691 LOGICAL :: lHasData, lHasChildren, lHasParents, lHasParentswData, lHasChildrenwData, lOverlapsAreOld 00692 00693 lHasData = (n >= 0) 00694 lHasParentswData = (n >= 1 .AND. .NOT. lRestart) 00695 lHasChildren = ((.NOT. lRestart .AND. n < MaxLevel) .OR. (lRestart .AND. n < RestartLevel)) 00696 lHasParents = (n > -2) 00697 nsteps=1 00698 step=1 00699 levels(n)%step=1 00700 levels(n+1:)%step=1 00701 IF (n > -1 .AND. MPI_ID == 0) THEN 00702 IF (lRestart) THEN 00703 write(*,*) 'Reloading Grids on level', n 00704 ELSE 00705 write(*,*) 'Initializing Grids on level', n 00706 END IF 00707 END IF 00708 ! ========================================================! 00709 ! Stage 1 ! 00710 ! Complete communication with parents and get level load ! 00711 ! ========================================================! 00712 00713 IF (lHasParents) THEN 00714 CALL CompRecvGridsFromParents(n) 00715 CALL CompSendGridsToChildren(n-1) 00716 CALL PostRecvOverlapsNeighbors(n) 00717 IF (lHasParentswData) CALL PostRecvParentsData(n) 00718 CALL SortNodes(n) 00719 00720 IF (lRestart) THEN 00721 IF (lHasData) CALL InitInfos(n) 00722 ELSE 00723 IF (lHasData) CALL InitialInitInfos(n) 00724 END IF 00725 00726 CALL CompRecvOverlapsNeighbors(n) 00727 CALL CompSendOverlapsNeighbors(n-1) 00728 00729 IF (lHasParentswData) THEN 00730 CALL CompRecvParentsData(n) 00731 CALL CompSendChildrenData(n-1) 00732 CALL ProlongateParentsData(n) 00733 END IF 00734 00735 IF (lHasData) CALL ChildMaskOverlaps(n) !need to know about neighbors for elliptic calls 00736 END IF 00737 IF (lHasData) CALL GetLevelLoad(n) 00738 00739 00740 ! ========================================================! 00741 ! Stage 2a ! 00742 ! Initialize grids and ghost data ! 00743 ! ========================================================! 00744 00745 IF (lRestart) THEN 00746 CALL IOReloadLevel(n) 00747 ELSE 00748 IF (lHasData) THEN 00749 IF (lHasParentswData) CALL AfterOverlaps(n) 00750 CALL InitGrids(n) 00751 CALL BeforeGlobalStep(n) 00752 IF (lNeedMeanDensity) CALL UpdateMeanDensity(n) 00753 levels(n)%step=levels(n)%steps 00754 CALL UpdateOverlaps(n) ! On last step - ghost zones are filled from neighbors - so overlaps=>neighbors 00755 CALL PostRecvOverlaps(n) 00756 CALL PostSendOverlaps(n,.true. ) 00757 CALL ApplyOverlaps(n,levels(n)%step) !Ghost zones are filled from neighbors 00758 CALL CompSendOverlaps(n) 00759 CALL CompRecvOverlaps(n) 00760 levels(n)%step=step 00761 CALL ApplyPhysicalBCs(n) 00762 CALL ParticlePreUpdate(n) 00763 CALL ParticlePostUpdate(n) 00764 # if defined HYPRE 00765 IF (lElliptic) THEN 00766 CALL ApplyEllipticBC(n) 00767 CALL InitialElliptic(n) 00768 END IF 00769 # endif 00770 IF (lExplicit) CALL InitialExplicit(n) 00771 END IF 00772 END IF 00773 00774 ! ========================================================! 00775 ! Stage 2b ! 00776 ! Create children ! 00777 ! ========================================================! 00778 00779 IF (lHasChildren) THEN 00780 IF (lRestart) THEN 00781 CALL DistributeChildrens(n, .false.) !no splitting allowed 00782 ELSE 00783 IF (lHasData) CALL SetErrFlags(n) 00784 CALL CreateChildrens(n) 00785 CALL DistributeChildrens(n) 00786 END IF 00787 00788 ! ========================================================! 00789 ! Stage 2c ! 00790 ! Communicate child information with neighbors ! 00791 ! Also post send/recv of child information to child procs ! 00792 ! ========================================================! 00793 00794 CALL PostRecvNeighboringChildren(n) 00795 CALL PostSendNeighboringChildren(n) 00796 CALL InheritNeighborsChildren(n) 00797 CALL CompRecvNeighboringChildren(n) 00798 CALL CompSendNeighboringChildren(n) 00799 00800 CALL PostSendGridsToChildren(n) 00801 CALL PostRecvGridsFromParents(n+1) 00802 00803 CALL PostSendOverlapsNeighbors(n) 00804 00805 IF (.NOT. lRestart) THEN 00806 IF (lHasData) CALL ClearChildFixups(n) 00807 IF (lHasData) CALL PostSendChildrenData(n) 00808 IF (lHasData) CALL PostRecvInitialChildrenData(n) 00809 END IF 00810 00811 ! ========================================================! 00812 ! Stage 2d ! 00813 ! CALL AMRStart ! 00814 ! ========================================================! 00815 00816 IF (lHasChildren) CALL AMRStart(n+1) 00817 00818 00819 ! ========================================================! 00820 ! Stage 2e ! 00821 ! Synchronize Data with children ! 00822 ! ========================================================! 00823 00824 IF (.NOT. lRestart) THEN 00825 IF (lHasChildren) THEN 00826 IF (lHasData) THEN 00827 CALL ApplyInitialChildrenData(n) 00828 CALL CompSendParentsInitialData(n+1) 00829 CALL CompRecvInitialChildrenData(n) 00830 END IF 00831 END IF 00832 END IF 00833 ! ========================================================! 00834 ! Stage 2e ! 00835 ! Synchronize Data with neighbors ! 00836 ! ========================================================! 00837 END IF 00838 IF (lHasData) CALL AfterFixups(n) 00839 ! IF (n > BaseLevel) CALL NullifyNeighbors(n) 00840 00841 ! ========================================================! 00842 ! Stage 2f ! 00843 ! Coarsen Data for parents ! 00844 ! ========================================================! 00845 00846 IF (lHasParentswData) THEN 00847 CALL CoarsenInitialDataForParents(n) 00848 CALL PostSendParentsInitialData(n) 00849 END IF 00850 00851 END SUBROUTINE AMRStart 00852 00853 00855 00856 00859 00862 SUBROUTINE ClearAllNodeLists(first_layer_gone) 00863 00864 INTEGER :: first_layer_gone 00865 INTEGER :: n 00866 00867 DO n = first_layer_gone, MaxLevel 00868 CALL AgeNodes(n) 00869 ! CALL NullifyNodes(n) 00870 END DO 00871 00872 END SUBROUTINE ClearAllNodeLists 00873 00875 SUBROUTINE DestroyAllNodes() 00876 INTEGER :: n 00877 DO n = MaxLevel, -2, -1 00878 CALL DestroyOldNodes(n) 00879 END DO 00880 DO n = MaxLevel, -2, -1 00881 CALL DestroyNodeList(BackupNodes(n)%p) 00882 CALL DestroyNodeList(BackupExternalNodes(n)%p) 00883 END DO 00884 00885 DO n = MaxLevel, -2, -1 00886 !IF (n <= BaseLevel) CALL NullifyNeighbors(n) 00887 CALL DestroyNodes(n) 00888 END DO 00889 END SUBROUTINE DestroyAllNodes 00890 00891 00892 00896 SUBROUTINE BackupData() 00897 INTEGER :: n 00898 TYPE(NodeDefList), POINTER :: nodelist 00899 CALL StartTimer(iBackupData, BaseLevel) 00900 ! First make backup copies of all nodes and their info structures... 00901 DO n=BaseLevel, MaxLevel 00902 !Destroy Previous Backups 00903 CALL DestroyNodeList(BackupNodes(n)%p) 00904 CALL DestroyNodeList(BackupExternalNodes(n)%p) 00905 !Make new Backups 00906 CALL BackupNodelist(Nodes(n)%p, BackupNodes(n)%p, .false.) 00907 CALL BackupNodelist(ExternalNodes(n)%p, BackupExternalNodes(n)%p, .false.) 00908 END DO 00909 ! Now make a backup of any necessary relationships. 00910 ! For static nodes - the parent and neighbor relationships don't need to be backed up 00911 ! However if the children are not static - these connections will be lost unless backed up. 00912 DO n=LastStaticLevel, MaxLevel 00913 nodelist=>BackupNodes(n)%p 00914 DO WHILE (ASSOCIATED(nodelist)) 00915 CALL BackupChildren(n, nodelist%self) 00916 IF (n > LastStaticLevel) CALL BackupParent(n, nodelist%self) 00917 nodelist=>nodelist%next 00918 END DO 00919 nodelist=>BackupExternalNodes(n)%p 00920 DO WHILE (ASSOCIATED(nodelist)) 00921 CALL BackupChildren(n, nodelist%self) 00922 nodelist=>nodelist%next 00923 END DO 00924 END DO 00925 CALL SinkParticleBackup() 00926 CALL StopTimer(iBackupData, BaseLevel) 00927 END SUBROUTINE BackupData 00928 00929 00930 00934 SUBROUTINE RestoreData() 00935 INTEGER :: n 00936 TYPE(NodeDefList), POINTER :: nodelist,backedupnodelist 00937 TYPE(NodeDef), POINTER :: node 00938 TYPE(InfoDef), POINTER :: tempinfo 00939 ! First restore info structures of permanent nodes 00940 DO n=0, min(LastStaticLevel, MaxLevel) 00941 nodelist=>Nodes(n)%p 00942 backedupnodelist=>BackupNodes(n)%p 00943 DO WHILE (ASSOCIATED(nodelist)) 00944 tempinfo=>nodelist%self%info 00945 NULLIFY(nodelist%self%info) 00946 CALL BackupInfo(backedupnodelist%Self%info, nodelist%self%info, .true.) 00947 IF (LastStaticLevel > n) THEN !Keep childfixups 00948 nodelist%self%info%childfixups=>tempinfo%childfixups 00949 nullify(tempinfo%childfixups) 00950 END IF 00951 CALL DestroyInfo(tempinfo) 00952 nodelist=>nodelist%next 00953 backedupnodelist=>backedupnodelist%next 00954 END DO 00955 END DO 00956 00957 ! Now restore old nodes 00958 DO n=LastStaticLevel+1, MaxLevel 00959 CALL DestroyNodeList(Nodes(n)%p) 00960 NULLIFY(Nodes(n)%p, LastLocalNode(n)%p) 00961 CALL BackupNodelist(BackupNodes(n)%p,Nodes(n)%p, .true., LastLocalNode(n)%p) 00962 CALL DestroyNodeList(ExternalNodes(n)%p) 00963 NULLIFY(ExternalNodes(n)%p, LastExternalNode(n)%p) 00964 CALL BackupNodelist(BackupExternalNodes(n)%p,ExternalNodes(n)%p, .true.,LastExternalNode(n)%p) 00965 END DO 00966 00967 ! First clear current child connections for last static level 00968 ! Now re-establish parent and child connections between restored nodes 00969 DO n=LastStaticLevel, MaxLevel 00970 nodelist=>Nodes(n)%p 00971 DO WHILE (ASSOCIATED(nodelist)) 00972 node=>nodelist%self 00973 IF (n == LastStaticLevel) THEN 00974 CALL ClearNodeList(node%children) 00975 NULLIFY(node%lastchild) 00976 END IF 00977 IF (n > LastStaticLevel) CALL RestoreParent(n, node) 00978 CALL RestoreChildren(n, node) 00979 nodelist=>nodelist%next 00980 END DO 00981 00982 ! Possible that an external node would have been added since backup because of a child being distributed locally 00983 nodelist=>ExternalNodes(n)%p 00984 DO WHILE (ASSOCIATED(nodelist)) 00985 node=>nodelist%self 00986 IF (n == LastStaticLevel) THEN 00987 CALL ClearNodeList(node%children) 00988 NULLIFY(node%lastchild) 00989 END IF 00990 CALL RestoreChildren(n, node, n /= LastStaticLevel) 00991 nodelist=>nodelist%next 00992 END DO 00993 END DO 00994 CALL SinkParticleRestore() 00995 END SUBROUTINE RestoreData 00996 00997 00998 00999 SUBROUTINE UpdateLevelStep(n, step, lAdvance) 01000 INTEGER :: n, step 01001 LOGICAL :: lAdvance 01002 IF (lAdvance) THEN 01003 IF (n > ROOT_LEVEL) THEN 01004 levels(n)%CurrentLevelStep=levels(n-1)%CoarsenRatio*(levels(n-1)%CurrentLevelStep-1)+step 01005 ELSE 01006 levels(n)%CurrentLevelStep=levels(n)%CurrentLevelStep+1 01007 END IF 01008 END IF 01009 levels(n)%step=step 01010 levels(n+1:)%step=1 01011 END SUBROUTINE UpdateLevelStep 01012 01014 SUBROUTINE MpiTest 01015 USE TreeLevelComms 01016 CALL PackTest 01017 END SUBROUTINE MpiTest 01018 01019 01020 01021 SUBROUTINE PrintAllocations 01022 REAL(KIND=qPREC), DIMENSION(3) :: totalallocators, maxallocators 01023 REAL(KIND=qPREC) :: totalcells(0:MaxDepth), cupspcpu, walltimeleft, cuprs 01024 INTEGER :: iErr, i 01025 REAL(KIND=qPREC), DIMENSION(:), ALLOCATABLE :: TotalWorkDone 01026 REAL(KIND=qPREC), DIMENSION(:), ALLOCATABLE :: TotalEllipticWorkDone 01027 REAL(KIND=qPREC), DIMENSION(:), ALLOCATABLE :: TotalInternalCellUpdates 01028 REAL(KIND=qPREC), DIMENSION(:), ALLOCATABLE :: TotalExternalCellUpdates 01029 ! CALL MPI_ALLREDUCE(maxallocation, maxallocators, 3, MPI_DOUBLE_PRECISION, MPI_MAX, MPI_COMM_WORLD, iErr) 01030 ! CALL MPI_ALLREDUCE(allocator, totalallocators, 3, MPI_DOUBLE_PRECISION, MPI_SUM, MPI_COMM_WORLD, iErr) 01031 ! CALL MPI_ALLREDUCE(NumCellsByLevel, totalcells, MaxDepth+1, MPI_DOUBLE_PRECISION, MPI_SUM, MPI_COMM_WORLD, iErr) 01032 ! IF (totalallocators(InfoAllocator) > 25.0*1024.0**3) THEN 01033 ! lPrintDebugFrame=.true. 01034 ! END IF 01035 ALLOCATE(TotalExternalCellUpdates(0:MaxLevel)) 01036 ALLOCATE(TotalInternalCellUpdates(0:MaxLevel)) 01037 ALLOCATE(TotalWorkDone(0:MaxLevel)) 01038 ALLOCATE(TotalEllipticWorkDone(0:MaxLevel)) 01039 01040 CALL MPI_REDUCE(maxallocation, maxallocators, 3, MPI_DOUBLE_PRECISION, MPI_MAX, 0, MPI_COMM_WORLD, iErr) 01041 CALL MPI_REDUCE(allocator, totalallocators, 3, MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, iErr) 01042 CALL MPI_REDUCE(NumCellsByLevel, totalcells, MaxDepth+1, MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, iErr) 01043 CALL MPI_REDUCE(REAL(InternalCellUpdates, KIND=qPREC), TotalInternalCellUpdates, MaxLevel+1, MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, iErr) 01044 CALL MPI_REDUCE(REAL(CellUpdates, KIND=qPREC), TotalExternalCellUpdates, MaxLevel+1, MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, iErr) 01045 CALL MPI_REDUCE(Timers(iAdvanceGrids)%Accumulator(0:MaxLevel), TotalWorkDone, MaxLevel+1, MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, iErr) 01046 IF (lElliptic) CALL MPI_REDUCE(Timers(iElliptic)%Accumulator(0:MaxLevel), TotalEllipticWorkDone, MaxLevel+1, MPI_DOUBLE_PRECISION, MPI_SUM, 0, MPI_COMM_WORLD, iErr) 01047 01048 IF (MPI_ID == 0) THEN 01049 cupspcpu=REAL(SUM(TotalInternalCellUpdates))/(MPI_NP*(mpi_wtime()-StartTime)) 01050 cuprs=REAL(SUM(totalcells(0:MaxLevel)*2**(/(i,i=0,MaxLevel)/))) !/cupspcpu/REAL(MPI_NP) 01051 walltimeleft=(final_time-levels(0)%tnow)/max(levels(0)%dt,olddt) * cuprs / (cupspcpu*MPI_NP) 01052 01053 write(*,'(A23,2A10)') 'Info allocations = ', printsize(totalallocators(1)), printsize(maxallocators(1)) 01054 write(*,'(A23,2A10)') 'message allocations = ', printsize(totalallocators(2)), printsize(maxallocators(2)) 01055 write(*,'(A23,2A10)') 'sweep allocations = ', printsize(totalallocators(3)), printsize(maxallocators(3)) 01056 write(*,'(A23,100F7.3)') 'filling fractions = ', (/(merge(real(totalcells(i))/real(totalcells(i-1)*levels(i-1)%CoarsenRatio**nDim),0.0,totalcells(i-1)>0), i=1, MaxLevel)/) 01057 IF (lElliptic) THEN 01058 write(*,'(A23,I3,A2,I3,A2,I3,A2)') 'Current efficiency = ', nint(SUM(TotalWorkDone)/(MPI_NP*(mpi_wtime()-StartTime))*100d0), '% ',nint(SUM(TotalEllipticWorkDone)/(MPI_NP*(mpi_wtime()-StartTime))*100d0), '% ',nint(SUM(TotalWorkDone+TotalEllipticWorkDone)/(MPI_NP*(mpi_wtime()-StartTime))*100d0), '% ' 01059 ELSE 01060 write(*,'(A23,I3,A2)') 'Current efficiency = ', nint(SUM(TotalWorkDone)/(MPI_NP*(mpi_wtime()-StartTime))*100d0), '% ' 01061 END IF 01062 write(*,'(A23,2I10,I4,A1)') 'Cell updates/second = ', nint(cupspcpu), nint(REAL(SUM(TotalExternalCellUpdates))/(MPI_NP*(mpi_wtime()-StartTime))), nint(100*SUM(TotalInternalCellUpdates)/SUM(TotalExternalCellUpdates)), '%' 01063 write(*,'(A23,A11,A10,F6.1,A4,I6)') 'Wall Time Remaining = ', printtime(walltimeleft), ' at frame ', (levels(0)%tnow-start_time)/(final_time-start_time)*REAL(final_frame-start_frame)+start_frame, ' of ',final_frame 01064 01065 write(*,'(A23,E16.4)') 'AMR Speed-Up Factor = ', .95*AdvanceCost(max(1,nint(REAL(levels(MaxLevel)%mX)/(REAL(MPI_NP)**(1d0/REAL(nDim))))))*2d0**MaxLevel/(cuprs/(cupspcpu*MPI_NP)) 01066 END IF 01067 01068 01069 DEALLOCATE(TotalWorkDone, TotalInternalCellUpdates, TotalEllipticWorkDone, TotalExternalCellUpdates) 01070 01071 END SUBROUTINE PrintAllocations 01072 01074 01075 01076 01077 END MODULE AmrControl 01078