Merge branch 'fix.SEE.particle.insert' into 'master.dev'

fixed SEE particles that are inserted into dielectric region when being...

See merge request piclas/piclas!623

src/particles/surfacemodel/surfacemodel_main.f90

@@ -109,8 +109,8 @@ SUBROUTINE SurfaceModel(PartID,SideID,GlobalElemID,n_Loc)
 ProductSpecNbr = 0
 
 ! Store info of impacting particle for possible surface charging
+PartPosImpact(1:3) = LastPartPos(1:3,PartID)+TrackInfo%PartTrajectory(1:3)*TrackInfo%alpha
 IF(DoDielectricSurfaceCharge.AND.PartBound%Dielectric(iBC)) THEN ! Surface charging active + dielectric surface contact
-  PartPosImpact(1:3) = LastPartPos(1:3,PartID)+TrackInfo%PartTrajectory(1:3)*TrackInfo%alpha
   IF(usevMPF)THEN
     MPF = PartMPF(PartID)
   ELSE
@@ -181,7 +181,7 @@ SUBROUTINE SurfaceModel(PartID,SideID,GlobalElemID,n_Loc)
   END IF
   ! Emit the secondary electrons
   IF (ProductSpec(2).GT.0) THEN
-    CALL SurfaceModel_ParticleEmission(n_loc, PartID, SideID, ProductSpec, ProductSpecNbr, TempErgy, GlobalElemID)
+    CALL SurfaceModel_ParticleEmission(n_loc, PartID, SideID, ProductSpec, ProductSpecNbr, TempErgy, GlobalElemID,PartPosImpact(1:3))
     ! Deposit opposite charge of SEE on node
     IF(DoDielectricSurfaceCharge.AND.PartBound%Dielectric(iBC)) THEN
       ! Get MPF

src/particles/surfacemodel/surfacemodel_tools.f90

@@ -33,7 +33,7 @@ MODULE MOD_SurfaceModel_Tools
 !===================================================================================================================================
 !> Routine for the particle emission at a surface
 !===================================================================================================================================
-SUBROUTINE SurfaceModel_ParticleEmission(n_loc, PartID, SideID, ProductSpec, ProductSpecNbr, TempErgy, GlobElemID)
+SUBROUTINE SurfaceModel_ParticleEmission(n_loc, PartID, SideID, ProductSpec, ProductSpecNbr, TempErgy, GlobElemID,POI_vec)
 ! MODULES
 USE MOD_Globals                   ,ONLY: OrthoNormVec
 USE MOD_Part_Tools                ,ONLY: VeloFromDistribution
