Merge branch 'improvement.rot.periodic.restart' into 'master.dev'

[improvement.rot.periodic.restart] Improved initialization performance for rotational periodic BC by splitting...

See merge request piclas/piclas!683

regressioncheck/NIG_DSMC/RotPeriodicBC/parameter.ini

@@ -6,9 +6,11 @@ IniExactFunc  = 0
 ! =============================================================================== !
 ! DISCRETIZATION
 ! =============================================================================== !
-N             = 1  ! Polynomial degree
-NAnalyze      = 1  ! Number of analyze points
-NVisu      = 0
+N                    = 1  ! Polynomial degree
+NAnalyze             = 1  ! Number of analyze points
+NVisu                = 1
+VisuParticles        = T
+!DisplayLostParticles = T
 ! =============================================================================== !
 ! MESH
 ! =============================================================================== !

regressioncheck/NIG_DSMC/RotPeriodicBC/pre-hopr-mortar/hopr_mortar.ini

@@ -1,14 +1,15 @@
 
-DEFVAR=(INT):    i01 = 6   ! no. elems in left and right block
-DEFVAR=(INT):    i02 = 6   ! no. elems in upper block (should be twice the value of i01)
+DEFVAR=(INT):    i01 = 6!2   ! no. elems in left and right block
+DEFVAR=(INT):    i02 = 6!6!2   ! no. elems in upper block (should be twice the value of i01)
 
-DEFVAR=(INT):    ir1 = 25   ! no. elems in r for first ring
+DEFVAR=(INT):    ir1 =  25!3   ! no. elems in r for first ring
 DEFVAR=(REAL):   r01 = 3.5 ! middle square dim
 DEFVAR=(REAL):   r02 = 7.0 ! middle square dim
 DEFVAR=(REAL):   s0  = 0.2857142857142857 ! middle square dim
 
 
-DEFVAR=(INT):    iz = 50    !
+DEFVAR=(INT):    iz1 = 50!3    !
+DEFVAR=(INT):    iz2 = 100!6    !
 DEFVAR=(REAL):   lz = 1    ! length of domain in z
 
 DEFVAR=(REAL):   f1 = 1.    ! stretching factor in first ring
@@ -33,7 +34,7 @@ nZones = 2                           ! number of boxes
 
 !right-lower (x+)
 Corner       =(/r01 , 0.  , -lz   ,,   r02 , 0.  , -lz   ,,   r02 , r02 , -lz   ,,   r01 , r01 , -lz  ,,   r01 , 0.  , lz   ,,   r02 , 0.  , lz   ,,   r02 , r02 , lz   ,,   r01 , r01 , lz /)
-nElems       =(/ir1,i01,iz/)                   ! number of elements in each direction
+nElems       =(/ir1,i01,iz1/)                   ! number of elements in each direction
 BCIndex      =(/1  , 5  , 4  , 0  , 3  , 2/)   ! Indices of Boundary Conditions for  six Boundary Faces (z- , y- , x+ , y+ , x- , z+)
 
 !            =(/z- , y- , x+ , y+ , x- , z+/)  ! Indices of Boundary Conditions
@@ -42,7 +43,7 @@ factor       =(/f1,1.,1./)                     ! stretching
 
 !right-upper (y+)
 Corner       =(/0.  , r01 , -lz   ,,   r01 , r01 , -lz   ,,   r02 , r02 , -lz   ,,   0.  , r02 , -lz   ,,   0.  , r01 , lz   ,,   r01 , r01 , lz   ,,   r02 , r02 , lz   ,,   0.  , r02 , lz /)
-nElems       =(/i02,50,100/)                   ! number of elements in each direction
+nElems       =(/i02,iz1,iz2/)                   ! number of elements in each direction
 BCIndex      =(/1  , 3  , 0  , 4  , 6  , 2/)   ! Indices of Boundary Conditions for  six Boundary Faces (z- , y- , x+ , y+ , x- , z+)
 !            =(/z- , y- , x+ , y+ , x- , z+/)  ! Indices of Boundary Conditions
 elemtype     =108                              ! element type (108: Hexahedral)

regressioncheck/NIG_DSMC/RotPeriodicBC/readme.md

@@ -7,4 +7,20 @@
   finding a corresponding rotationally periodic side)
   - Debugging information should be written to TestRotatingWall_LostRotPeriodicSides_000.00000000000000000.h5 file
     (cannot be tested automatically unfortunately) when setting CalcMeshInfo=T
