JGonzalez/fpakc
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Common scatter subroutine

Browse source 
All subroutines of scattering particle properties to the nodes of a
volume have been converged into one in moduleMesh.
This commit is contained in:
Jorge Gonzalez 2021-04-21 23:40:58 +02:00
commit 35bd61fda9
9 changed files with 39 additions and 408 deletions
Download patch file Download diff file Expand all files Collapse all files

25
src/modules/mesh/0D/moduleMesh0D.f90
View file

@ -18,7 +18,6 @@ MODULE moduleMesh0D
PROCEDURE, PASS:: getNodes => getNodes0D
PROCEDURE, PASS:: randPos => randPos0D
PROCEDURE, NOPASS:: fPsi => fPsi0D
PROCEDURE, PASS:: scatter => scatter0D
PROCEDURE, PASS:: gatherEF => gatherEF0D
PROCEDURE, PASS:: elemK => elemK0D
PROCEDURE, PASS:: elemF => elemF0D
@ -113,30 +112,6 @@ MODULE moduleMesh0D
END FUNCTION fPsi0D
SUBROUTINE scatter0D(self, part)
USE moduleMath
USE moduleSpecies
USE OMP_LIB
IMPLICIT NONE
CLASS(meshVol0D), INTENT(in):: self
CLASS(particle), INTENT(in):: part
REAL(8):: tensorS(1:3,1:3)
CLASS(meshNode), POINTER:: node
INTEGER:: sp
tensorS = outerProduct(part%v, part%v)
node => self%n1
sp = part%species%n
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*tensorS
CALL OMP_UNSET_LOCK(node%lock)
END SUBROUTINE scatter0D
PURE FUNCTION gatherEF0D(self, xi) RESULT(EF)
IMPLICIT NONE

45
src/modules/mesh/1DCart/moduleMesh1DCart.f90
View file

@ -75,9 +75,7 @@ MODULE moduleMesh1DCart
PROCEDURE, PASS:: partialDer => partialDerSegm
PROCEDURE, PASS:: elemK => elemKSegm
PROCEDURE, PASS:: elemF => elemFSegm
PROCEDURE, NOPASS:: weight => weightSegm
PROCEDURE, NOPASS:: inside => insideSegm
PROCEDURE, PASS:: scatter => scatterSegm
PROCEDURE, PASS:: gatherEF => gatherEFSegm
PROCEDURE, PASS:: getNodes => getNodesSegm
PROCEDURE, PASS:: phy2log => phy2logSegm
@ -357,16 +355,6 @@ MODULE moduleMesh1DCart
END FUNCTION elemFSegm
PURE FUNCTION weightSegm(xi) RESULT(w)
IMPLICIT NONE
REAL(8), INTENT(in):: xi(1:3)
REAL(8):: w(1:2)
w = fPsiSegm(xi)
END FUNCTION weightSegm
PURE FUNCTION insideSegm(xi) RESULT(ins)
IMPLICIT NONE
@ -378,39 +366,6 @@ MODULE moduleMesh1DCart
END FUNCTION insideSegm
SUBROUTINE scatterSegm(self, part)
USE moduleMath
USE moduleSpecies
USE OMP_LIB
IMPLICIT NONE
CLASS(meshVol1DCartSegm), INTENT(in):: self
CLASS(particle), INTENT(in):: part
REAL(8):: w_p(1:2)
REAL(8):: tensorS(1:3,1:3)
CLASS(meshNode), POINTER:: node
INTEGER:: sp
w_p = self%weight(part%xi)
tensorS = outerProduct(part%v, part%v)
sp = part%species%n
node => self%n1
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(1)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(1)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(1)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n2
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(1)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(1)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(1)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
END SUBROUTINE scatterSegm
!Gathers EF at position Xii
PURE FUNCTION gatherEFSegm(self, xi) RESULT(EF)
IMPLICIT NONE

45
src/modules/mesh/1DRad/moduleMesh1DRad.f90
View file

