Scrambler  1
tree_declarations.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 !    tree_declarations.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 
00029 
00039 MODULE TreeDeclarations
00040    USE DataDeclarations
00041    USE GlobalDeclarations
00042 
00043    IMPLICIT NONE
00044    PRIVATE
00045 
00046    ! Public Data types
00047    PUBLIC :: NodeBox, NodeBoxList, NodeDef, NodeDefList, pNodeDef, pNodeDefList
00048 
00049    !Routines for creating nodes
00050    PUBLIC :: AddNode, AddOldNode, AddFindNode, AddFindOldNode
00051 
00052    !Routines for destroying nodes
00053    PUBLIC :: NullifyNodeFields, DestroyNode, BackupNode
00054 
00055    !Routines for adding connections to tree
00056    PUBLIC :: AddParent, AddOverlap, AddFindOverlap, AddNeighbor, AddChild, AddFindChild
00057 
00058    !Routines for managing nodelists
00059    PUBLIC :: AddNodeToList, ClearNodeList, DestroyNodeList, BackupNodeList
00060 
00061    !Routines for finding nodes
00062    PUBLIC :: FindNode, StrictFindNode, FindOldNode, StrictFindOldNode, GetChildID, FindNodeInList, FindBackupNode, StrictFindBackupNode, FindParent, StrictFindParent, StrictFindAnyNode
00063 
00064    !Miscellaneous Routines
00065    PUBLIC :: isAncestor, GetFinestLevel, NodeCount, traverse, CellCount
00066 
00067    !Routines for determining necessary connections in the tree
00068    PUBLIC :: ChildCanBeNeighbor, NephewCanBeNeighbor, NephewCanBeOverlap, Overlaps, Neighbors
00069 
00070    !Routines for manipulating node box indices
00071    PUBLIC :: GetChildMGlobal, GetChildMBounds
00072 
00073    !Routines for manipulating node boxes
00074    PUBLIC :: CreateNodeBox, DestroyNodeBox, MatchBox
00075 
00076    !Routines for manipulating nodeboxlists
00077    PUBLIC :: AddNodeBoxToList, ClearNodeBoxList, NodeBoxCount
00078 
00079    SAVE
00080 
00082    TYPE NodeBox
00083       SEQUENCE
00084       INTEGER :: mGlobal(3,2)
00085       INTEGER :: MPI_ID
00086    End TYPE NodeBox
00087 
00089    TYPE NodeBoxList
00090       SEQUENCE
00091       TYPE(NodeBox) :: self
00092       TYPE(NodeBoxList), POINTER :: next
00093    End TYPE NodeBoxList
00094 
00096    TYPE NodeDef
00097       SEQUENCE
00098       TYPE(NodeBox) :: box
00099       TYPE(NodeDef), POINTER :: parent
00100       TYPE(NodeDefList), POINTER :: children
00101       TYPE(NodeDefList), POINTER :: oldchildren
00102       TYPE(NodeDefList), POINTER :: neighbors
00103       TYPE(NodeDefList), POINTER :: overlaps
00104       TYPE(NodeDefList), POINTER :: lastchild
00105       TYPE(NodeDefList), POINTER :: lastneighbor
00106       TYPE(NodeDefList), POINTER :: lastoverlap
00107       REAL, POINTER :: ProcTime(:)
00108       INTEGER, POINTER :: ProcList(:)
00109       TYPE(InfoDef), POINTER :: Info
00110       REAL(KIND=qPrec) :: stamp
00111       INTEGER :: ID
00112    End TYPE NodeDef
00113 
00115    TYPE NodeDefList
00116       SEQUENCE
00117       TYPE(NodeDef), POINTER :: self
00118       TYPE(NodeDefList), POINTER :: next
00119    End TYPE NodeDefList
00120 
00122    TYPE pNodeDefList
00123       SEQUENCE
00124       TYPE(NodeDefList), POINTER :: p
00125    End TYPE pNodeDefList
00126 
00128    TYPE pNodeDef
00129       SEQUENCE
00130       TYPE(NodeDef), POINTER :: p
00131    End TYPE pNodeDef
00132 
00133    TYPE(pNodeDefList), DIMENSION(:), POINTER, PUBLIC :: Nodes, OldNodes, ExternalNodes, OldExternalNodes, BackupNodes, BackupExternalNodes
00134    TYPE(pNodeDefList), DIMENSION(:), POINTER, PUBLIC :: LastLocalNode, LastExternalNode, LastOldLocalNode, LastOldExternalNode
00135 
00136    INTEGER :: nNodes=0
00137    !=========
00138 CONTAINS 
00139    !=========
00140 
00143 
00150    SUBROUTINE AddNode(level, box, node, proclist, proctime)
00151 
00152       INTEGER :: level
00153       TYPE(NodeBox) :: box
00154       INTEGER, OPTIONAL, DIMENSION(:), POINTER :: proclist
00155       REAL, OPTIONAL, DIMENSION(:), POINTER :: proctime
00156       TYPE(NodeDef), POINTER :: node
00157       TYPE(NodeDefList), POINTER :: LastNode, nodelist
00158       INTEGER :: iErr
00159 
00160       IF (ASSOCIATED(node)) THEN
00161          PRINT *, "node already associated!!!"
00162          STOP
00163       END IF
00164 
00165       !      PRINT "('AddNode(', i2, ') = [', 6i4, ']')", level, box%mGlobal
00166       ! Nodes where the data is going to reside on another processor should go
00167       ! in the external nodes list.
00168       IF (box%MPI_ID == MPI_ID) THEN
00169          CALL AddNodeToList(node,LastLocalNode(level)%p, Nodes(level)%p)
00170       ELSE
00171          CALL AddNodeToList(node,LastExternalNode(level)%p, ExternalNodes(level)%p)
00172       END IF
00173       node%box=box
00174       !  Assign a processor pool and timeslot to the node if they are supplied.
00175       IF (PRESENT(proclist)) THEN
00176          node%proclist=>proclist
00177          IF (PRESENT(proctime)) THEN
00178             node%proctime=>proctime
00179          ELSE
00180             IF (ASSOCIATED(node%proctime)) THEN
00181                PRINT*, "node%proctime already associated"
00182                STOP
00183             END IF
00184             NULLIFY(node%proctime)
00185             ALLOCATE(node%proctime(SIZE(proclist)), STAT=iErr)
00186             IF (iErr /= 0) THEN
00187                PRINT *, "AddNode error: unable to allocate proctime list."
00188                STOP
00189             END IF
00190             node%proctime = 0.0
00191          END IF
00192       END IF
00193 
00194       ! If box data is to reside on this processor, then create a InfoDef structure.
00195 
00196       ! Let's do this later after we've built the tree so we can pass in parents? jjc
00197       !    IF (box%MPI_ID == MPI_ID) THEN
00198       !        ALLOCATE(node%Info)
00199       !        CALL InitInfo(node%info)
00200       !    END IF
00201       nNodes=nNodes+1
00202       node%iD=nNodes
00203 !      IF (box%MPI_id == MPI_id) THEN
00204 !         PRINT "('***AddNode(', i2, '), [', 6i4, ']) success = ', l2, ' ID=', i4,i4,i4' .')", level, &
00205 !              box%mGlobal, ASSOCIATED(Nodes(level)%p), node%id, node%box%MPI_ID
00206 !      ELSE
00207 !         PRINT "('***AddExternalNode(', i2, '), [', 6i4, ']) success = ', l2, ' ID=', i4,i4,i4' .')", level, &
00208 !              box%mGlobal, ASSOCIATED(ExternalNodes(level)%p), node%id, node%box%MPI_ID
00209 !      END IF
00210       CALL CPU_TIME(node%stamp)
00211    END SUBROUTINE AddNode
00212 
00213 
00220    SUBROUTINE AddOldNode(level, box, node, proclist, proctime)
00221 
00222       INTEGER :: level
00223       TYPE(NodeBox) :: box
00224       INTEGER, OPTIONAL, DIMENSION(:), POINTER :: proclist
00225       REAL, OPTIONAL, DIMENSION(:), POINTER :: proctime
00226       TYPE(NodeDef), POINTER :: node
00227       TYPE(NodeDefList), POINTER :: LastNode, nodelist
00228       INTEGER :: iErr
00229 
00230       IF (ASSOCIATED(node)) THEN
00231          PRINT *, "node already associated!!!"
00232          STOP
00233       END IF
00234 
00235       !      PRINT "('AddNode(', i2, ') = [', 6i4, ']')", level, box%mGlobal
00236       ! Nodes where the data is going to reside on another processor should go
00237       ! in the external nodes list.
00238       IF (box%MPI_ID == MPI_ID) THEN
00239          CALL AddNodeToList(node,LastOldLocalNode(level)%p, OldNodes(level)%p)
00240          !         LastNode=>LastLocalNode(level)%p
00241       ELSE
00242          CALL AddNodeToList(node,LastOldExternalNode(level)%p, OldExternalNodes(level)%p)
00243          !         LastNode=>LastExternalNode(level)%p
00244       END IF
00245       node%box=box
00246       !  Assign a processor pool and timeslot to the node if they are supplied.
00247       IF (PRESENT(proclist)) THEN
00248          node%proclist=>proclist
00249          IF (PRESENT(proctime)) THEN
00250             node%proctime=>proctime
00251          ELSE
00252             IF (ASSOCIATED(node%proctime)) THEN
00253                PRINT*, "node%proctime already associated"
00254                STOP
00255             END IF
00256             NULLIFY(node%proctime)
00257             ALLOCATE(node%proctime(SIZE(proclist)), STAT=iErr)
00258             IF (iErr /= 0) THEN
00259                PRINT *, "AddNode error: unable to allocate proctime list."
00260                STOP
00261             END IF
00262             node%proctime = 0.0
00263          END IF
00264       END IF
00265 
00266       ! If box data is to reside on this processor, then create a InfoDef structure.
00267 
00268       ! Let's do this later after we've built the tree so we can pass in parents? jjc
00269       !    IF (box%MPI_ID == MPI_ID) THEN
00270       !        ALLOCATE(node%Info)
00271       !        CALL InitInfo(node%info)
00272       !    END IF
00273       nNodes=nNodes+1
00274       node%iD=nNodes
00275 !            IF (box%MPI_id == MPI_id) THEN
00276 !               PRINT "('***AddOldNode(', i2, '), [', 6i4, ']) success = ', l2, ' ID=', i4, i4, ' .')", level, &
00277 !                      box%mGlobal, ASSOCIATED(OldNodes(level)%p), node%id, node%box%mpi_id
00278 !            ELSE
00279 !               PRINT "('***AddOldExternalNode(', i2, '), [', 6i4, ']) success = ', l2, ' ID=', i4,i4, ' .')", level, &
00280 !                      box%mGlobal, ASSOCIATED(OldExternalNodes(level)%p), node%id, node%box%mpi_id
00281 !            END IF
00282 
00283    END SUBROUTINE AddOldNode
00284 
00285 
00286 
00293 
00294    SUBROUTINE AddFindNode(level, box, node, proclist, proctime)
00295 
00296       INTEGER :: level
00297       TYPE(NodeBox) :: box
00298       TYPE(NodeDef), POINTER :: node
00299       !TYPE(NodeDef), POINTER, INTENT(OUT) :: node
00300       INTEGER, OPTIONAL, DIMENSION(:), POINTER :: proclist
00301       REAL, OPTIONAL, DIMENSION(:), POINTER :: proctime
00302 
00303 
00304       CALL FindNode(level, box, node)
00305 
00306       IF (ASSOCIATED(node)) THEN
00307          ! If node is found and there is a processor pool/time list supplied,
00308          ! then assign the optional lists.
00309 
00310          IF (present(proclist)) THEN
00311             node%proclist=>proclist
00312             IF (present(proctime))  node%proctime=>proctime
00313          END IF
00314       ELSE
00315          ! If node is not found, then create it.
00316          IF(present(proclist)) THEN
00317             IF (present(proctime)) THEN
00318                CALL AddNode(level,box,node,proclist,proctime)
00319             ELSE
00320                CALL AddNode(level,box,node,proclist)
00321             END IF
00322          ELSE
00323             CALL AddNode(level,box,node)
00324          END IF
00325       END IF
00326 
00327    END SUBROUTINE AddFindNode
00328 
00335    SUBROUTINE AddFindOldNode(level, box, node, proclist, proctime)
00336 
00337       INTEGER :: level
00338       TYPE(NodeBox) :: box
00339       TYPE(NodeDef), POINTER, INTENT(OUT) :: node
00340       INTEGER, OPTIONAL, DIMENSION(:), POINTER :: proclist
00341       REAL, OPTIONAL, DIMENSION(:), POINTER :: proctime
00342 
00343 
00344       CALL FindOldNode(level, box, node)
00345 
00346       IF (ASSOCIATED(node)) THEN
00347          ! If node is found and there is a processor pool/time list supplied,
00348          ! then assign the optional lists.
00349          IF (present(proclist)) THEN
00350             node%proclist=>proclist
00351             IF (present(proctime))  node%proctime=>proctime
00352          END IF
00353       ELSE
00354          ! If node is not found, then create it.
00355          IF(present(proclist)) THEN
00356             IF (present(proctime)) THEN
00357                CALL AddOldNode(level,box,node,proclist,proctime)
00358             ELSE
00359                CALL AddOldNode(level,box,node,proclist)
00360             END IF
00361          ELSE
00362             CALL AddOldNode(level,box,node)
00363          END IF
00364       END IF
00365 
00366    END SUBROUTINE AddFindOldNode
00367 
00369 
00370 
00373 
00376    SUBROUTINE NullifyNodeFields(node)
00377 
00378       TYPE(NodeDef), POINTER :: node
00379 
00380       NULLIFY(node%proclist)
00381       NULLIFY(node%proctime)
00382       NULLIFY(node%parent)
00383       NULLIFY(node%Info)
00384       NULLIFY(node%children, node%lastchild)
00385       NULLIFY(node%oldchildren)
00386       NULLIFY(node%neighbors, node%lastneighbor)
00387       NULLIFY(node%overlaps, node%lastoverlap)
00388 
00389    END SUBROUTINE NullifyNodeFields
00390 
00391 
00394    SUBROUTINE DestroyNode(node)
00395 
00396       TYPE(NodeDef), POINTER :: node
00397       TYPE(InfoDef), POINTER :: Info
00398       
00399       IF (ASSOCIATED(node%children))  CALL ClearNodeList(node%children)
00400       IF (ASSOCIATED(node%oldchildren)) CALL ClearNodeList(node%oldchildren)
00401       !    CALL ClearNodeList(node%parent)
00402       IF (ASSOCIATED(node%overlaps, target=node%neighbors)) THEN
00403          NULLIFY(node%neighbors)
00404       ELSE
00405          CALL ClearNodeList(node%neighbors)
00406       END IF
00407       CALL ClearNodeList(node%overlaps)
00408       NULLIFY(node%lastchild)
00409       NULLIFY(node%lastneighbor)
00410       NULLIFY(node%lastoverlap)
00411       NULLIFY(node%parent)
00412 
00413       IF (ASSOCIATED(node%Info)) CALL DestroyInfo(node%Info)     
00414       IF (ASSOCIATED(node%ProcList)) THEN
00415          DEALLOCATE(node%ProcList)
00416          NULLIFY(node%ProcList)
00417       END IF
00418       IF (ASSOCIATED(node%proctime)) THEN
00419          DEALLOCATE(node%proctime)
00420          NULLIFY(node%proctime)
00421       END IF
00422       DEALLOCATE(node)
00423       NULLIFY(node)
00424 
00425    END SUBROUTINE DestroyNode
00426 
00427 
00431    SUBROUTINE BackupNode(original, node, lRestore)
00432       TYPE(NodeDef), POINTER :: original, node
00433       LOGICAL :: lRestore
00434       node%box=original%box
00435       node%stamp=original%stamp
00436       node%id=original%id
00437       IF (ASSOCIATED(original%proctime)) THEN
00438          ALLOCATE(node%proctime(size(original%proctime)))
00439          node%proctime=original%proctime
00440       END IF
00441       IF (ASSOCIATED(original%proclist)) THEN
00442          ALLOCATE(node%proclist(size(original%proclist)))
00443          node%proclist=original%proclist
00444       END IF
00445       IF (ASSOCIATED(original%info)) THEN
00446          CALL BackupInfo(original%info, node%info, lRestore)
00447       END IF                 
00448    END SUBROUTINE BackupNode
00449 
00450 
00452 
00453 
00454 
00455 
00458 
00462    SUBROUTINE AddParent(node, parent)
00463 
00464       TYPE(NodeDef), POINTER :: node
00465       TYPE(NodeDef), POINTER :: parent
00466 
00467 !          PRINT "('***AddParent():  node ID = ', I4, '; parent ID = ', I4, '.')", node%ID, parent%ID
00468       node%Parent=>parent
00469    END SUBROUTINE AddParent
00470 
00471 
00475    SUBROUTINE AddOverlap(node,overlap)
00476 
00477       TYPE(NodeDef), POINTER :: node, overlap
00478       INTEGER :: iErr
00479 !      PRINT "('***AddOverlap:   Node: ', 2I4, '[', 6I4, '], Overlap: ', 2I4, ', [', 6I4, ']')", node%box%mpi_ID, node%id, node%box%mGlobal, overlap%box%mpi_ID, overlap%id, overlap%box%mGlobal
00480       CALL AddNodeToList(overlap, node%lastoverlap, node%overlaps)
00481 
00482    END SUBROUTINE AddOverlap
00483 
00487    SUBROUTINE AddFindOverlap(node, overlap)
00488       TYPE(NodeDefList), POINTER :: overlaplist
00489       TYPE(NodeDef), POINTER :: node, overlap, temp
00490 
00491       overlaplist=>node%overlaps
00492 
00493       ! Iterate over list to find a box that matches the target box dimensions.
00494       DO WHILE (associated(overlaplist))
00495          temp=>overlaplist%self
00496          IF (MatchBox(temp%box,overlap%box)) EXIT
00497          overlaplist=>overlaplist%next
00498       END DO
00499 
00500       IF (.NOT. ASSOCIATED(overlaplist))  CALL AddOverlap(node,overlap)
00501 
00502    END SUBROUTINE AddFindOverlap
00503 
00504 
00508    SUBROUTINE AddNeighbor(node,neighbor)
00509 
00510       TYPE(NodeDef), POINTER :: node, neighbor
00511       INTEGER :: iErr
00512 !      PRINT "('***AddNeighbor:  Node: ', 2I4, '[', 6I4, '], Neighbor: ', 2I4, ', [', 6I4, ']')", node%box%MPI_ID, node%ID, node%box%mGlobal, neighbor%box%MPI_ID, neighbor%ID, neighbor%box%mGlobal
00513       CALL AddNodeToList(neighbor, node%lastneighbor, node%neighbors)
00514    END SUBROUTINE AddNeighbor
00515 
00516 
00520    SUBROUTINE AddChild(node, child)
00521       TYPE(NodeDef), POINTER :: node, child
00522       INTEGER :: iErr
00523 !         PRINT "('***AddChild:  Node: ', I4, '[', 6I4, '], Child: ', I4, ', [', 6I4, ']')", node%ID, node%box%mGlobal, child%ID, child%box%mGlobal
00524       CALL AddNodeToList(child, node%lastchild, node%children)
00525 
00526    END SUBROUTINE AddChild
00527 
00531    SUBROUTINE AddFindChild(node, child)
00532       TYPE(NodeDefList), POINTER :: childlist
00533       TYPE(NodeDef), POINTER :: node, child, temp
00534 
00535       childlist=>node%children
00536 
00537       ! Iterate over list to find a box that matches the target box dimensions.
00538       DO WHILE (associated(childlist))
00539          temp=>childlist%self
00540          IF (MatchBox(temp%box,child%box)) EXIT
00541          childlist=>childlist%next
00542       END DO
00543 
00544       IF (.NOT. Associated(childlist))  CALL AddChild(node,child)
00545 
00546    END SUBROUTINE AddFindChild
00547 
00549 
00550 
00553 
00558    SUBROUTINE AddNodeToList(node, lastnode, nodelist)
00559       TYPE(NodeDefList), POINTER :: lastnode
00560       TYPE(NodeDefList), POINTER, OPTIONAL :: nodelist
00561       TYPE(NodeDef), POINTER :: node
00562       INTEGER :: iErr
00563       integer :: ierr2
00564 
00565       IF (.NOT. ASSOCIATED(LastNode)) THEN     
00566          ALLOCATE(LastNode, STAT=iErr)
00567          NULLIFY(LastNode%next)
00568          NULLIFY(LastNode%self)
00569          IF (iErr /= 0) THEN
00570             PRINT *, "AddNode() error: unable to allocate LastNode list object."
00571             STOP
00572          END IF
00573          IF (present(nodelist)) nodelist=>lastnode
00574       ELSE
00575          IF (ASSOCIATED(LastNode%next)) THEN
00576             PRINT *, "Error - last node next allocated"
00577             STOP
00578          END IF
00579          ALLOCATE(LastNode%next, STAT=iErr)
00580          IF (iErr /= 0) THEN
00581             PRINT *, "AddNode() error: unable to allocate LastNode%next list object."
00582             STOP
00583          END IF
00584          LastNode=>LastNode%next         
00585          NULLIFY(LastNode%next)
00586          NULLIFY(LastNode%self)
00587 
00588       END IF
00589 
00590       IF (ASSOCIATED(node)) THEN
00591          LastNode%self=>node
00592       ELSE
00593          ALLOCATE(LastNode%self, STAT=iErr)
00594 
00595          IF (iErr /= 0) THEN
00596             PRINT *, "AddNode() error: unable to allocate LastNode%self object."
00597             STOP
00598          END IF
00599          node=>LastNode%self
00600          NULLIFY(node%parent, node%children, node%oldchildren, node%neighbors, node%overlaps, node%lastchild, &
00601               node%lastneighbor, node%lastoverlap, node%proctime, node%proclist, node%info)         
00602       END IF
00603 
00604    END SUBROUTINE AddNodeToList
00605 
00606 
00609    RECURSIVE SUBROUTINE ClearNodeList(nodelist)
00610       TYPE(NodeDefList), POINTER :: nodelist
00611       IF (.NOT. ASSOCIATED(nodelist))  RETURN
00612       IF (ASSOCIATED(nodelist%next)) CALL ClearNodeList(nodelist%next)
00613       DEALLOCATE(nodelist)
00614       NULLIFY(nodelist)
00615    END SUBROUTINE ClearNodeList
00616 
00620    RECURSIVE SUBROUTINE DestroyNodeList(nodelist)
00621 
00622       TYPE(NodeDefList), POINTER :: nodelist
00623       IF (.NOT. ASSOCIATED(nodelist))  RETURN
00624       IF (ASSOCIATED(nodelist%next))  CALL DestroyNodeList(nodelist%next)
00625       IF (ASSOCIATED(nodelist%self))  CALL DestroyNode(nodelist%self)
00626       DEALLOCATE(nodelist)
00627       NULLIFY(nodelist)
00628 
00629    END SUBROUTINE  DestroyNodeList
00630 
00634    SUBROUTINE BackupNodelist(original, backup, lRestore, lastnode)
00635       TYPE(NodeDefList), POINTER :: original, backup, lastbackup, nodelist
00636       TYPE(NodeDef), POINTER :: node
00637       TYPE(NodeDefList), POINTER, OPTIONAL :: lastnode
00638       LOGICAL :: lRestore
00639       !      IF (ASSOCIATED(original)) ALLOCATE(backup)
00640       NULLIFY(lastbackup, backup)
00641       nodelist=>original
00642       DO WHILE (ASSOCIATED(nodelist))
00643          ALLOCATE(node)
00644          CALL NullifyNodeFields(node)
00645          CALL BackupNode(nodelist%self, node, lRestore)
00646          CALL AddNodeToList(node, lastbackup, backup)
00647          nodelist=>nodelist%next
00648       END DO
00649       IF (PRESENT(LastNode)) LastNode=>lastbackup
00650 
00651    END SUBROUTINE BackupNodelist
00652 
00653 
00655 
00658 
00659 
00664    SUBROUTINE FindNode(level,box,node)
00665 
00666       INTEGER :: level
00667       TYPE(NodeBox) :: box
00668       TYPE(NodeDef), POINTER :: node
00669 
00670       IF (box%MPI_ID==MPI_ID) THEN
00671          CALL FindNodeInList(box,Nodes(level)%p, node)
00672       ELSE
00673          CALL FindNodeInList(box,ExternalNodes(level)%p, node)
00674       END IF
00675 
00676    END SUBROUTINE FindNode
00677 
00683    SUBROUTINE StrictFindNode(level, box, node, caller)
00684 
00685       INTEGER :: level
00686       TYPE(NodeBox) :: box
00687       TYPE(NodeDef), POINTER :: node
00688       CHARACTER(*) :: caller
00689 
00690 
00691       CALL FindNode(level, box, node)
00692 
00693       IF (.NOT. ASSOCIATED(node)) THEN
00694          PRINT *, caller, "::StrictFindNode() error: node not associated."
00695          PRINT "('StrictFindNode(box%MPI_id = ', i4, ', level = ', i2, ') failed on box [', 6I5, ']')", box%MPI_id, level, box%mGlobal
00696          STOP
00697       END IF
00698 
00699    END SUBROUTINE StrictFindNode
00700 
00701 
00707 
00708    SUBROUTINE StrictFindAnyNode(level, box, node, caller)
00709 
00710       INTEGER :: level
00711       TYPE(NodeBox) :: box
00712       TYPE(NodeDef), POINTER :: node
00713       CHARACTER(*) :: caller
00714       TYPE(NodeDefList), POINTER :: nodelist
00715       NULLIFY(node)
00716       nodelist=>Nodes(level)%p
00717       DO WHILE (ASSOCIATED(nodelist))
00718          IF (ALL(nodelist%self%box%mGlobal==box%mGlobal)) THEN
00719             node=>nodelist%self
00720             EXIT
00721          END IF
00722          nodelist=>nodelist%next
00723       END DO
00724       IF (.NOT. ASSOCIATED(node)) THEN
00725          
00726          nodelist=>ExternalNodes(level)%p
00727          DO WHILE (ASSOCIATED(nodelist))
00728             IF (ALL(nodelist%self%box%mGlobal==box%mGlobal)) THEN
00729                node=>nodelist%self
00730                EXIT
00731             END IF
00732             nodelist=>nodelist%next
00733          END DO
00734       END IF
00735 
00736       IF (.NOT. ASSOCIATED(node)) THEN
00737          PRINT *, caller, "::StrictFindAnyNode() error: node not associated."
00738          PRINT "('StrictFindAnyNode(box%MPI_id = ', i4, ', level = ', i2, ') failed on box [', 6I5, ']')", box%MPI_id, level, box%mGlobal
00739          STOP
00740       END IF
00741    END SUBROUTINE StrictFindAnyNode
00742 
00743 
00749    SUBROUTINE StrictFindBackupNode(level, box, node, caller)
00750 
00751       INTEGER :: level
00752       TYPE(NodeBox) :: box
00753       TYPE(NodeDef), POINTER :: node
00754       CHARACTER(*) :: caller
00755 
00756 
00757       CALL FindBackupNode(level, box, node)
00758 
00759       IF (.NOT. ASSOCIATED(node)) THEN
00760          PRINT *, caller, "::StrictFindBackupNode() error: node not associated."
00761          PRINT "('StrictFindBackupNode(box%MPI_id = ', i4, ', level = ', i2, ') failed on box [', 6I5, '] on processor ', i5)", box%MPI_id, level, box%mGlobal, MPI_ID
00762          STOP
00763       END IF
00764 
00765    END SUBROUTINE StrictFindBackupNode
00766 
00771    SUBROUTINE FindOldNode(level,box,node)
00772 
00773       INTEGER :: level
00774       TYPE(NodeBox) :: box
00775       TYPE(NodeDef), POINTER :: node
00776 
00777       ! Nodes where the data is going to reside on another processor should go
00778       ! in the external nodes list.
00779       IF (box%MPI_ID==MPI_ID) THEN
00780          CALL FindNodeInList(box,OldNodes(level)%p, node)
00781       ELSE
00782          CALL FindNodeInList(box,OldExternalNodes(level)%p, node)
00783       END IF
00784    END SUBROUTINE FindOldNode
00785 
00786 
00791    SUBROUTINE FindBackupNode(level,box,node)
00792 
00793       INTEGER :: level
00794       TYPE(NodeBox) :: box
00795       TYPE(NodeDef), POINTER :: node
00796 
00797       ! Nodes where the data is going to reside on another processor should go
00798       ! in the external nodes list.
00799       IF (box%MPI_ID==MPI_ID) THEN
00800          CALL FindNodeInList(box,BackupNodes(level)%p, node)
00801       ELSE
00802          CALL FindNodeInList(box,BackupExternalNodes(level)%p, node)
00803       END IF
00804    END SUBROUTINE FindBackupNode
00805 
00806 
00807 
00808    SUBROUTINE FindNodeInList(box,nodelist,node)
00809       TYPE(Nodebox) :: box
00810       TYPE(NodeDef), POINTER :: node, temp
00811       TYPE(NodeDefList), POINTER :: nodelist, traverse
00812 
00813       ! Make sure that the pointer is null by default.
00814       NULLIFY(node)
00815       traverse=>nodelist
00816       ! Iterate over list to find a box that matches the target box dimensions.
00817       DO WHILE (associated(traverse))
00818          temp=>traverse%self
00819 
00820          IF (MatchBox(temp%box,box)) THEN
00821             node=>temp
00822             EXIT
00823          END IF
00824 
00825          traverse=>traverse%next
00826       END DO
00827 
00828    END SUBROUTINE FindNodeInList
00829 
00835    SUBROUTINE StrictFindOldNode(level, box, node, caller)
00836 
00837       INTEGER :: level
00838       TYPE(NodeBox) :: box
00839       TYPE(NodeDef), POINTER :: node
00840       CHARACTER(*) :: caller
00841 
00842       CALL FindOldNode(level, box, node)
00843 
00844       IF (.NOT. ASSOCIATED(node)) THEN
00845          PRINT *, caller, "::StrictFindOldNode() error: node not associated."
00846          PRINT "('StrictFindNode(box%MPI_id = ', i4, ') failed on box [', 6I5, ']')", box%MPI_id, box%mGlobal
00847          STOP
00848       END IF
00849 
00850    END SUBROUTINE StrictFindOldNode
00851 
00855    INTEGER FUNCTION GetChildID(node,child)
00856       TYPE(NodeDefList), POINTER :: childlist
00857       TYPE(NodeDef), POINTER :: node, child
00858       INTEGER :: i
00859       i=0
00860       childlist=>node%children
00861       DO WHILE (ASSOCIATED(childlist))
00862          i=i+1
00863          IF (MatchBox(child%box,childlist%self%box)) THEN
00864             GetChildID=i
00865             RETURN
00866          END IF
00867          childlist=>childlist%next
00868       END DO
00869       PRINT*, "Amr_node.f90: Couldn't find child! Stopping..."
00870       STOP
00871    END FUNCTION GetChildID
00872 
00874 
00875 
00878 
00881    LOGICAL FUNCTION isAncestor(childmask)
00882       INTEGER :: childmask
00883       isAncestor = (childmask < 0 .AND. childmask /= NEIGHBORCHILD)
00884    END FUNCTION isAncestor
00885 
00888    FUNCTION GetFinestLevel() RESULT(level_buf)
00889       INTEGER :: level
00890       INTEGER :: finest_level
00891       INTEGER :: level_buf
00892       INTEGER :: iErr
00893 
00894 
00895       finest_level = Maxlevel
00896 
00897       ! Loop over the levels and retrieve the finest level on this processor.
00898       DO level=0, MaxLevel
00899          IF (NodeCount(Nodes(level)%p) == 0) THEN
00900             finest_level = level-1 
00901             EXIT
00902          END IF
00903       END DO
00904 
00905       ! Gathers the FinestLevel values from all the processors, finds the maximum, and 
00906       ! redistributes the maximum to all processors.
00907       CALL MPI_ALLREDUCE(finest_level, level_buf, 1, MPI_INTEGER, MPI_MAX, MPI_COMM_WORLD, iErr)
00908 
00909    END FUNCTION GetFinestLevel
00910 
00913    INTEGER FUNCTION NodeCount(node_list)
00914 
00915       TYPE(NodeDefList), POINTER :: node_list, temp
00916 
00917       INTEGER :: iCount
00918       iCount=0
00919 
00920       temp=>node_list
00921       DO WHILE (ASSOCIATED(temp))
00922          IF (ASSOCIATED(temp%self))  iCount = iCount + 1
00923          temp => temp%next
00924       END DO
00925 
00926       NodeCount = iCount
00927    END FUNCTION NodeCount
00928 
00929 
00932    INTEGER FUNCTION CellCount(node_list)
00933 
00934       TYPE(NodeDefList), POINTER :: node_list, temp
00935       CellCount=0
00936 
00937       temp=>node_list
00938       DO WHILE (ASSOCIATED(temp))
00939         ! IF (ASSOCIATED(temp%self%box%MPI_ID == MPI_ID))  
00940          CellCount = CellCount + product(temp%self%box%mGlobal(:,2)-temp%self%box%mGlobal(:,1)+1)
00941          temp => temp%next
00942       END DO
00943    END FUNCTION CellCount
00944 
00947    SUBROUTINE traverse(node_list)
00948 
00949       TYPE(NodeDefList), POINTER :: node_list, temp
00950       INTEGER :: ierr 
00951       write(*,*) "Traversing NodeList"
00952       temp=>node_list
00953       DO WHILE (ASSOCIATED(temp))
00954          IF (ASSOCIATED(temp%self)) THEN
00955             !            write(*,'(7I4)') temp%self%box%mGlobal, temp%self%id
00956             PRINT "('Node ID ', I5, '; box = ', 7I4, '.')", temp%self%id, temp%self%box%mGlobal,temp%self%box%mpi_id
00957             IF (.NOT. ASSOCIATED(temp%self%info)) THEN
00958                PRINT *, "WARNING:  No Info associated."
00959             END IF
00960          ELSE
00961             write(*,*) "node not associated"
00962             write(*,*) ierr
00963          END IF
00964 
00965          temp => temp%next
00966       END DO
00967    END SUBROUTINE traverse
00968 
00970 
00973 
00978    LOGICAL FUNCTION ChildCanBeNeighbor(node,child,level)
00979       ! Neighbors need to know about children that will have common fluxes
00980       ! and/or data for setting ghost for the second step on the child level
00981       ! (ie within mbc cells).  Children that are not well within their parents
00982       ! can potentially neighbor neighbor's children.
00983 
00984       TYPE(NodeDef) :: node, child
00985       INTEGER, DIMENSION(3,2) :: mGlobal
00986       INTEGER :: level, mbc
00987 
00988       !    mbc=hyperbolic_mbc  !<-this works         !levels(level)%gmbc(1)
00989       mbc=levels(level+1)%nmbc
00990 
00991       mGlobal=LevelUp(node%box%mGlobal,level)
00992       !      mGlobal(1:nDim,2)=node%box%mGlobal(1:nDim,2)*2
00993       !      mGlobal(1:nDim,1)=(node%box%mGlobal(1:nDim,1)-1)*2
00994 
00995       IF (MaintainAuxArrays) THEN
00996          mGlobal(1:nDim,1)=mGlobal(1:nDim,1)     + mbc + 1
00997          mGlobal(1:nDim,2)=mGlobal(1:nDim,2)     - mbc - 1
00998       ELSE
00999          mGlobal(1:nDim,1)=mGlobal(1:nDim,1)     + mbc
01000          mGlobal(1:nDim,2)=mGlobal(1:nDim,2)     - mbc
01001       END IF
01002 
01003       ChildCanBeNeighbor = .NOT. Within(child%box%mGlobal,mGlobal)
01004 
01005    END FUNCTION ChildCanBeNeighbor
01006 
01007 
01011    LOGICAL FUNCTION Within(mGlobal1,mGlobal2)
01012       INTEGER, DIMENSION(3,2) :: mGlobal1, mGlobal2
01013       Within= ALL(mGlobal1(1:nDim,1) >= mGlobal2(1:nDim,1)) &
01014            .AND. ALL(mGlobal1(1:nDim,2) <= mGlobal2(1:nDim,2))
01015    END FUNCTION Within
01016 
01021    LOGICAL FUNCTION NephewCanBeNeighbor(neighbor,child,neighborlevel)
01022       ! Neighbors need to know about children that will have common fluxes
01023       ! and/or data for setting ghost for the second step on the child level
01024       ! (ie within mbc cells).  This routine checks for children within mbc cells 
01025       ! of potential neighbor children.
01026 
01027       TYPE(NodeDef) :: neighbor, child
01028       INTEGER :: neighborlevel,i,j,k,mbc
01029       INTEGER, DIMENSION(3,2) :: mGlobal, ioffset
01030 
01031 
01032       !    mbc=hyperbolic_mbc   !<- This works   !levels(neighborlevel)%gmbc(1)
01033       mbc=levels(neighborlevel+1)%nmbc
01034       NephewCanBeNeighbor=.false.
01035       ioffset=0
01036       WHERE(lAnyPeriodic(1:nDim)) ioffset(1:nDim,2)=nperiodic_overlaps(1:nDim)
01037       ioffset(1:nDim,1)=-ioffset(1:nDim,2)
01038       DO i=ioffset(1,1),ioffset(1,2)
01039          DO j=ioffset(2,1),ioffset(2,2)
01040             DO k=ioffset(3,1),ioffset(3,2)
01041 
01042                mGlobal(:,:)=LevelUp(neighbor%box%mGlobal(:,:)+SPREAD((/i,j,k/)*levels(neighborlevel)%mx(:),2,2),neighborlevel)
01043 
01044                IF (MaintainAuxArrays) THEN
01045                   mGlobal(1:nDim,1)=mGlobal(1:nDim,1)     - mbc - 1
01046                   mGlobal(1:nDim,2)=mGlobal(1:nDim,2)     + mbc + 1
01047                ELSE
01048                   mGlobal(1:nDim,1)=mGlobal(1:nDim,1)     - mbc
01049                   mGlobal(1:nDim,2)=mGlobal(1:nDim,2)     + mbc
01050                END IF
01051                NephewCanBeNeighbor =  NephewCanBeNeighbor .OR. BoxOverlap(mGlobal,child%box%mglobal)
01052             END DO
01053          END DO
01054       END DO
01055    END FUNCTION NephewCanBeNeighbor
01056 
01061    LOGICAL FUNCTION NephewCanBeOverlap(overlap,child,level)
01062       ! Overlaps need to know about children that will need to send or receive data
01063       ! for the extended ghost zones of new boxs.  This routine checks whether or not
01064       ! the given child can overlap within rmbc of the overlap's potential children.
01065 
01066       TYPE(NodeDef) :: overlap, child
01067       INTEGER, DIMENSION(3,2) :: mGlobal,ioffset
01068       INTEGER :: level,i,j,k,rmbc
01069 
01070       !    NephewCanBeOverlap = .TRUE.
01071       !    RETURN
01072 
01073       !    rmbc=levels(level)%gmbc(levels(level)%step)
01074       !    rmbc=levels(level+1)%gmbc(1) !<- This works
01075       rmbc=levels(level+1)%ombc(1)
01076 
01077       !      rmbc=mbc*levels(level)%CoarsenRatio
01078       NephewCanBeOverlap=.false.
01079       ioffset=0
01080       WHERE(lAnyPeriodic(1:nDim)) ioffset(1:nDim,2)=nperiodic_overlaps(1:nDim)
01081       ioffset(1:nDim,1)=-ioffset(1:nDim,2)
01082       DO i=ioffset(1,1),ioffset(1,2)
01083          DO j=ioffset(2,1),ioffset(2,2)
01084             DO k=ioffset(3,1),ioffset(3,2)
01085 
01086                mGlobal(:,:)=LevelUp(overlap%box%mGlobal(:,:)+SPREAD((/i,j,k/)*levels(level)%mx(:),2,2),level)
01087 
01088                IF (MaintainAuxArrays) THEN
01089                   mGlobal(1:nDim,1)=mGlobal(1:nDim,1)     - rmbc - 1
01090                   mGlobal(1:nDim,2)=mGlobal(1:nDim,2)     + rmbc + 1
01091                ELSE
01092                   mGlobal(1:nDim,1)=mGlobal(1:nDim,1)     - rmbc
01093                   mGlobal(1:nDim,2)=mGlobal(1:nDim,2)     + rmbc
01094                END IF
01095                NephewCanBeOverlap = NephewCanBeOverlap .OR. BoxOverlap(child%box%mGlobal,mGlobal)
01096             END DO
01097          END DO
01098       END DO
01099    END FUNCTION NephewCanBeOverlap
01100 
01105    LOGICAL FUNCTION Overlaps(child,overlapchild,level)
01106       ! This routine checks whether or not the extended region of one box overlaps with the 
01107       ! non-extended region of the other.  Will determine whether or not children need to share
01108       ! data with eachother before or after their two steps.
01109 
01110       TYPE(NodeDef) :: child, overlapchild
01111       INTEGER :: rmbc,level,i,j,k
01112 
01113       INTEGER, DIMENSION(3,2) :: mGlobal,ioffset
01114       !    rmbc=levels(level)%gmbc(levels(level)%step)  !<- Not sure why this worked... Reduced size of ghosted region to not include extended zones......
01115       rmbc=levels(level)%ombc(1)
01116       Overlaps=.false.
01117       ioffset=0
01118       WHERE(lAnyPeriodic(1:nDim)) ioffset(1:nDim,2)=nperiodic_overlaps(1:nDim)
01119       ioffset(1:nDim,1)=-ioffset(1:nDim,2)
01120       DO i=ioffset(1,1),ioffset(1,2)
01121          DO j=ioffset(2,1),ioffset(2,2)
01122             DO k=ioffset(3,1),ioffset(3,2)
01123 
01124                mGlobal(:,:)=overlapchild%box%mGlobal(:,:)+SPREAD((/i,j,k/)*levels(level)%mX(:),2,2)
01125 
01126                IF (MaintainAuxArrays) THEN
01127                   mGlobal(1:nDim,1)=mGlobal(1:nDim,1) - rmbc - 1
01128                   mGlobal(1:nDim,2)=mGlobal(1:nDim,2) + rmbc + 1
01129                ELSE
01130                   mGlobal(1:nDim,1)=mGlobal(1:nDim,1) - rmbc
01131                   mGlobal(1:nDim,2)=mGlobal(1:nDim,2) + rmbc
01132                END IF
01133 
01134 
01135                Overlaps=Overlaps .OR. BoxOverlap(child%box%mGlobal, mGlobal)
01136             END DO
01137          END DO
01138       END DO
01139    END FUNCTION Overlaps
01140 
01145    LOGICAL FUNCTION Neighbors(child,neighborchild,level)
01146       ! This routine checks to see if two boxs are within mbc of each other.  This determins
01147       ! whether boxs will need to synchronize fluxes/emf's and share data in between the
01148       ! first and second step.
01149 
01150       TYPE(NodeDef) :: child, neighborchild
01151 
01152       INTEGER, DIMENSION(3,2) :: mGlobal,ioffset
01153       INTEGER :: i,j,k,level,mbc
01154 
01155       !    mbc=hyperbolic_mbc    !<- This worked     !levels(level)%gmbc(1)
01156       mbc=levels(level)%nmbc
01157       Neighbors=.false.
01158       ioffset=0
01159       WHERE(lAnyPeriodic(1:nDim)) ioffset(1:nDim,2)=nperiodic_overlaps(1:nDim)
01160       ioffset(1:nDim,1)=-ioffset(1:nDim,2)
01161       DO i=ioffset(1,1),ioffset(1,2)
01162          DO j=ioffset(2,1),ioffset(2,2)
01163             DO k=ioffset(3,1),ioffset(3,2)
01164 
01165                mGlobal(:,:)=neighborchild%box%mGlobal(:,:)+SPREAD((/i,j,k/)*levels(level)%mX(:),2,2)
01166 
01167                IF (MaintainAuxArrays) THEN
01168                   mGlobal(1:nDim,1)=mGlobal(1:nDim,1) - mbc - 1
01169                   mGlobal(1:nDim,2)=mGlobal(1:nDim,2) + mbc + 1
01170                ELSE
01171                   mGlobal(1:nDim,1)=mGlobal(1:nDim,1) - mbc
01172                   mGlobal(1:nDim,2)=mGlobal(1:nDim,2) + mbc
01173                END IF
01174                Neighbors=Neighbors .OR. BoxOverlap(child%box%mGlobal, mGlobal)
01175             END DO
01176          END DO
01177       END DO
01178 
01179    END FUNCTION Neighbors
01180 
01182 
01183 
01186 
01191    FUNCTION GetChildMGlobal(node,mB,level)
01192       TYPE(NodeDef), POINTER :: node
01193       INTEGER, DIMENSION(3,2) :: GetChildMGlobal, mB
01194       INTEGER :: level
01195       GetChildmGlobal=1
01196 !      GetChildmGlobal(nDim+1:3,:)=1
01197       GetChildMGlobal(1:nDim,1)=(node%box%mGlobal(1:nDim,1)-2+mB(1:nDim,1))*levels(level)%CoarsenRatio+1
01198       GetChildMGlobal(1:nDim,2)=(node%box%mGlobal(1:nDim,1)-1+mB(1:nDim,2))*levels(level)%CoarsenRatio
01199    END FUNCTION GetChildMGlobal
01200 
01205    FUNCTION GetChildMBounds(node,childmGlobal,level)
01206       TYPE(NodeDef), POINTER :: node
01207       INTEGER, DIMENSION(3,2) :: GetChildMBounds, childmGlobal
01208       INTEGER :: level
01209       INTEGER, DIMENSION(3,2) :: work_array
01210 
01211 !      GetChildmBounds(nDim+1:3,:)=1
01212       GetChildmBounds=1
01213 !PRINT "('GetChildMBounds(', i1 ,')::ChildmGlobal = [', 6i4, '].')", MPI_id, ChildmGlobal
01214 !PRINT "('GetChildMBounds(', i1 ,')::node%box%mGlobal = [', 6i4, '].')", MPI_id, node%box%mGlobal
01215 !PRINT "('GetChildMBounds(', i1 ,')::CoarsenRatio = ', i4, '.')", MPI_id, levels(level)%CoarsenRatio
01216 !      work_array = 1
01217 !      work_array(1:nDim,1)=(childmGlobal(1:nDim,1)-1)/levels(level)%CoarsenRatio - node%box%mGlobal(1:nDim,1)+2
01218 !      work_array(1:nDim,2)=(childmGlobal(1:nDim,2)  )/levels(level)%CoarsenRatio - node%box%mGlobal(1:nDim,1)+1
01219 !PRINT "('GetChildMBounds(', i1 ,')::work_array = [', 6i4, '].')", MPI_id, work_array
01220 !PRINT *
01221 
01222       GetChildMBounds(1:nDim,1)=(childmGlobal(1:nDim,1)-1)/levels(level)%CoarsenRatio - node%box%mGlobal(1:nDim,1)+2
01223       GetChildMBounds(1:nDim,2)=(childmGlobal(1:nDim,2)  )/levels(level)%CoarsenRatio - node%box%mGlobal(1:nDim,1)+1
01224    END FUNCTION GetChildMBounds
01226 
01229 
01234    SUBROUTINE CreateNodeBox(box_array, node_box, proc_id)
01235 
01236       INTEGER, DIMENSION(3,2) :: box_array
01237       INTEGER, OPTIONAL :: proc_id
01238       TYPE(NodeBox), POINTER :: node_box
01239 
01240       INTEGER :: iErr
01241 
01242 
01243 !      IF (ASSOCIATED(node_box)) THEN
01244 !         DEALLOCATE(node_box)
01245 !         NULLIFY(node_box)
01246 !      END IF
01247 
01248       ALLOCATE(node_box, STAT=iErr)
01249 
01250       IF (iErr /= 0) THEN
01251          PRINT *, "AppendNodeBox() error:  unable to allocate new NodeBox structure."
01252          STOP
01253       END IF
01254 
01255       node_box%mGlobal = box_array
01256       IF (PRESENT(proc_id)) THEN
01257          node_box%MPI_id = proc_id
01258          !          PRINT "('CreateNodeBox([', 6i4, '], ', i2, ') done.')", box_array, proc_id
01259       ELSE
01260          !          PRINT "('CreateNodeBox([', 6i4, '], ', i2, ') done.')", box_array, MPI_id
01261       END IF
01262 
01263    END SUBROUTINE CreateNodeBox
01264 
01267    SUBROUTINE DestroyNodeBox(node_box)
01268 
01269       TYPE(NodeBox), POINTER :: node_box
01270 
01271       IF (ASSOCIATED(node_box)) THEN
01272          DEALLOCATE(node_box)
01273          NULLIFY(node_box)
01274       END IF
01275 
01276    END SUBROUTINE DestroyNodeBox
01277 
01281    FUNCTION MatchBox(box1,box2)
01282       TYPE(NodeBox) :: box1,box2
01283       LOGICAL :: MatchBox
01284       MatchBox=box1%MPI_ID==box2%MPI_ID .AND. ALL(box1%mGlobal==box2%mGlobal)
01285    END FUNCTION MatchBox
01286 
01288 
01291 
01296    SUBROUTINE AddNodeBoxToList(last_box, first_box, current_box)
01297 
01298       TYPE(NodeBoxList), POINTER :: last_box
01299       TYPE(NodeBoxList), POINTER, OPTIONAL :: first_box
01300       TYPE(NodeBoxList), POINTER, OPTIONAL :: current_box
01301 
01302       INTEGER :: iErr
01303 
01304 
01305       iErr = 0
01306 
01307       IF (.NOT. ASSOCIATED(last_box)) THEN
01308          ALLOCATE(last_box, STAT=iErr)
01309          IF (PRESENT(first_box))  first_box => last_box
01310          IF (PRESENT(current_box))  current_box => first_box
01311 
01312       ELSE
01313 
01314          IF (ASSOCIATED(last_box%next)) THEN
01315             PRINT *, "AddNodeBox error: last_box%next associated."
01316             STOP
01317          END IF
01318 
01319          ALLOCATE(last_box%next, STAT=iErr)
01320          last_box => last_box%next
01321 
01322       END IF
01323 
01324       IF (iErr /= 0) THEN
01325          PRINT *, "AddNodeBox() error: unable to allocate new node box."
01326          STOP
01327       END IF
01328 
01329       NULLIFY(last_box%next)
01330 
01331    END SUBROUTINE AddNodeBoxToList
01332 
01335    RECURSIVE SUBROUTINE ClearNodeBoxList(box_list)
01336 
01337       TYPE(NodeBoxList), POINTER :: box_list
01338 
01339 
01340       IF (.NOT. ASSOCIATED(box_list))  RETURN
01341 
01342       IF (ASSOCIATED(box_list%next))  CALL ClearNodeBoxList(box_list%next)
01343 
01344       DEALLOCATE(box_list)
01345       NULLIFY(box_list)
01346 
01347    END SUBROUTINE ClearNodeBoxList
01348 
01349    !! Counts the number of elements in a NodeBox list.
01350    !! @param box_list The list whose elements are to be counted.
01351    INTEGER FUNCTION NodeBoxCount(box_list)
01352 
01353       TYPE(NodeBoxList), POINTER :: box_list
01354 
01355       INTEGER :: work_count
01356       TYPE(NodeBoxList), POINTER :: iter
01357 
01358 
01359       iter => box_list
01360       work_count = 0
01361 
01362       DO WHILE (ASSOCIATED(iter))
01363          work_count = work_count + 1
01364          iter => iter%next
01365       END DO
01366 
01367       NodeBoxCount = work_count
01368 
01369    END FUNCTION NodeBoxCount
01370 
01375    SUBROUTINE FindParent(level, mGlobal, parent)
01376        INTEGER :: level
01377        INTEGER, DIMENSION(MAX_DIMS, 2) :: mGlobal
01378        TYPE(NodeDef), POINTER :: parent
01379 
01380        TYPE(NodeDefList), POINTER :: iter
01381 
01382 
01383        NULLIFY(parent)
01384 
01385        ! Search through local nodes for parent.
01386        iter => Nodes(level - 1)%p
01387 
01388        DO WHILE (ASSOCIATED(iter))
01389 
01390           IF (ALL(mGlobal(:,1) >= iter%self%box%mGlobal(:,1)) .AND. &
01391               ALL(mGlobal(:,2) <= iter%self%box%mGlobal(:,2))) THEN
01392               parent => iter%self
01393               RETURN
01394           END IF
01395 
01396        END DO
01397 
01398        ! search through external nodes for parent.
01399        iter => ExternalNodes(level - 1)%p
01400 
01401        DO WHILE (ASSOCIATED(iter))
01402 
01403           IF (ALL(mGlobal(:,1) >= iter%self%box%mGlobal(:,1)) .AND. &
01404               ALL(mGlobal(:,2) <= iter%self%box%mGlobal(:,2))) THEN
01405               parent => iter%self
01406               RETURN
01407           END IF
01408 
01409        END DO
01410 
01411    END SUBROUTINE FindParent
01412 
01418    SUBROUTINE StrictFindParent(level, mGlobal, parent, caller)
01419 
01420        INTEGER :: level
01421        INTEGER, DIMENSION(MAX_DIMS, 2) :: mGlobal
01422        TYPE(NodeDef), POINTER :: parent
01423        CHARACTER(LEN=*) :: caller
01424 
01425 
01426        ! Search for parent node.
01427        CALL FindParent(level, mGlobal, parent)
01428 
01429        ! Kill program with detailed error info if node is not found.
01430        IF (.NOT. ASSOCIATED(parent)) THEN
01431            PRINT *, caller, "::StrictFindNode() error: parent not found."
01432            PRINT "('StrictFindParent(level = ', i3, ') failed on box [', 6I5, ']')", level, mGlobal
01433            STOP
01434        END IF
01435 ! Comment test
01436    END SUBROUTINE StrictFindParent
01437 
01439 
01440 END MODULE TreeDeclarations
01441 
 All Classes Files Functions Variables