Merge branch 'issue/triaInitial' into 'development'

Issue with random position in volumes

See merge request JorgeGonz/fpakc!21
This commit is contained in:
Jorge Gonzalez 2021-05-24 10:41:28 +00:00
commit c6e3238810
7 changed files with 95 additions and 35 deletions

View file

@ -0,0 +1,64 @@
reset
set macros
term = 'postscript eps color enhanced'
font = 'font "CMUSerif-Roman,20"'
linew = 'lw 2.0'
times='"{/Symbol \264}"'
size_x=8.5
size_y=5.2
tmar=0.998
bmar=0.145
lmarLabel = 0.070
lmar=lmarLabel+0.085
rmar=0.97
space_x=(rmar-lmar)
xt_off_0=0.5
yt_off_0=0.6
set xtics nomirror offset 0.0,xt_off_0
set mxtics 2
set ytics nomirror offset yt_off_0,0.0
set mytics 2
set pointsize 1.5
set style line 1 pt 4 lc rgb "#B50427" #Squares red
set style line 2 pt 6 lc rgb "#3B4CC1" #Circles blue
set style line 3 pt 1 lc rgb "#2CA02C" #Crosses green
set style line 4 pt 2 lc rgb "#FE7F0E" #Exes orange
set style line 5 pt 8 lc rgb "#D6696B" #Triangles light red
set style line 10 lt 1 lw 2.0 lc rgb "black" #Black solid line
set terminal @term size size_x cm, size_y cm @linew @font
set output "comp_temp.eps"
#files
filename1 = "OUTPUT_Argon.dat"
filename2 = "OUTPUT_Argon+.dat"
set lmargin at screen lmar
set rmargin at screen rmar
set xrange [-0.01:1.01]
set xtics 0.2
set xlabel "Time (s)" offset 0.0,1.2
set yrange [250:3550]
set format y "%3.0f"
set ytics 500
set ylabel "Temperature (K)" offset 1.4,0.0
set key width 0.5 height 0.1 spacing 1.3 samplen 0.2 box opaque font ",16"
set key at graph 0.95, graph 0.9 right top
plot filename1 u ($1):($7) t "Ar" ls 1 with lines, \
filename2 u ($1):($7) t "Ar^{+}" ls 2 with lines, \
'< paste OUTPUT_Argon.dat OUTPUT_Argon+.dat' u ($1):($7 + $14) t "Sum" ls 3 with lines

View file

@ -622,7 +622,7 @@ MODULE moduleMesh2DCart
REAL(8), ALLOCATABLE:: fPsi(:)
xii(1) = random( 0.D0, 1.D0)
xii(2) = random( 0.D0, 1.D0)
xii(2) = random( 0.D0, 1.D0 - xii(1))
xii(3) = 0.D0
fPsi = self%fPsi(xii)

View file

@ -643,10 +643,9 @@ MODULE moduleMesh2DCyl
REAL(8), ALLOCATABLE:: fPsi(:)
xii(1) = random( 0.D0, 1.D0)
xii(2) = random( 0.D0, 1.D0)
xii(2) = random( 0.D0, 1.D0 - xii(1))
xii(3) = 0.D0
ALLOCATE(fPsi(1:3))
fPsi = self%fPsi(xii)
r(1) = DOT_PRODUCT(fPsi, self%z)

View file