@ -76,9 +76,7 @@ MODULE moduleMesh1DRad
PROCEDURE, PASS:: partialDer => partialDerRad
PROCEDURE, PASS:: elemK => elemKRad
PROCEDURE, PASS:: elemF => elemFRad
PROCEDURE, NOPASS:: weight => weightRad
PROCEDURE, NOPASS:: inside => insideRad
PROCEDURE, PASS:: scatter => scatterRad
PROCEDURE, PASS:: gatherEF => gatherEFRad
PROCEDURE, PASS:: getNodes => getNodesRad
PROCEDURE, PASS:: phy2log => phy2logRad
@ -369,16 +367,6 @@ MODULE moduleMesh1DRad
END FUNCTION elemFRad
PURE FUNCTION weightRad(xi) RESULT(w)
IMPLICIT NONE
REAL(8), INTENT(in):: xi(1:3)
REAL(8):: w(1:2)
w = fPsiRad(xi)
END FUNCTION weightRad
PURE FUNCTION insideRad(xi) RESULT(ins)
IMPLICIT NONE
@ -390,39 +378,6 @@ MODULE moduleMesh1DRad
END FUNCTION insideRad
SUBROUTINE scatterRad(self, part)
USE moduleMath
USE moduleSpecies
USE OMP_LIB
IMPLICIT NONE
CLASS(meshVol1DRadSegm), INTENT(in):: self
CLASS(particle), INTENT(in):: part
REAL(8):: w_p(1:2)
REAL(8):: tensorS(1:3,1:3)
CLASS(meshNode), POINTER:: node
INTEGER:: sp
w_p = self%weight(part%xi)
tensorS = outerProduct(part%v, part%v)
sp = part%species%n
node => self%n1
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(1)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(1)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(1)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n2
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(2)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(2)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(2)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
END SUBROUTINE scatterRad
!Gathers EF at position Xii
PURE FUNCTION gatherEFRad(self, xi) RESULT(EF)
IMPLICIT NONE

115
src/modules/mesh/2DCart/moduleMesh2DCart.f90
View file

@ -83,9 +83,7 @@ MODULE moduleMesh2DCart
PROCEDURE, PASS:: partialDer => partialDerQuad
PROCEDURE, PASS:: elemK => elemKQuad
PROCEDURE, PASS:: elemF => elemFQuad
PROCEDURE, NOPASS:: weight => weightQuad
PROCEDURE, NOPASS:: inside => insideQuad
PROCEDURE, PASS:: scatter => scatterQuad
PROCEDURE, PASS:: gatherEF => gatherEFQuad
PROCEDURE, PASS:: getNodes => getNodesQuad
PROCEDURE, PASS:: phy2log => phy2logQuad
@ -114,9 +112,7 @@ MODULE moduleMesh2DCart
PROCEDURE, PASS:: partialDer => partialDerTria
PROCEDURE, PASS:: elemK => elemKTria
PROCEDURE, PASS:: elemF => elemFTria
PROCEDURE, NOPASS:: weight => weightTria
PROCEDURE, NOPASS:: inside => insideTria
PROCEDURE, PASS:: scatter => scatterTria
PROCEDURE, PASS:: gatherEF => gatherEFTria
PROCEDURE, PASS:: getNodes => getNodesTria
PROCEDURE, PASS:: phy2log => phy2logTria
@ -470,17 +466,6 @@ MODULE moduleMesh2DCart
END FUNCTION elemFQuad
!Computes weights in the element nodes
PURE FUNCTION weightQuad(xi) RESULT(w)
IMPLICIT NONE
REAL(8), INTENT(in):: xi(1:3)
REAL(8):: w(1:4)
w = fPsiQuad(xi)
END FUNCTION weightQuad
!Checks if a particle is inside a quad element
PURE FUNCTION insideQuad(xi) RESULT(ins)
IMPLICIT NONE
@ -493,54 +478,6 @@ MODULE moduleMesh2DCart
END FUNCTION insideQuad
!Scatter properties of particle into element nodes
SUBROUTINE scatterQuad(self, part)
USE moduleMath
USE moduleSpecies
USE OMP_LIB
IMPLICIT NONE
CLASS(meshVol2DCartQuad), INTENT(in):: self
CLASS(particle), INTENT(in):: part
REAL(8):: w_p(1:4)
REAL(8):: tensorS(1:3,1:3)
CLASS(meshNode), POINTER:: node
INTEGER:: sp
w_p = self%weight(part%xi)
tensorS = outerProduct(part%v, part%v)
sp = part%species%n
node => self%n1
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(1)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(1)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(1)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n2
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(2)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(2)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(2)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n3
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(3)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(3)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(3)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n4
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(4)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(4)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(4)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
END SUBROUTINE scatterQuad
!Gathers the electric field at position xi
PURE FUNCTION gatherEFQuad(self,xi) RESULT(EF)
IMPLICIT NONE
@ -841,17 +778,6 @@ MODULE moduleMesh2DCart
END FUNCTION elemFTria
!Computes weights in the element nodes
PURE FUNCTION weightTria(xi) RESULT(w)
IMPLICIT NONE
REAL(8), INTENT(in):: xi(1:3)
REAL(8), ALLOCATABLE:: w(:)
w = fPsiTria(xi)
END FUNCTION weightTria
PURE FUNCTION insideTria(xi) RESULT(ins)
IMPLICIT NONE
@ -864,47 +790,6 @@ MODULE moduleMesh2DCart
END FUNCTION insideTria
!Scatter properties of particles into element
SUBROUTINE scatterTria(self, part)
USE moduleMath
USE moduleSpecies
USE OMP_LIB
IMPLICIT NONE
CLASS(meshVol2DCartTria), INTENT(in):: self
CLASS(particle), INTENT(in):: part
REAL(8):: w_p(1:3)
REAL(8):: tensorS(1:3,1:3)
CLASS(meshNode), POINTER:: node
INTEGER:: sp
w_p = self%weight(part%xi)
tensorS = outerProduct(part%v, part%v)
sp = part%species%n
node => self%n1
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(1)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(1)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(1)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n2
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(2)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(2)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(2)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n3
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(3)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(3)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(3)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
END SUBROUTINE scatterTria
!Gathers the electric field at position xi
PURE FUNCTION gatherEFTria(self,xi) RESULT(EF)
IMPLICIT NONE

