Merge branch 'master.dev' into 'master'

Release 3.3.1

See merge request piclas/piclas!972

.githooks/pre-commit

@@ -171,19 +171,19 @@ if [ -n "$CHANGED" ]; then
   fi
 
   if [ "$EXEERR" -ne 0 ]; then
-    printf "${RED}The following file(s) are executable, which is not allowed. Remove the execute permission via 'chmod -x filename' and try again.${NC}\n"
+    printf "${RED}The following file(s) are executable, which is not allowed. Remove the execute permission via 'chmod -x filename' to fix the files. Run 'git add' on the fixed files and then 'git commit' afterwards.${NC}\n"
     printf "$EXELINE"
   fi
 
   if [ "$TABERR" -ne 0 ]; then
-    printf "${RED}The following file(s) contain tabs, which is not allowed. Remove the tab stops and try again.${NC}\n"
+    printf "${RED}The following file(s) contain tabs (see https://en.wikipedia.org/wiki/Tab_key), which is not allowed. Remove the tab stops from all files. Run 'git add' on the fixed files and 'git commit' afterwards.${NC}\n"
     printf "$TABLINE"
   fi
 
   # Abort if hook encountered any error
   if [ "$ERROR" -ne 0 ]; then
-    echo '------------------------------------------------------------------------------------------------------------------------------------'
-    echo 'Commit rejected! You can override this check by passing "GIT_OVERRIDE_LIMITS=1" to git, e.g., run "GIT_OVERRIDE_LIMITS=1 git commit".'
+    echo "------------------------------------------------------------------------------------------------------------------------------------"
+    echo "Commit rejected! If you still want to commit your work by ignoring the errors, you can override this check by passing 'GIT_OVERRIDE_LIMITS=1' to git, e.g., run 'GIT_OVERRIDE_LIMITS=1 git commit'."
     exit 1
   fi
 

.github/workflows/cmake-ninja.yml

@@ -411,7 +411,7 @@ jobs:
     - name: Upload artifacts
       uses: actions/upload-artifact@v4
       with:
-        name: piclas-binaries-v3.3.0
+        name: piclas-binaries-v3.3.1
         path: artifacts
 
     - name: Upload release asset

docs/documentation/userguide/features-and-models/particle-initialization-and-emission.md

@@ -531,6 +531,13 @@ simulation is then determined by
 $$\dot{m} = \frac{QM}{1000RT},$$
 
 where $R=8.314$ J mol$^{-1}$K$^{-1}$ is the gas constant, $M$ the molar mass in [g mol$^{-1}$] and $T$ is the gas temperature [K].
+
+In some cases it might be useful to utilize averaged values across the complete BC by enabling
+
+    AdaptiveBC-AverageValuesOverBC = T
+
+Here, the cell-local values are averaged and the resulting macroscopic values are utilized for the particle emission.
+
 It should be noted that while multiple adaptive boundaries are possible, adjacent boundaries that share a mesh element should be avoided or treated carefully.
 Examples are given as part of the regression tests in `regressioncheck/CHE_DSMC/SurfFlux_Tria_Adaptive_ConstMassflow` and `SurfFlux_Tria_Adaptive_ConstPressure`.
 

regressioncheck/CHE_BGK/2D_VTS_Insert_CellLocal/parameter.ini

@@ -31,6 +31,7 @@ CFLscale   = 0.2  ! Scaling of theoretical CFL number
 DoParticleLatencyHiding=T
 DoLoadBalance=T
 DoInitialAutoRestart=T
+UseH5IOLoadBalance = T
 ! =============================================================================== !
 ! PARTICLES
 ! =============================================================================== !

regressioncheck/CHE_BGK/builds.ini

@@ -2,11 +2,10 @@
 binary=./bin/piclas
 
 ! fixed compiler flags
-CMAKE_BUILD_TYPE         = RELEASE
+CMAKE_BUILD_TYPE         = Debug
 LIBS_BUILD_HDF5          = OFF
 PICLAS_POLYNOMIAL_DEGREE = N
 PICLAS_EQNSYSNAME        = maxwell
 PICLAS_TIMEDISCMETHOD    = BGK-Flow
 LIBS_USE_MPI             = ON
 PICLAS_NODETYPE          = GAUSS
-PICLAS_LOADBALANCE       = OFF

src/globals/globals_vars.f90

@@ -24,7 +24,7 @@ MODULE MOD_Globals_Vars
 CHARACTER(LEN=6),PARAMETER :: ProgramName  = 'PICLas'              !> name of this program
 INTEGER,PARAMETER          :: MajorVersion = 3                     !> FileVersion number saved in each hdf5 file with hdf5 header
 INTEGER,PARAMETER          :: MinorVersion = 3                     !> FileVersion number saved in each hdf5 file with hdf5 header
