Adding a time step for collisions
A new option has been added in which MCC are computed with its own time step. If no time is provided, then the minimum time step of the simulation is employed.
This commit is contained in:
parent
c6e3238810
commit
b6a7eb9ced
9 changed files with 88 additions and 53 deletions
|
|
@ -42,7 +42,7 @@ MODULE moduleMeshOutput0D
|
|||
USE moduleOutput
|
||||
IMPLICIT NONE
|
||||
|
||||
CLASS(meshGeneric), INTENT(in):: self
|
||||
CLASS(meshGeneric), INTENT(inout):: self
|
||||
INTEGER, INTENT(in):: t
|
||||
CHARACTER(:), ALLOCATABLE:: fileName
|
||||
|
||||
|
|
@ -56,7 +56,7 @@ MODULE moduleMeshOutput0D
|
|||
END IF
|
||||
|
||||
OPEN(20, file = path // folder // '/' // fileName, position = 'append', action = 'write')
|
||||
WRITE(20, "(ES20.6E3, I20)") REAL(t)*tauMin*ti_ref, mesh%vols(1)%obj%nColl
|
||||
WRITE(20, "(ES20.6E3, I20)") REAL(t)*tauMin*ti_ref, self%vols(1)%obj%nColl
|
||||
CLOSE(20)
|
||||
|
||||
END SUBROUTINE printColl0D
|
||||
|
|
|
|||
|
|
@ -95,7 +95,7 @@ MODULE moduleMeshOutputGmsh2
|
|||
USE moduleOutput
|
||||
IMPLICIT NONE
|
||||
|
||||
CLASS(meshGeneric), INTENT(in):: self
|
||||
CLASS(meshGeneric), INTENT(inout):: self
|
||||
INTEGER, INTENT(in):: t
|
||||
INTEGER:: numEdges
|
||||
INTEGER:: n
|
||||
|
|
|
|||
|
|
@ -310,7 +310,7 @@ MODULE moduleMesh
|
|||
SUBROUTINE printColl_interface(self, t)
|
||||
IMPORT meshGeneric
|
||||
|
||||
CLASS(meshGeneric), INTENT(in):: self
|
||||
CLASS(meshGeneric), INTENT(inout):: self
|
||||
INTEGER, INTENT(in):: t
|
||||
|
||||
END SUBROUTINE printColl_interface
|
||||
|
|
@ -637,7 +637,7 @@ MODULE moduleMesh
|
|||
END FUNCTION findCellBrute
|
||||
|
||||
!Computes collisions in element
|
||||
SUBROUTINE doCollisions(self)
|
||||
SUBROUTINE doCollisions(self, t)
|
||||
USE moduleCollisions
|
||||
USE moduleSpecies
|
||||
USE moduleList
|
||||
|
|
@ -646,6 +646,7 @@ MODULE moduleMesh
|
|||
IMPLICIT NONE
|
||||
|
||||
CLASS(meshGeneric), INTENT(inout), TARGET:: self
|
||||
INTEGER, INTENT(in):: t
|
||||
INTEGER:: e
|
||||
CLASS(meshVol), POINTER:: vol
|
||||
INTEGER:: nPart !Number of particles inside the cell
|
||||
|
|
@ -657,49 +658,57 @@ MODULE moduleMesh
|
|||
REAL(8):: sigmaVrelMaxNew
|
||||
TYPE(pointerArray), ALLOCATABLE:: partTemp(:)
|
||||
|
||||
!$OMP DO SCHEDULE(DYNAMIC)
|
||||
DO e=1, self%numVols
|
||||
vol => self%vols(e)%obj
|
||||
nPart = vol%listPart_in%amount
|
||||
!Calculates number of collisions if there is more than one particle in the cell
|
||||
IF (nPart > 1) THEN
|
||||
!Probability of collision
|
||||
pMax = vol%totalWeight*vol%sigmaVrelMax*tauMin/vol%volume
|
||||
IF (MOD(t, everyColl) == 0) THEN
|
||||
!Collisions need to be performed in this iteration
|
||||
!$OMP DO SCHEDULE(DYNAMIC)
|
||||
DO e=1, self%numVols
|
||||
vol => self%vols(e)%obj
|
||||
nPart = vol%listPart_in%amount
|
||||
|
||||
!Number of collisions in the cell
|
||||
vol%nColl = NINT(REAL(nPart)*pMax*0.5D0)
|
||||
!Resets the number of collisions
|
||||
vol%nColl = 0
|
||||
|
||||
IF (vol%nColl > 0) THEN
|
||||
!Converts the list of particles to an array for easy access
|
||||
partTemp = vol%listPart_in%convert2Array()
|
||||
!Calculates number of collisions if there is more than one particle in the cell
|
||||
IF (nPart > 1) THEN
|
||||
!Probability of collision
|
||||
pMax = vol%totalWeight*vol%sigmaVrelMax*tauColl/vol%volume
|
||||
|
||||
END IF
|
||||
!Number of collisions in the cell
|
||||
vol%nColl = NINT(REAL(nPart)*pMax*0.5D0)
|
||||
|
||||
DO n = 1, vol%nColl
|
||||
!Select random numbers
|
||||
rnd = random(1, nPart)
|
||||
part_i => partTemp(rnd)%part
|
||||
rnd = random(1, nPart)
|
||||
part_j => partTemp(rnd)%part
|
||||
ij = interactionIndex(part_i%species%n, part_j%species%n)
|
||||
sigmaVrelMaxNew = 0.D0
|
||||
DO k = 1, interactionMatrix(ij)%amount
|
||||
CALL interactionMatrix(ij)%collisions(k)%obj%collide(vol%sigmaVrelMax, sigmaVrelMaxNew, part_i, part_j)
|
||||
|
||||
END DO
|
||||
|
||||
!Update maximum cross section*v_rel per each collision
|
||||
IF (sigmaVrelMaxNew > vol%sigmaVrelMax) THEN
|
||||
vol%sigmaVrelMax = sigmaVrelMaxNew
|
||||
IF (vol%nColl > 0) THEN
|
||||
!Converts the list of particles to an array for easy access
|
||||
partTemp = vol%listPart_in%convert2Array()
|
||||
|
||||
END IF
|
||||
|
||||
END DO
|
||||
DO n = 1, vol%nColl
|
||||
!Select random numbers
|
||||
rnd = random(1, nPart)
|
||||
part_i => partTemp(rnd)%part
|
||||
rnd = random(1, nPart)
|
||||
part_j => partTemp(rnd)%part
|
||||
ij = interactionIndex(part_i%species%n, part_j%species%n)
|
||||
sigmaVrelMaxNew = 0.D0
|
||||
DO k = 1, interactionMatrix(ij)%amount
|
||||
CALL interactionMatrix(ij)%collisions(k)%obj%collide(vol%sigmaVrelMax, sigmaVrelMaxNew, part_i, part_j)
|
||||
|
||||
END IF
|
||||
END DO
|
||||
|
||||
END DO
|
||||
!$OMP END DO
|
||||
!Update maximum cross section*v_rel per each collision
|
||||
IF (sigmaVrelMaxNew > vol%sigmaVrelMax) THEN
|
||||
vol%sigmaVrelMax = sigmaVrelMaxNew
|
||||
|
||||
END IF
|
||||
|
||||
END DO
|
||||
|
||||
END IF
|
||||
|
||||
END DO
|
||||
!$OMP END DO
|
||||
|
||||
END IF
|
||||
|
||||
END SUBROUTINE doCollisions
|
||||
|
||||
|
|
|
|||
|
|
@ -4,6 +4,7 @@ MODULE moduleCaseParam
|
|||
INTEGER:: tFinal, tInitial = 0
|
||||
REAL(8), ALLOCATABLE:: tau(:)
|
||||
REAL(8):: tauMin
|
||||
REAL(8):: tauColl
|
||||
|
||||
END MODULE moduleCaseParam
|
||||
|
||||
|
|
|
|||
|
|
@ -2,6 +2,9 @@ MODULE moduleCollisions
|
|||
USE moduleSpecies
|
||||
USE moduleTable
|
||||
|
||||
!Integer for when collisions are computed
|
||||
INTEGER:: everyColl
|
||||
|
||||
!Abstract type for collision between two particles
|
||||
TYPE, ABSTRACT:: collisionBinary
|
||||
REAL(8):: rMass !Reduced mass
|
||||
|
|
|
|||
|
|
@ -40,6 +40,11 @@ MODULE moduleInput
|
|||
CALL readSpecies(config)
|
||||
CALL checkStatus(config, "readSpecies")
|
||||
|
||||
!Reads case parameters
|
||||
CALL verboseError('Reading Case parameters...')
|
||||
CALL readCase(config)
|
||||
CALL checkStatus(config, "readCase")
|
||||
|
||||
!Read interactions between species
|
||||
CALL verboseError('Reading interaction between species...')
|
||||
CALL readInteractions(config)
|
||||
|
|
@ -55,10 +60,10 @@ MODULE moduleInput
|
|||
CALL readGeometry(config)
|
||||
CALL checkStatus(config, "readGeometry")
|
||||
|
||||
!Reads case parameters
|
||||
CALL verboseError('Reading Case parameters...')
|
||||
CALL readCase(config)
|
||||
CALL checkStatus(config, "readCase")
|
||||
!Read initial state for species
|
||||
CALL verboseError('Reading Initial state...')
|
||||
CALL readInitial(config)
|
||||
CALL checkStatus(config, "readInitial")
|
||||
|
||||
!Read injection of particles
|
||||
CALL verboseError('Reading injection of particles from boundaries...')
|
||||
|
|
@ -233,11 +238,6 @@ MODULE moduleInput
|
|||
WRITE(tString, '(I1)') iterationDigits
|
||||
iterationFormat = "(I" // tString // "." // tString // ")"
|
||||
|
||||
!Read initial state for species
|
||||
CALL verboseError('Reading Initial state...')
|
||||
CALL readInitial(config)
|
||||
CALL checkStatus(config, "readInitial")
|
||||
|
||||
END SUBROUTINE readCase
|
||||
|
||||
!Reads the initial information for the species
|
||||
|
|
@ -558,6 +558,8 @@ MODULE moduleInput
|
|||
USE moduleCollisions
|
||||
USE moduleErrors
|
||||
USE moduleMesh
|
||||
USE moduleCaseParam
|
||||
USE moduleRefParam
|
||||
USE OMP_LIB
|
||||
USE json_module
|
||||
IMPLICIT NONE
|
||||
|
|
@ -595,6 +597,24 @@ MODULE moduleInput
|
|||
|
||||
END IF
|
||||
|
||||
!Reads collision time step
|
||||
CALL config%info('interactions.timeStep', found)
|
||||
IF (found) THEN
|
||||
CALL config%get('interactions.timeStep', tauColl, found)
|
||||
tauColl = tauColl / ti_ref
|
||||
|
||||
ELSE
|
||||
tauColl = tauMin
|
||||
|
||||
END IF
|
||||
|
||||
IF (tauColl < tauMin) THEN
|
||||
CALL criticalError('Collisional time step below minimum time step.', 'readInteractions')
|
||||
|
||||
END IF
|
||||
|
||||
everyColl = NINT(tauColl / tauMin)
|
||||
|
||||
!Inits the MCC matrix
|
||||
CALL initInteractionMatrix(interactionMatrix)
|
||||
|
||||
|
|
|
|||
Loading…
Add table
Add a link
Reference in a new issue