+- for debugging, use 2 procs and smaller mesh via
+
+        DEFVAR=(INT):    i01 = 2   ! no. elems in left and right block
+        DEFVAR=(INT):    i02 = 2   ! no. elems in upper block (should be twice the value of i01)
+
+        DEFVAR=(INT):    ir1 = 3   ! no. elems in r for first ring
+        DEFVAR=(REAL):   r01 = 3.5 ! middle square dim
+        DEFVAR=(REAL):   r02 = 7.0 ! middle square dim
+        DEFVAR=(REAL):   s0  = 0.2857142857142857 ! middle square dim
+
+
+        DEFVAR=(INT):    iz1 = 3
+        DEFVAR=(INT):    iz2 = 6
+        DEFVAR=(REAL):   lz = 1    ! length of domain in z
+
+        DEFVAR=(REAL):   f1 = 1.    ! stretching factor in first ring
 

src/io_hdf5/hdf5_output_elemdata.f90

@@ -28,10 +28,10 @@ MODULE MOD_HDF5_Output_ElemData
 
 #if USE_MPI
 PUBLIC :: WriteMyInvisibleRankToHDF5
+#endif /*USE_MPI*/
 #if defined(PARTICLES)
 PUBLIC :: WriteLostRotPeriodicSidesToHDF5
 #endif /*defined(PARTICLES)*/
-#endif /*USE_MPI*/
 
 PUBLIC :: WriteAdditionalElemData
 !===================================================================================================================================
@@ -191,6 +191,7 @@ SUBROUTINE WriteMyInvisibleRankToHDF5()
 WRITE(UNIT_stdOut,'(a)',ADVANCE='YES')'DONE'
 #endif
 END SUBROUTINE WriteMyInvisibleRankToHDF5
+#endif /*USE_MPI*/
 
 
 #if defined(PARTICLES)
@@ -247,7 +248,6 @@ SUBROUTINE WriteLostRotPeriodicSidesToHDF5()
 #endif
 END SUBROUTINE WriteLostRotPeriodicSidesToHDF5
 #endif /*defined(PARTICLES)*/
-#endif /*USE_MPI*/
 
 
 END MODULE MOD_HDF5_Output_ElemData

src/mesh/mesh_vars.f90

@@ -174,10 +174,11 @@ MODULE MOD_Mesh_Vars
                                                                      !< [1:4] (mortar) neighbors
                                                                      !< [1:6] local sides
                                                                      !< [OffSetElem+1:OffsetElem+PP_nElems]