115
src/modules/mesh/2DCyl/moduleMesh2DCyl.f90
View file

@ -84,9 +84,7 @@ MODULE moduleMesh2DCyl
PROCEDURE, PASS:: partialDer => partialDerQuad
PROCEDURE, PASS:: elemK => elemKQuad
PROCEDURE, PASS:: elemF => elemFQuad
PROCEDURE, NOPASS:: weight => weightQuad
PROCEDURE, NOPASS:: inside => insideQuad
PROCEDURE, PASS:: scatter => scatterQuad
PROCEDURE, PASS:: gatherEF => gatherEFQuad
PROCEDURE, PASS:: getNodes => getNodesQuad
PROCEDURE, PASS:: phy2log => phy2logQuad
@ -115,9 +113,7 @@ MODULE moduleMesh2DCyl
PROCEDURE, PASS:: partialDer => partialDerTria
PROCEDURE, PASS:: elemK => elemKTria
PROCEDURE, PASS:: elemF => elemFTria
PROCEDURE, NOPASS:: weight => weightTria
PROCEDURE, NOPASS:: inside => insideTria
PROCEDURE, PASS:: scatter => scatterTria
PROCEDURE, PASS:: gatherEF => gatherEFTria
PROCEDURE, PASS:: getNodes => getNodesTria
PROCEDURE, PASS:: phy2log => phy2logTria
@ -491,17 +487,6 @@ MODULE moduleMesh2DCyl
END FUNCTION elemFQuad
!Computes weights in the element nodes
PURE FUNCTION weightQuad(xi) RESULT(w)
IMPLICIT NONE
REAL(8), INTENT(in):: xi(1:3)
REAL(8):: w(1:4)
w = fPsiQuad(xi)
END FUNCTION weightQuad
!Checks if a particle is inside a quad element
PURE FUNCTION insideQuad(xi) RESULT(ins)
IMPLICIT NONE
@ -514,54 +499,6 @@ MODULE moduleMesh2DCyl
END FUNCTION insideQuad
!Scatter properties of particle into element nodes
SUBROUTINE scatterQuad(self, part)
USE moduleMath
USE moduleSpecies
USE OMP_LIB
IMPLICIT NONE
CLASS(meshVol2DCylQuad), INTENT(in):: self
CLASS(particle), INTENT(in):: part
REAL(8):: w_p(1:4)
REAL(8):: tensorS(1:3,1:3)
CLASS(meshNode), POINTER:: node
INTEGER:: sp
w_p = self%weight(part%xi)
tensorS = outerProduct(part%v, part%v)
sp = part%species%n
node => self%n1
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(1)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(1)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(1)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n2
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(2)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(2)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(2)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n3
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(3)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(3)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(3)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n4
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(4)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(4)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(4)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
END SUBROUTINE scatterQuad
!Gathers the electric field at position xi
PURE FUNCTION gatherEFQuad(self,xi) RESULT(EF)
IMPLICIT NONE
@ -871,17 +808,6 @@ MODULE moduleMesh2DCyl
END FUNCTION elemFTria
!Computes weights in the element nodes
PURE FUNCTION weightTria(xi) RESULT(w)
IMPLICIT NONE
REAL(8), INTENT(in):: xi(1:3)
REAL(8), ALLOCATABLE:: w(:)
w = fPsiTria(xi)
END FUNCTION weightTria
PURE FUNCTION insideTria(xi) RESULT(ins)
IMPLICIT NONE
@ -894,47 +820,6 @@ MODULE moduleMesh2DCyl
END FUNCTION insideTria
!Scatter properties of particles into element
SUBROUTINE scatterTria(self, part)
USE moduleMath
USE moduleSpecies
USE OMP_LIB
IMPLICIT NONE
CLASS(meshVol2DCylTria), INTENT(in):: self
CLASS(particle), INTENT(in):: part
REAL(8):: w_p(1:3)
REAL(8):: tensorS(1:3,1:3)
CLASS(meshNode), POINTER:: node
INTEGER:: sp
w_p = self%weight(part%xi)
tensorS = outerProduct(part%v, part%v)
sp = part%species%n
node => self%n1
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(1)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(1)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(1)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n2
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(2)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(2)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(2)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n3
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(3)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(3)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(3)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
END SUBROUTINE scatterTria
!Gathers the electric field at position xi
PURE FUNCTION gatherEFTria(self,xi) RESULT(EF)
IMPLICIT NONE