@@ -47,11 +47,13 @@ SUBROUTINE SurfaceModel_ParticleEmission(n_loc, PartID, SideID, ProductSpec, Pro
 USE MOD_Particle_Tracking_Vars    ,ONLY: TrackInfo
 USE MOD_Particle_Vars             ,ONLY: usevMPF,PartMPF
 USE MOD_part_tools                ,ONLY: CalcRadWeightMPF
+USE MOD_Particle_Mesh_Vars        ,ONLY: BoundsOfElem_Shared
 ! IMPLICIT VARIABLE HANDLING
 IMPLICIT NONE
 !-----------------------------------------------------------------------------------------------------------------------------------
 ! INPUT VARIABLES
 REAL,INTENT(IN)    :: n_loc(1:3)       !< normal vector of the surface
+REAL,INTENT(IN)    :: POI_vec(1:3)     !< Point Of Intersection
 INTEGER,INTENT(IN) :: PartID, SideID   !< Particle index and side index
 INTEGER,INTENT(IN) :: GlobElemID       !< global element ID of the impacting particle (used for creating a new particle)
 INTEGER,INTENT(IN) :: ProductSpec(2)   !< 1: product species of incident particle (also used for simple reflection)
@@ -66,36 +68,43 @@ SUBROUTINE SurfaceModel_ParticleEmission(n_loc, PartID, SideID, ProductSpec, Pro
 !-----------------------------------------------------------------------------------------------------------------------------------
 ! LOCAL VARIABLES
 INTEGER            :: iNewPart, NewPartID, locBCID, SurfSideID
-REAL               :: tang1(1:3), tang2(1:3), WallVelo(1:3), WallTemp, NewVelo(3), POI_vec(3), OldMPF
+REAL               :: tang1(1:3), tang2(1:3), WallVelo(1:3), WallTemp, NewVelo(3), OldMPF, BoundsOfElemCenter(1:3),NewPos(1:3)
+REAL,PARAMETER     :: eps=1e-6
+REAL,PARAMETER     :: eps2=1.0-eps
 !===================================================================================================================================
 locBCID    = PartBound%MapToPartBC(SideInfo_Shared(SIDE_BCID,SideID))
 SurfSideID = GlobalSide2SurfSide(SURF_SIDEID,SideID)
 WallTemp   = PartBound%WallTemp(locBCID)
 WallVelo   = PartBound%WallVelo(1:3,locBCID)
 
 IF(PartBound%RotVelo(locBCID)) THEN
-  POI_vec(1:3) = LastPartPos(1:3,PartID) + TrackInfo%PartTrajectory(1:3)*TrackInfo%alpha
   CALL CalcRotWallVelo(locBCID,POI_vec,WallVelo)
 END IF
 
 CALL OrthoNormVec(n_loc,tang1,tang2)
 
+! Get Elem Center
+BoundsOfElemCenter(1:3) = (/SUM(BoundsOfElem_Shared(1:2,1,GlobElemID)), &
+                            SUM(BoundsOfElem_Shared(1:2,2,GlobElemID)), &
+                            SUM(BoundsOfElem_Shared(1:2,3,GlobElemID)) /) / 2.
+
 ! Create new particles
 DO iNewPart = 1, ProductSpecNbr
   ! create new particle and assign correct energies
   ! sample newly created velocity
   NewVelo(1:3) = VeloFromDistribution(SurfModEnergyDistribution(locBCID),TempErgy,iNewPart,ProductSpecNbr)
   ! Rotate velocity vector from global coordinate system into the surface local coordinates (important: n_loc points outwards)
   NewVelo(1:3) = tang1(1:3)*NewVelo(1) + tang2(1:3)*NewVelo(2) - n_Loc(1:3)*NewVelo(3) + WallVelo(1:3)
-  ! Create new particle and get a free particle index
+  ! Create new position by using POI and moving the particle by eps in the direction of the element center
+  NewPos(1:3) = eps*BoundsOfElemCenter(1:3) + eps2*POI_vec(1:3)
   IF(usevMPF)THEN
     ! Get MPF of old particle
     OldMPF = PartMPF(PartID)
     ! New particle acquires the MPF of the impacting particle (not necessarily the MPF of the newly created particle species)
-    CALL CreateParticle(ProductSpec(2),LastPartPos(1:3,PartID),GlobElemID,NewVelo(1:3),0.,0.,0.,NewPartID=NewPartID, NewMPF=OldMPF)
+    CALL CreateParticle(ProductSpec(2),NewPos(1:3),GlobElemID,NewVelo(1:3),0.,0.,0.,NewPartID=NewPartID, NewMPF=OldMPF)
   ELSE
     ! New particle acquires the MPF of the new particle species
-    CALL CreateParticle(ProductSpec(2),LastPartPos(1:3,PartID),GlobElemID,NewVelo(1:3),0.,0.,0.,NewPartID=NewPartID)
+    CALL CreateParticle(ProductSpec(2),NewPos(1:3),GlobElemID,NewVelo(1:3),0.,0.,0.,NewPartID=NewPartID)
   END IF ! usevMPF
   ! Adding the energy that is transferred from the surface onto the internal energies of the particle
   CALL SurfaceModel_EnergyAccommodation(NewPartID,locBCID,WallTemp)