-INTEGER(KIND=8),ALLOCATABLE     ::  ElemGlobalID(:)                  !< global element id of each element
-INTEGER(KIND=8),ALLOCATABLE     ::  myInvisibleRank(:)               !< global proc ID which the current proc cannot see (particle
+INTEGER(KIND=8),ALLOCATABLE     :: ElemGlobalID(:)                   !< global element id of each element
+INTEGER(KIND=8),ALLOCATABLE     :: myInvisibleRank(:)                !< global proc ID which the current proc cannot see (particle
                                                                      !< communication)
-INTEGER(KIND=8),ALLOCATABLE     ::  LostRotPeriodicSides(:)          !< Number of lost sides during rotational periodic search
+INTEGER(KIND=8),ALLOCATABLE     :: LostRotPeriodicSides(:)           !< Number of lost sides during rotational periodic search
+!LOGICAL                         :: RotPeriodicReBuild                !< Force re-building of mapping (might already exist)
 !-----------------------------------------------------------------------------------------------------------------------------------
 CHARACTER(LEN=255),ALLOCATABLE   :: BoundaryName(:)
 CHARACTER(LEN=255)               :: MeshFile        !< name of hdf5 meshfile (write with ending .h5!)

src/particles/boundary/particle_boundary_condition.f90

@@ -398,7 +398,7 @@ SUBROUTINE RotPeriodicBC(PartID,SideID,ElemID)
 USE MOD_Particle_Mesh_Tools    ,ONLY: ParticleInsideQuad3D
 USE MOD_part_operations        ,ONLY: RemoveParticle
 USE MOD_part_tools             ,ONLY: StoreLostParticleProperties
-USE MOD_Particle_Tracking_Vars ,ONLY: NbrOfLostParticles, TrackInfo, CountNbrOfLostParts
+USE MOD_Particle_Tracking_Vars ,ONLY: NbrOfLostParticles, TrackInfo, CountNbrOfLostParts,DisplayLostParticles
 USE MOD_DSMC_Vars              ,ONLY: DSMC, AmbipolElecVelo
 #ifdef CODE_ANALYZE
 USE MOD_Particle_Tracking_Vars ,ONLY: PartOut,MPIRankOut
@@ -413,7 +413,7 @@ SUBROUTINE RotPeriodicBC(PartID,SideID,ElemID)
 INTEGER,INTENT(INOUT),OPTIONAL    :: ElemID
 !-----------------------------------------------------------------------------------------------------------------------------------
 ! LOCAL VARIABLES
-INTEGER                              :: SideID2, ElemID2, iNeigh, RotSideID
+INTEGER                              :: iNeigh, RotSideID, GlobalElemID
 REAL                                 :: adaptTimeStep
 LOGICAL                              :: FoundInElem
 REAL                                 :: LastPartPos_old(1:3),Velo_old(1:3), Velo_oldAmbi(1:3)
@@ -480,7 +480,7 @@ SUBROUTINE RotPeriodicBC(PartID,SideID,ElemID)
       END IF
   END SELECT
 END ASSOCIATE
-! (2) update particle positon after periodic BC
+! (2) update particle position after periodic BC
 PartState(1:3,PartID)   = LastPartPos(1:3,PartID) + (1.0 - TrackInfo%alpha/TrackInfo%lengthPartTrajectory) * dt*RKdtFrac &
                         * PartState(4:6,PartID) * adaptTimeStep
 ! compute moved particle || rest of movement
@@ -496,26 +496,31 @@ SUBROUTINE RotPeriodicBC(PartID,SideID,ElemID)
   IF(PartID.EQ.PARTOUT)THEN
     IPWRITE(UNIT_stdout,'(I0,A)') '     RotPeriodicBC: '
     IPWRITE(UNIT_stdout,'(I0,A,3(1X,G0))') ' ParticlePosition-pp: ',PartState(1:3,PartID)
-    IPWRITE(UNIT_stdout,'(I0,A,3(1X,G0))') ' LastPartPo-pp:       ',LastPartPos(1:3,PartID)
+    IPWRITE(UNIT_stdout,'(I0,A,3(1X,G0))') ' LastPartPosition-pp: ',LastPartPos(1:3,PartID)
   END IF
 END IF
 #endif /*CODE_ANALYZE*/
 ! (3) move particle from old element to new element
-RotSideID = SurfSide2RotPeriodicSide((GlobalSide2SurfSide(SURF_SIDEID,SideID)))
+RotSideID = SurfSide2RotPeriodicSide(GlobalSide2SurfSide(SURF_SIDEID,SideID))
 DO iNeigh=1,NumRotPeriodicNeigh(RotSideID)
-  SideID2 = RotPeriodicSideMapping(RotSideID,iNeigh)
-  IF(SideID2.EQ.-1) THEN
-    CALL abort(__STAMP__,' ERROR: Halo-rot-periodic side has no corresponding side.')
-  END IF
-  ElemID2 = SideInfo_Shared(SIDE_ELEMID,SideID2)
-  ! find rotational periodic SideID2 through localization in all potentional rotational periodic sides
-  CALL ParticleInsideQuad3D(LastPartPos(1:3,PartID),ElemID2,FoundInElem)
+  ! find rotational periodic elem through localization in all potential rotational periodic sides
+  GlobalElemID = RotPeriodicSideMapping(RotSideID,iNeigh)
+  IF(GlobalElemID.EQ.-1) CALL abort(__STAMP__,' ERROR: Halo-rot-periodic side has no corresponding element.')
+  CALL ParticleInsideQuad3D(LastPartPos(1:3,PartID),GlobalElemID,FoundInElem)
   IF(FoundInElem) THEN
-    ElemID = ElemID2
+    ElemID = GlobalElemID
     EXIT
   END IF
 END DO
 IF(.NOT.FoundInElem) THEN
+  ! Particle appears to have not crossed any of the checked sides. Deleted!
+  IF(DisplayLostParticles)THEN
+    IPWRITE(*,*) 'Error in Particle TriaTracking! Particle Number',PartID,'lost. Element:', ElemID,'(species:',PartSpecies(PartID),')'
+    IPWRITE(*,*) 'LastPos: ', LastPartPos(1:3,PartID)
+    IPWRITE(*,*) 'Pos:     ', PartState(1:3,PartID)
+    IPWRITE(*,*) 'Velo:    ', PartState(4:6,PartID)
+    IPWRITE(*,*) 'Particle deleted!'
+  END IF ! DisplayLostParticles
   IF(CountNbrOfLostParts)THEN
     CALL StoreLostParticleProperties(PartID,ElemID)
     NbrOfLostParticles=NbrOfLostParticles+1