@ -215,7 +215,9 @@ MODULE moduleMesh3DCart
REAL(8):: xii(1:3)
REAL(8):: fPsi(1:3)
xii = (/random(), random(), 0.D0 /)
xii(1) = random( 0.D0, 1.D0)
xii(2) = random( 0.D0, 1.D0 - xii(1))
xii(3) = 0.D0
fPsi = self%fPsi(xii)
r = (/DOT_PRODUCT(fPsi, self%x), &
@ -294,9 +296,9 @@ MODULE moduleMesh3DCart
REAL(8):: xii(1:3)
REAL(8), ALLOCATABLE:: fPsi(:)
xii(1) = random(0.D0, 1.D0)
xii(2) = random(0.D0, 1.D0)
xii(3) = random(0.D0, 1.D0)
xii(1) = random( 0.D0, 1.D0)
xii(2) = random( 0.D0, 1.D0 - xii(1))
xii(3) = random( 0.D0, 1.D0 - xii(1) - xii(2))
ALLOCATE(fPsi(1:4))
fPsi = self%fPsi(xii)

View file

@ -43,7 +43,7 @@ MODULE moduleMeshOutputGmsh2
WRITE(60, "(A)") '$EndNodeData'
WRITE(60, "(A)") '$NodeData'
WRITE(60, "(A)") '1'
WRITE(60, "(A)") '"' // species(i)%obj%name // ' velocity (m/s)"'
WRITE(60, "(A)") '"' // species(i)%obj%name // ' velocity (m s^-1)"'
WRITE(60, *) 1
WRITE(60, *) time
WRITE(60, *) 3
@ -189,7 +189,7 @@ MODULE moduleMeshOutputGmsh2
WRITE(20, "(A)") '$ElementData'
WRITE(20, "(A)") '1'
WRITE(20, "(A)") '"Electric Field (V/m)"'
WRITE(20, "(A)") '"Electric Field (V m^-1)"'
WRITE(20, *) 1
WRITE(20, *) time
WRITE(20, *) 3

View file

@ -186,15 +186,6 @@ MODULE moduleMesh
END SUBROUTINE initVol_interface
SUBROUTINE scatter_interface(self, part)
USE moduleSpecies
IMPORT:: meshVol
CLASS(meshVol), INTENT(in):: self
CLASS(particle), INTENT(in):: part
END SUBROUTINE scatter_interface
PURE FUNCTION gatherEF_interface(self, xi) RESULT(EF)
IMPORT:: meshVol
CLASS(meshVol), INTENT(in):: self
@ -574,6 +565,7 @@ MODULE moduleMesh
!TODO: try to combine this with the findCell for a regular mesh
!Find the volume in which particle reside in the mesh for collisions
!No boundary interaction taken into account
SUBROUTINE findCellCollMesh(part)
USE moduleSpecies
IMPLICIT NONE
@ -669,7 +661,7 @@ MODULE moduleMesh
DO e=1, self%numVols
vol => self%vols(e)%obj
nPart = vol%listPart_in%amount
!Computes iterations if there is more than one particle in the cell
!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

View file

@ -5,7 +5,7 @@ MODULE moduleCollisions
!Abstract type for collision between two particles
TYPE, ABSTRACT:: collisionBinary
REAL(8):: rMass !Reduced mass
REAL(8):: sMass !Summed mass
REAL(8):: sMassInv !Summed mass
TYPE(table1D):: crossSec !cross section of collision
CONTAINS
PROCEDURE(collideBinary_interface), PASS, DEFERRED:: collide
@ -169,8 +169,8 @@ MODULE moduleCollisions
CALL collision%crossSec%convert(eV2J/(m_ref*v_ref**2), 1.D0/L_ref**2)
!Calculates reduced mass
collision%sMass = mass_i+mass_j
collision%rMass = (mass_i*mass_j)/collision%sMass
collision%sMassInv = 1.D0/(mass_i+mass_j)
collision%rMass = (mass_i*mass_j)*collision%sMassInv
END SUBROUTINE initBinaryElastic
@ -180,6 +180,7 @@ MODULE moduleCollisions
USE moduleSpecies
USE moduleConstParam
USE moduleRandom
USE moduleMath
IMPLICIT NONE
CLASS(collisionBinaryElastic), INTENT(in):: self
@ -193,7 +194,7 @@ MODULE moduleCollisions
REAL(8), DIMENSION(1:3):: vCM
REAL(8):: vp(1:3)
vRel = SUM(DABS(part_i%v-part_j%v)) !TODO make function of norm1
vRel = NORM2(part_i%v-part_j%v)
eRel = self%rMass*vRel**2
sigmaVrel = self%crossSec%get(eRel)*vRel
sigmaVrelMaxNew = sigmaVrelMaxNew + sigmaVrel
@ -205,8 +206,8 @@ MODULE moduleCollisions
vp = vRel*randomDirectionVHS()
!Assign velocities to particles
part_i%v = vCM + m_j*vp/self%sMass
part_j%v = vCM - m_i*vp/self%sMass
part_i%v = vCM + m_j*vp*self%sMassInv
part_j%v = vCM - m_i*vp*self%sMassInv
END IF
@ -238,8 +239,8 @@ MODULE moduleCollisions
CALL collision%crossSec%convert(eV2J/(m_ref*v_ref**2), 1.D0/L_ref**2)
!Calculates reduced mass
collision%sMass = mass_i+mass_j
collision%rMass = (mass_i*mass_j)/collision%sMass
collision%sMassInv = 1.D0/(mass_i+mass_j)
collision%rMass = (mass_i*mass_j)*collision%sMassInv
!Specific parameters for ionization collision
SELECT TYPE(collision)
@ -269,6 +270,7 @@ MODULE moduleCollisions
USE moduleErrors
USE moduleList
USE moduleRandom
USE moduleMath
USE OMP_LIB
IMPLICIT NONE
@ -283,7 +285,7 @@ MODULE moduleCollisions
REAL(8), DIMENSION(1:3):: vp_e, vp_n
!eRel (in units of [m][L]^2[t]^-2
vRel = SUM(DABS(part_i%v-part_j%v)) !TODO make function of norm1
vRel = NORM2(part_i%v-part_j%v)
eRel = self%rMass*vRel**2
!Relative energy must be higher than threshold
IF (eRel > self%eThreshold) THEN
@ -370,8 +372,8 @@ MODULE moduleCollisions
CALL collision%crossSec%convert(eV2J/(m_ref*v_ref**2), 1.D0/L_ref**2)
!Calculates reduced mass
collision%sMass = mass_i+mass_j
collision%rMass = (mass_i*mass_j)/collision%sMass
collision%sMassInv = 1.D0/(mass_i+mass_j)
collision%rMass = (mass_i*mass_j)*collision%sMassInv
!Specific parameters for ionization collision
SELECT TYPE(collision)
@ -401,7 +403,7 @@ MODULE moduleCollisions
USE moduleErrors
USE moduleList
USE moduleRandom
USE OMP_LIB
USE moduleMath
IMPLICIT NONE
CLASS(collisionBinaryRecombination), INTENT(in):: self
@ -414,7 +416,7 @@ MODULE moduleCollisions
REAL(8), DIMENSION(1:3):: vp_i
!eRel (in units of [m][L]^2[t]^-2
vRel = SUM(DABS(part_i%v-part_j%v)) !TODO make function of norm1
vRel = NORM2(part_i%v-part_j%v)
eRel = self%rMass*vRel**2
!Relative energy must be higher than threshold
sigmaVrel = self%crossSec%get(eRel)*vRel
@ -476,8 +478,8 @@ MODULE moduleCollisions
CALL collision%crossSec%convert(eV2J/(m_ref*v_ref**2), 1.D0/L_ref**2)
!Calculates reduced mass
collision%sMass = mass_i+mass_j
collision%rMass = (mass_i*mass_j)/collision%sMass
collision%sMassInv = 1.D0/(mass_i+mass_j)
collision%rMass = (mass_i*mass_j)/collision%sMassInv
END SUBROUTINE initBinaryChargeExchange
@ -485,6 +487,7 @@ MODULE moduleCollisions
part_i, part_j)
USE moduleSpecies
USE moduleRandom
USE moduleMath
IMPLICIT NONE
CLASS(collisionBinaryChargeExchange), INTENT(in):: self
@ -496,7 +499,7 @@ MODULE moduleCollisions
REAL(8):: eRel !relative energy
!eRel (in units of [m][L]^2[t]^-2
vRel = SUM(DABS(part_i%v-part_j%v)) !TODO make function of norm1
vRel = NORM2(part_i%v-part_j%v)
eRel = self%rMass*vRel**2
sigmaVrel = self%crossSec%get(eRel)*vRel
sigmaVrelMaxNew = sigmaVrelMaxNew + sigmaVrel