-INTEGER,PARAMETER          :: PatchVersion = 0                     !> FileVersion number saved in each hdf5 file with hdf5 header
+INTEGER,PARAMETER          :: PatchVersion = 1                     !> FileVersion number saved in each hdf5 file with hdf5 header
 REAL,PARAMETER             :: FileVersionReal  = REAL(MajorVersion,8)+REAL(MinorVersion,8)/10.+REAL(PatchVersion,8)/100.
                                                                    !> OLD number saved in each hdf5 file with hdf5 header
 INTEGER,PARAMETER          :: FileVersionInt = PatchVersion+MinorVersion*100+MajorVersion*10000

src/particles/bgk/bgk_colloperator.f90

@@ -185,6 +185,12 @@ SUBROUTINE BGK_CollisionOperator(iPartIndx_Node, nPart, NodeVolume, AveragingVal
   IF(DSMC%CalcQualityFactors) THEN
     BGK_MaxRotRelaxFactor          = MAX(BGK_MaxRotRelaxFactor,MAXVAL(rotrelaxfreqSpec(:))*dtCell)
   END IF
+ELSE
+  rotrelaxfreqSpec(:) = 0.0
+  vibrelaxfreqSpec(:) = 0.0
+  IF(DSMC%CalcQualityFactors) THEN
+    BGK_MaxRotRelaxFactor          = 0.0
+  END IF
 END IF
 
 ! 4.) Determine the relaxation temperatures and the number of particles undergoing a relaxation (including vibration + rotation)
