Scrambler  1
amr_control.f90
Go to the documentation of this file.
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 
 All Classes Files Functions Variables