58
src/modules/mesh/3DCart/moduleMesh3DCart.f90
View file

@ -76,9 +76,7 @@ MODULE moduleMesh3DCart
PROCEDURE, PASS:: partialDer => partialDerTetra
PROCEDURE, PASS:: elemK => elemKTetra
PROCEDURE, PASS:: elemF => elemFTetra
PROCEDURE, NOPASS:: weight => weightTetra
PROCEDURE, NOPASS:: inside => insideTetra
PROCEDURE, PASS:: scatter => scatterTetra
PROCEDURE, PASS:: gatherEF => gatherEFTetra
PROCEDURE, PASS:: getNodes => getNodesTetra
PROCEDURE, PASS:: phy2log => phy2logTetra
@ -459,15 +457,6 @@ MODULE moduleMesh3DCart
END FUNCTION elemFTetra
PURE FUNCTION weightTetra(xii) RESULT(w)
IMPLICIT NONE
REAL(8), INTENT(in):: xii(1:3)
REAL(8), ALLOCATABLE:: w(:)
w = fPsiTetra(xii)
END FUNCTION weightTetra
PURE FUNCTION insideTetra(xi) RESULT(ins)
IMPLICIT NONE
@ -481,53 +470,6 @@ MODULE moduleMesh3DCart
END FUNCTION insideTetra
SUBROUTINE scatterTetra(self, part)
USE moduleMath
USE moduleSpecies
USE OMP_LIB
IMPLICIT NONE
CLASS(meshVol3DCartTetra), INTENT(in):: self
CLASS(particle), INTENT(in):: part
REAL(8):: w_p(1:4)
REAL(8):: tensorS(1:3, 1:3)
CLASS(meshNode), POINTER:: node
INTEGER:: sp
w_p = self%weight(part%xi)
tensorS = outerProduct(part%v, part%v)
sp = part%species%n
node => self%n1
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(1)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(1)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(1)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n2
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(2)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(2)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(2)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n3
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(3)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(3)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(3)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
node => self%n4
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*w_p(4)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*w_p(4)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*w_p(4)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
END SUBROUTINE scatterTetra
PURE FUNCTION gatherEFTetra(self, xi) RESULT(EF)
IMPLICIT NONE