@@ -299,10 +305,12 @@ SUBROUTINE BGK_CollisionOperator(iPartIndx_Node, nPart, NodeVolume, AveragingVal
 DO iSpec = 1, nSpecies
   ! Calculate scaling factor alpha per species, see F. Hild, M. Pfeiffer, "Multi-species modeling in the particle-based ellipsoidal
   ! statistical Bhatnagar-Gross-Krook method including internal degrees of freedom", subitted to Phys. Fluids, August 2023
-  IF (NewEnRot(iSpec).GT.0.0) THEN
-    alphaRot(iSpec) = OldEn/NewEnRot(iSpec)*(Xi_RotSpec(iSpec)*RotRelaxWeightSpec(iSpec)/(Xi_RotTotal+3.*(totalWeight-1.)))
-  ELSE
-    alphaRot(iSpec) = 0.0
+  IF((Species(iSpec)%InterID.EQ.2).OR.(Species(iSpec)%InterID.EQ.20)) THEN
+    IF (NewEnRot(iSpec).GT.0.0) THEN
+      alphaRot(iSpec) = OldEn/NewEnRot(iSpec)*(Xi_RotSpec(iSpec)*RotRelaxWeightSpec(iSpec)/(Xi_RotTotal+3.*(totalWeight-1.)))
+    ELSE
+      alphaRot(iSpec) = 0.0
+    END IF
   END IF
 END DO
 DO iLoop = 1, nRotRelax

src/particles/bgk/bgk_init.f90

@@ -121,7 +121,7 @@ SUBROUTINE InitBGK()
 DO iSpec=1, nSpecies
   IF ((Species(iSpec)%InterID.EQ.2).OR.(Species(iSpec)%InterID.EQ.20)) MoleculePresent = .TRUE.
   ALLOCATE(SpecBGK(iSpec)%CollFreqPreFactor(nSpecies))
-  ! Calculation of the prefacor of the collision frequency per species
+  ! Calculation of the prefactor of the collision frequency per species
   ! S. Chapman and T.G. Cowling, "The mathematical Theory of Non-Uniform Gases", Cambridge University Press, 1970, S. 87f
   DO iSpec2=1, nSpecies
     SpecBGK(iSpec)%CollFreqPreFactor(iSpec2)= 4.*CollInf%dref(iSpec,iSpec2)**2.0 &
@@ -144,6 +144,9 @@ SUBROUTINE InitBGK()
   IF (DoVirtualCellMerge) THEN
     CALL abort(__STAMP__,'Virtual cell merge not implemented for coupled DSMC-BGK simulations!')
   END IF
+  IF(DSMC%RotRelaxProb.GT.1.0.OR.DSMC%VibRelaxProb.GT.1.0) THEN
+    CALL abort(__STAMP__,'Variable relaxation probabilities not implemented for coupled DSMC-BGK simulations!')
+  END IF
 #if USE_MPI
   IF (DoParticleLatencyHiding) THEN
     CALL abort(__STAMP__,'Particle latency hiding not implemented for coupled DSMC-BGK simulations!')

src/particles/dsmc/dsmc_relaxation.f90

@@ -311,7 +311,7 @@ SUBROUTINE InitCalcVibRelaxProb()
 ! Initialize the calculation of the variable vibrational relaxation probability in the cell for each iteration
 !===================================================================================================================================
 ! MODULES
-USE MOD_DSMC_Vars          ,ONLY: DSMC, VarVibRelaxProb 
+USE MOD_DSMC_Vars          ,ONLY: DSMC, VarVibRelaxProb
 USE MOD_Particle_Vars      ,ONLY: nSpecies
 
 ! IMPLICIT VARIABLE HANDLING
@@ -321,7 +321,7 @@ SUBROUTINE InitCalcVibRelaxProb()
 INTEGER                   :: iSpec
 !===================================================================================================================================
 
-IF(DSMC%VibRelaxProb.EQ.2.0) THEN ! Set summs for variable vibrational relaxation to zero
+IF(DSMC%VibRelaxProb.EQ.2.0) THEN ! Set sums for variable vibrational relaxation to zero
   DO iSpec=1,nSpecies
     VarVibRelaxProb%ProbVibAvNew(iSpec) = 0
     VarVibRelaxProb%nCollis(iSpec) = 0
@@ -336,13 +336,15 @@ SUBROUTINE DSMC_calc_P_rot(iSpec1, iSpec2, iPair, iPart, Xi_rel, ProbRot, ProbRo
 ! Calculation of probability for rotational relaxation. Different Models implemented:
 ! 0 - Constant Probability
 ! 1 - No rotational relaxation. RotRelaxProb = 0
-! 2 - Boyd
+! 2 - Boyd for N2
 ! 3 - Zhang (Nonequilibrium Direction Dependent)
 !===================================================================================================================================
 ! MODULES
-USE MOD_Globals            ,ONLY : Abort
-USE MOD_Globals_Vars       ,ONLY : Pi, BoltzmannConst
-USE MOD_DSMC_Vars          ,ONLY : SpecDSMC, Coll_pData, PartStateIntEn, DSMC, useRelaxProbCorrFactor, CollInf
+USE MOD_Globals            ,ONLY: Abort
+USE MOD_Globals_Vars       ,ONLY: Pi, BoltzmannConst
+USE MOD_Particle_Vars      ,ONLY: UseVarTimeStep, usevMPF
+USE MOD_part_tools         ,ONLY: GetParticleWeight
+USE MOD_DSMC_Vars          ,ONLY: SpecDSMC, Coll_pData, PartStateIntEn, DSMC, useRelaxProbCorrFactor, CollInf, RadialWeighting
 ! IMPLICIT VARIABLE HANDLING
 IMPLICIT NONE
 !-----------------------------------------------------------------------------------------------------------------------------------
@@ -357,8 +359,13 @@ SUBROUTINE DSMC_calc_P_rot(iSpec1, iSpec2, iPair, iPart, Xi_rel, ProbRot, ProbRo
 REAL                      :: TransEn, RotEn, RotDOF, CorrFact           ! CorrFact: To correct sample Bias
                                                                         ! (fewer DSMC particles than natural ones)
 !===================================================================================================================================
+! Note that during probability calculation, collision energy only contains translational part
+IF (usevMPF.OR.RadialWeighting%DoRadialWeighting.OR.UseVarTimeStep) THEN
+  TransEn = Coll_pData(iPair)%Ec / GetParticleWeight(iPart)
+ELSE
+  TransEn = Coll_pData(iPair)%Ec
+END IF
 
-TransEn    = Coll_pData(iPair)%Ec ! notice that during probability calculation,Collision energy only contains translational part
 RotDOF     = SpecDSMC(iSpec1)%Xi_Rot
 RotEn      = PartStateIntEn(2,iPart)
 ProbRot    = 0.
@@ -405,7 +412,7 @@ END SUBROUTINE DSMC_calc_P_rot
 
 SUBROUTINE DSMC_calc_P_elec(iSpec1, iSpec2, ProbElec)
 !===================================================================================================================================
-! Calculation of probability for electronic relaxation. 
+! Calculation of probability for electronic relaxation.
 !===================================================================================================================================
 ! MODULES
 USE MOD_DSMC_Vars          ,ONLY : SpecDSMC, useRelaxProbCorrFactor, DSMC
@@ -516,7 +523,7 @@ END SUBROUTINE DSMC_calc_P_vib
 SUBROUTINE DSMC_calc_var_P_vib(iSpec, jSpec, iPair, ProbVib)
 !===================================================================================================================================
 ! Calculation of probability for vibrational relaxation for variable relaxation rates. This has to average over all collisions!
-! No instantanious variable probability calculateable
+! No instantaneous variable probability calculateable
 !===================================================================================================================================
 ! MODULES
 USE MOD_Globals      ,ONLY: Abort
@@ -607,10 +614,10 @@ END SUBROUTINE SumVibRelaxProb
 
 SUBROUTINE FinalizeCalcVibRelaxProb(iElem)
 !===================================================================================================================================
-  ! Finalize the calculation of the variable vibrational relaxation probability in the cell for each iteration
+! Finalize the calculation of the variable vibrational relaxation probability in the cell for each iteration
 !===================================================================================================================================
 ! MODULES
-USE MOD_DSMC_Vars          ,ONLY: DSMC, VarVibRelaxProb 
+USE MOD_DSMC_Vars          ,ONLY: DSMC, VarVibRelaxProb
 USE MOD_Particle_Vars      ,ONLY: nSpecies
 
 ! IMPLICIT VARIABLE HANDLING

src/particles/emission/particle_surface_flux.f90

@@ -1205,6 +1205,8 @@ END SUBROUTINE CalcConstMassflowWeight
 
 !===================================================================================================================================
 !> SurfaceFlux: Determine the particle velocity of each inserted particle when inserted at a surface
+!> Based on: Garcia, A. L., & Wagner, W. (2006). Generation of the Maxwellian inflow distribution. Journal of Computational Physics,
+!>           217(2), 693–708. https://doi.org/10.1016/j.jcp.2006.01.025
 !===================================================================================================================================
 SUBROUTINE SetSurfacefluxVelocities(Mode,iSpec,iSF,iSample,jSample,iSide,BCSideID,SideID,NbrOfParticle,PartIns,particle_xis)
 ! MODULES

src/particles/restart/particle_restart.f90

@@ -626,14 +626,12 @@ SUBROUTINE ParticleRestart()
       CounterElec = 0
       CounterAmbi = 0
 
+      ! Increase the array size by TotalNbrOfMissingParticlesSum if needed
+      IF(PDM%ParticleVecLength+TotalNbrOfMissingParticlesSum.GT.PDM%maxParticleNumber) &
+        CALL IncreaseMaxParticleNumber(TotalNbrOfMissingParticlesSum)
+
       DO iPart = 1, TotalNbrOfMissingParticlesSum
         ! Sanity check
-        IF(CurrentPartNum.GT.PDM%maxParticleNumber)THEn
-          IPWRITE(UNIT_StdOut,'(I0,A,I0)') " CurrentPartNum        = ",  CurrentPartNum
-          IPWRITE(UNIT_StdOut,'(I0,A,I0)') " PDM%maxParticleNumber = ",  PDM%maxParticleNumber
-          CALL abort(__STAMP__,'Missing particle ID > PDM%maxParticleNumber. Increase Part-MaxParticleNumber!')
-        END IF !CurrentPartNum.GT.PDM%maxParticleNumber
-
         ! Do not search particles twice: Skip my own particles, because these have already been searched for before they are
         ! sent to all other procs
         ASSOCIATE( myFirst => OffsetTotalNbrOfMissingParticles(myRank) + 1 ,&
@@ -745,7 +743,6 @@ SUBROUTINE ParticleRestart()
         END IF
       END DO
 #endif
-      ! IF(PDM%ParticleVecLength.GT.PDM%maxParticleNumber) CALL IncreaseMaxParticleNumber(PDM%ParticleVecLength*CEILING(1+0.5*PDM%MaxPartNumIncrease)-PDM%maxParticleNumber)
 
       ! Combine number of found particles to make sure none are lost completely or found twice
       IF(MPIroot)THEN

src/particles/sampling/particle_sampling_vars.f90

@@ -27,7 +27,7 @@ MODULE MOD_Particle_Sampling_Vars
 !-----------------------------------------------------------------------------------------------------------------------------------
 ! Sampling of elements with a boundary for adaptive surface flux and porous BC
 LOGICAL                                 :: UseAdaptiveBC                  ! Flag is set if an adaptive boundary is present
-LOGICAL                                 :: AdaptBCAverageValBC            ! Flag to enable/disable averaging accross the whole BC
+LOGICAL                                 :: AdaptBCAverageValBC            ! Flag to enable/disable averaging across the whole BC
 REAL, ALLOCATABLE                       :: AdaptBCAverageMacroVal(:,:,:)  ! Macroscopic values averaged over BC
                                                                           ! (1:3, 1:nSpecies, 1:nSurfaceFluxBCs)
                                                                           !  1: Number density