36
src/modules/mesh/moduleMesh.f90
View file

@ -164,7 +164,7 @@ MODULE moduleMesh
PROCEDURE(getNodesVol_interface), DEFERRED, PASS:: getNodes
PROCEDURE(randPosVol_interface), DEFERRED, PASS:: randPos
PROCEDURE(fPsi_interface), DEFERRED, NOPASS:: fPsi
PROCEDURE(scatter_interface), DEFERRED, PASS:: scatter
PROCEDURE, PASS:: scatter
PROCEDURE(gatherEF_interface), DEFERRED, PASS:: gatherEF
PROCEDURE(elemK_interface), DEFERRED, PASS:: elemK
PROCEDURE(elemF_interface), DEFERRED, PASS:: elemF
@ -468,6 +468,40 @@ MODULE moduleMesh
END SUBROUTINE resetOutput
!Scatters particle properties into vol nodes
SUBROUTINE scatter(self, part)
USE moduleMath
USE moduleSpecies
USE OMP_LIB
IMPLICIT NONE
CLASS(meshVol), INTENT(inout):: self
CLASS(particle), INTENT(in):: part
REAL(8), ALLOCATABLE:: fPsi(:)
INTEGER, ALLOCATABLE:: volNodes(:)
REAL(8):: tensorS(1:3, 1:3)
INTEGER:: sp
INTEGER:: i, nNodes
CLASS(meshNode), POINTER:: node
fPsi = self%fPsi(part%xi)
tensorS = outerProduct(part%v, part%v)
sp = part%species%n
volNodes = self%getNodes()
nNodes = SIZE(volNodes)
DO i = 1, nNodes
node => mesh%nodes(volNodes(i))%obj
CALL OMP_SET_LOCK(node%lock)
node%output(sp)%den = node%output(sp)%den + part%weight*fPsi(i)
node%output(sp)%mom(:) = node%output(sp)%mom(:) + part%weight*fPsi(i)*part%v(:)
node%output(sp)%tensorS(:,:) = node%output(sp)%tensorS(:,:) + part%weight*fPsi(i)*tensorS
CALL OMP_UNSET_LOCK(node%lock)
END DO
END SUBROUTINE scatter
!Find next cell for particle
RECURSIVE SUBROUTINE findCell(self, part, oldCell)
USE moduleSpecies

2
src/modules/moduleInject.f90
View file

@ -127,7 +127,7 @@ MODULE moduleInject
et = et + 1
self%edges(et) = mesh%edges(e)%obj%n
!Assign connectivity between injection edge and meshColl volume
IF (ASSOCIATED(meshForMCC, meshColl)) THEN
IF (doubleMesh) THEN
nVolColl = findCellBrute(meshColl, mesh%edges(e)%obj%randPos())
IF (nVolColl > 0) THEN
mesh%edges(e)%obj%eColl => meshColl%vols(nVolColl)%obj

6
src/modules/moduleInput.f90
View file

@ -345,11 +345,11 @@ MODULE moduleInput
partNew%v(2) = velocityXi(2) + vTh*randomMaxwellian()
partNew%v(3) = velocityXi(3) + vTh*randomMaxwellian()
partNew%vol = e
IF (ASSOCIATED(meshForMCC, mesh)) THEN
partNew%volColl = partNew%vol
IF (doubleMesh) THEN
partNew%volColl = findCellBrute(meshColl, partNew%r)
ELSE
partNew%volColl = findCellBrute(meshColl, partNew%r)
partNew%volColl = partNew%vol
END IF
partNew%n_in = .TRUE.