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

Modification in boundary conditions:

Browse source 
- Now absorption scatte the particle  properties into the edge nodes.
  - New boundary condition 'transparent' subsitute old absorption.
This commit is contained in:
Jorge Gonzalez 2020-12-17 21:16:35 +01:00
commit 9e3a1a771b
15 changed files with 212 additions and 148 deletions
Download patch file Download diff file Expand all files Collapse all files

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

@ -47,6 +47,15 @@ MODULE moduleMesh1DCart
END SUBROUTINE absorption
MODULE SUBROUTINE transparent(edge, part)
USE moduleSpecies
IMPLICIT NONE
CLASS(meshEdge), INTENT(inout):: edge
CLASS(particle), INTENT(inout):: part
END SUBROUTINE transparent
END INTERFACE
TYPE, PUBLIC, ABSTRACT, EXTENDS(meshVol):: meshVol1DCart
@ -105,7 +114,6 @@ MODULE moduleMesh1DCart
PROCEDURE, PASS:: getNodes => getNodesSegm
PROCEDURE, PASS:: phy2log => phy2logSegm
PROCEDURE, PASS:: nextElement => nextElementSegm
PROCEDURE, PASS:: resetOutput => resetOutputSegm
END TYPE meshVol1DCartSegm
@ -178,6 +186,9 @@ MODULE moduleMesh1DCart
TYPE IS(boundaryReflection)
self%fBoundary(s)%apply => reflection
TYPE IS(boundaryTransparent)
self%fBoundary(s)%apply => transparent
CLASS DEFAULT
CALL criticalError("Boundary type not defined in this geometry", 'initEdge1DCart')
@ -455,29 +466,6 @@ MODULE moduleMesh1DCart
END SUBROUTINE nextElementSegm
!Reset the output of nodes in element
PURE SUBROUTINE resetOutputSegm(self)
USE moduleSpecies
USE moduleOutput
IMPLICIT NONE
CLASS(meshVol1DCartSegm), INTENT(inout):: self
INTEGER:: k
DO k = 1, nSpecies
self%n1%output(k)%den = 0.D0
self%n1%output(k)%mom = 0.D0
self%n1%output(k)%tensorS = 0.D0
self%n2%output(k)%den = 0.D0
self%n2%output(k)%mom = 0.D0
self%n2%output(k)%tensorS = 0.D0
END DO
END SUBROUTINE resetOutputSegm
!COMMON FUNCTIONS FOR 1D VOLUME ELEMENTS
!Computes the element Jacobian determinant
PURE FUNCTION detJ1DCart(self, xi, dPsi_in) RESULT(dJ)

34
src/modules/mesh/1DCart/moduleMesh1DCartBoundary.f90
View file

@ -24,9 +24,43 @@ SUBMODULE (moduleMesh1DCart) moduleMesh1DCartBoundary
CLASS(meshEdge), INTENT(inout):: edge
CLASS(particle), INTENT(inout):: part
REAL(8):: rEdge(1) !Position of particle in the edge
REAL(8):: d !Distance from particle to edge
SELECT TYPE(edge)
TYPE IS(meshEdge1DCart)
rEdge(1) = edge%x
d = DABS(part%r(1) - rEdge(1))
IF (d > 0.D0) THEN
part%weight = part%weight / d
part%r(1) = rEdge(1)
END IF
IF (ASSOCIATED(edge%e1)) THEN
CALL edge%e1%scatter(part)
ELSE
CALL edge%e2%scatter(part)
END IF
END SELECT
part%n_in = .FALSE.
END SUBROUTINE absorption
SUBROUTINE transparent(edge, part)
USE moduleSpecies
IMPLICIT NONE
CLASS(meshEdge), INTENT(inout):: edge
CLASS(particle), INTENT(inout):: part
part%n_in = .FALSE.
END SUBROUTINE transparent
END SUBMODULE moduleMesh1DCartBoundary

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

@ -47,6 +47,15 @@ MODULE moduleMesh1DRad
END SUBROUTINE absorption
MODULE SUBROUTINE transparent(edge, part)
USE moduleSpecies
IMPLICIT NONE
CLASS(meshEdge), INTENT(inout):: edge
CLASS(particle), INTENT(inout):: part
END SUBROUTINE transparent
END INTERFACE
TYPE, PUBLIC, ABSTRACT, EXTENDS(meshVol):: meshVol1DRad
@ -106,7 +115,6 @@ MODULE moduleMesh1DRad
PROCEDURE, PASS:: getNodes => getNodesRad
PROCEDURE, PASS:: phy2log => phy2logRad
PROCEDURE, PASS:: nextElement => nextElementRad
PROCEDURE, PASS:: resetOutput => resetOutputRad
END TYPE meshVol1DRadSegm
@ -179,6 +187,9 @@ MODULE moduleMesh1DRad
TYPE IS(boundaryReflection)
self%fBoundary(s)%apply => reflection
TYPE IS(boundaryTransparent)
self%fBoundary(s)%apply => transparent
CLASS DEFAULT
CALL criticalError("Boundary type not defined in this geometry", 'initEdge1DRad')
@ -467,29 +478,6 @@ MODULE moduleMesh1DRad
END SUBROUTINE nextElementRad
!Reset the output of nodes in element
PURE SUBROUTINE resetOutputRad(self)
USE moduleSpecies
USE moduleOutput
IMPLICIT NONE
CLASS(meshVol1DRadSegm), INTENT(inout):: self
INTEGER:: k
DO k = 1, nSpecies
self%n1%output(k)%den = 0.D0
self%n1%output(k)%mom = 0.D0
self%n1%output(k)%tensorS = 0.D0
self%n2%output(k)%den = 0.D0
self%n2%output(k)%mom = 0.D0
self%n2%output(k)%tensorS = 0.D0
END DO
END SUBROUTINE resetOutputRad
!COMMON FUNCTIONS FOR 1D VOLUME ELEMENTS
!Computes the element Jacobian determinant
PURE FUNCTION detJ1DRad(self, xi, dPsi_in) RESULT(dJ)

36
src/modules/mesh/1DRad/moduleMesh1DRadBoundary.f90
View file

@ -16,6 +16,8 @@ SUBMODULE (moduleMesh1DRad) moduleMesh1DRadBoundary
END SELECT
part%n_in = .TRUE.
END SUBROUTINE reflection
SUBROUTINE absorption(edge, part)
@ -24,9 +26,43 @@ SUBMODULE (moduleMesh1DRad) moduleMesh1DRadBoundary
CLASS(meshEdge), INTENT(inout):: edge
CLASS(particle), INTENT(inout):: part
REAL(8):: rEdge(1) !Position of particle in the edge
REAL(8):: d !Distance from particle to edge
SELECT TYPE(edge)
TYPE IS(meshEdge1DRad)
rEdge(1) = edge%r
d = DABS(part%r(1) - rEdge(1))
IF (d > 0.D0) THEN
part%weight = part%weight / d
part%r(1) = rEdge(1)
END IF
IF (ASSOCIATED(edge%e1)) THEN
CALL edge%e1%scatter(part)
ELSE
CALL edge%e2%scatter(part)
END IF
END SELECT
part%n_in = .FALSE.
END SUBROUTINE absorption
SUBROUTINE transparent(edge, part)
USE moduleSpecies
IMPLICIT NONE
CLASS(meshEdge), INTENT(inout):: edge
CLASS(particle), INTENT(inout):: part
part%n_in = .FALSE.
END SUBROUTINE transparent
END SUBMODULE moduleMesh1DRadBoundary

71
src/modules/mesh/2DCyl/moduleMeshCyl.f90
View file

@ -55,6 +55,15 @@ MODULE moduleMeshCyl
END SUBROUTINE absorption
MODULE SUBROUTINE transparent(edge, part)
USE moduleSpecies
IMPLICIT NONE
CLASS(meshEdge), INTENT(inout):: edge
CLASS(particle), INTENT(inout):: part
END SUBROUTINE transparent
MODULE SUBROUTINE symmetryAxis(edge, part)
USE moduleSpecies
IMPLICIT NONE
@ -126,7 +135,6 @@ MODULE moduleMeshCyl
PROCEDURE, PASS:: getNodes => getNodesQuad
PROCEDURE, PASS:: phy2log => phy2logQuad
PROCEDURE, PASS:: nextElement => nextElementQuad
PROCEDURE, PASS:: resetOutput => resetOutputQuad
END TYPE meshVolCylQuad
@ -159,7 +167,6 @@ MODULE moduleMeshCyl
PROCEDURE, PASS:: getNodes => getNodesTria
PROCEDURE, PASS:: phy2log => phy2logTria
PROCEDURE, PASS:: nextElement => nextElementTria
PROCEDURE, PASS:: resetOutput => resetOutputTria
END TYPE meshVolCylTria
@ -237,6 +244,9 @@ MODULE moduleMeshCyl
TYPE IS(boundaryReflection)
self%fBoundary(s)%apply => reflection
TYPE IS(boundaryTransparent)
self%fBoundary(s)%apply => transparent
TYPE IS(boundaryAxis)
self%fBoundary(s)%apply => symmetryAxis
@ -636,36 +646,6 @@ MODULE moduleMeshCyl
END SUBROUTINE nextElementQuad
!Reset the output of nodes in quad element
PURE SUBROUTINE resetOutputQuad(self)
USE moduleSpecies
USE moduleOutput
IMPLICIT NONE
CLASS(meshVolCylQuad), INTENT(inout):: self
INTEGER:: k
DO k = 1, nSpecies
self%n1%output(k)%den = 0.D0
self%n1%output(k)%mom = 0.D0
self%n1%output(k)%tensorS = 0.D0
self%n2%output(k)%den = 0.D0
self%n2%output(k)%mom = 0.D0
self%n2%output(k)%tensorS = 0.D0
self%n3%output(k)%den = 0.D0
self%n3%output(k)%mom = 0.D0
self%n3%output(k)%tensorS = 0.D0
self%n4%output(k)%den = 0.D0
self%n4%output(k)%mom = 0.D0
self%n4%output(k)%tensorS = 0.D0
END DO
END SUBROUTINE resetOutputQuad
!TRIA ELEMENT
!Init tria element
SUBROUTINE initVolTriaCyl(self, n, p)
@ -999,33 +979,6 @@ MODULE moduleMeshCyl
END SUBROUTINE nextElementTria
!Reset the output of nodes in tria element
PURE SUBROUTINE resetOutputTria(self)
USE moduleSpecies
USE moduleOutput
IMPLICIT NONE
CLASS(meshVolCylTria), INTENT(inout):: self
INTEGER:: k
DO k = 1, nSpecies
self%n1%output(k)%den = 0.D0
self%n1%output(k)%mom = 0.D0
self%n1%output(k)%tensorS = 0.D0
self%n2%output(k)%den = 0.D0
self%n2%output(k)%mom = 0.D0
self%n2%output(k)%tensorS = 0.D0
self%n3%output(k)%den = 0.D0
self%n3%output(k)%mom = 0.D0
self%n3%output(k)%tensorS = 0.D0
END DO
END SUBROUTINE resetOutputTria
!COMMON FUNCTIONS FOR CYLINDRICAL VOLUME ELEMENTS
!Computes element Jacobian determinant
PURE FUNCTION detJCyl(self, xi, dPsi_in) RESULT(dJ)

53
src/modules/mesh/2DCyl/moduleMeshCylBoundary.f90
View file

@ -34,10 +34,10 @@ SUBMODULE (moduleMeshCyl) moduleMeshCylBoundary
part%v(1) = cosT*vpp(1) + sinT*vpp(2)
part%v(2) = -sinT*vpp(1) + cosT*vpp(2)
part%n_in = .TRUE.
END SELECT
part%n_in = .TRUE.
END SUBROUTINE reflection
!Absoption in a surface
@ -47,13 +47,60 @@ SUBMODULE (moduleMeshCyl) moduleMeshCylBoundary
CLASS(meshEdge), INTENT(inout):: edge
CLASS(particle), INTENT(inout):: part
REAL(8):: rEdge(1:2) !Position of particle projected to the edge
REAL(8):: a, b, c
REAL(8):: a2b2
REAL(8):: d !Distance from particle to edge
SELECT TYPE(edge)
TYPE IS(meshEdgeCyl)
a = (edge%z(1) - edge%z(2))
b = (edge%r(1) - edge%r(2))
c = edge%z(1)*edge%r(2) - edge%z(2)*edge%r(1)
!TODO: Add scatter to mesh nodes
a2b2 = a**2 + b**2
rEdge(1) = (b*( b*part%r(1) - a*part%r(2)) - a*c)/a2b2
rEdge(2) = (a*(-b*part%r(1) + a*part%r(2)) - b*c)/a2b2
d = NORM2(rEdge - part%r(1:2))
!Reduce weight of particle by the distance to the edge and move it to the edge
IF (d > 0.D0) THEN
part%weight = part%weight / d
part%r(1:2) = rEdge
END IF
!Scatter particle in associated volume
IF (ASSOCIATED(edge%e1)) THEN
CALL edge%e1%scatter(part)
ELSE
CALL edge%e2%scatter(part)
END IF
END SELECT
!Remove particle from the domain
part%n_in = .FALSE.
END SUBROUTINE absorption
!Transparent boundary condition
SUBROUTINE transparent(edge, part)
USE moduleSpecies
IMPLICIT NONE
CLASS(meshEdge), INTENT(inout):: edge
CLASS(particle), INTENT(inout):: part
!Removes particle from domain
part%n_in = .FALSE.
END SUBROUTINE transparent
!Symmetry axis. Dummy function
SUBROUTINE symmetryAxis(edge, part)
USE moduleSpecies
IMPLICIT NONE

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

@ -15,8 +15,9 @@ MODULE moduleMesh
TYPE(outputNode), ALLOCATABLE:: output(:)
TYPE(emNode):: emData
CONTAINS
PROCEDURE(initNode_interface), DEFERRED, PASS:: init
PROCEDURE(getCoord_interface), DEFERRED, PASS:: getCoordinates
PROCEDURE(initNode_interface), DEFERRED, PASS:: init
PROCEDURE(getCoord_interface), DEFERRED, PASS:: getCoordinates
PROCEDURE, PASS:: resetOutput
END TYPE meshNode
@ -146,7 +147,6 @@ MODULE moduleMesh
PROCEDURE(inside_interface), DEFERRED, NOPASS:: inside
PROCEDURE(nextElement_interface), DEFERRED, PASS:: nextElement
PROCEDURE, PASS:: collision
PROCEDURE(resetOutput_interface), DEFERRED, PASS:: resetOutput
END TYPE meshVol
@ -215,18 +215,6 @@ MODULE moduleMesh
END FUNCTION inside_interface
SUBROUTINE collision_interface(self)
IMPORT:: meshVol
CLASS(meshVol), INTENT(inout):: self
END SUBROUTINE collision_interface
PURE SUBROUTINE resetOutput_interface(self)
IMPORT:: meshVol
CLASS(meshVol), INTENT(inout):: self
END SUBROUTINE resetOutput_interface
END INTERFACE
!Containers for volumes in the mesh
@ -309,6 +297,24 @@ MODULE moduleMesh
CLASS(meshGeneric), ALLOCATABLE, TARGET:: mesh
CONTAINS
!Reset the output of node
PURE SUBROUTINE resetOutput(self)
USE moduleSpecies
USE moduleOutput
IMPLICIT NONE
CLASS(meshNode), INTENT(inout):: self
INTEGER:: k
DO k = 1, nSpecies
self%output(k)%den = 0.D0
self%output(k)%mom = 0.D0
self%output(k)%tensorS = 0.D0
END DO
END SUBROUTINE resetOutput
!Find next cell for particle
RECURSIVE SUBROUTINE findCell(self, part, oldCell)
USE moduleSpecies
@ -422,9 +428,6 @@ MODULE moduleMesh
self%totalWeight = 0.D0
!Reset output in nodes
CALL self%resetOutput()
!Erase the list of particles inside the cell
CALL self%listPart_in%erase()

6
src/modules/moduleInput.f90
View file

@ -369,9 +369,15 @@ MODULE moduleInput
CASE('absorption')
ALLOCATE(boundaryAbsorption:: boundary(i)%bTypes(s)%obj)
CASE('transparent')
ALLOCATE(boundaryTransparent:: boundary(i)%bTypes(s)%obj)
CASE('axis')
ALLOCATE(boundaryAxis:: boundary(i)%bTypes(s)%obj)
CASE DEFAULT
CALL criticalError('Boundary type ' // bType // ' undefined', 'readBoundary')
END SELECT
END DO

2
src/modules/moduleOutput.f90
View file

@ -3,7 +3,7 @@ MODULE moduleOutput
IMPLICIT NONE
!Output for each node
TYPE outputNode
REAL(8):: den, mom(1:3), tensorS(1:3,1:3)
REAL(8):: den = 0.D0, mom(1:3) = 0.D0, tensorS(1:3,1:3) = 0.D0
END TYPE

13
src/modules/moduleSolver.f90
View file

@ -311,10 +311,11 @@ MODULE moduleSolver
SUBROUTINE doReset()
USE moduleSpecies
USE moduleMesh
USE moduleList
IMPLICIT NONE
INTEGER:: nn, n
INTEGER:: nn, n, e
INTEGER, SAVE:: nPartNew
INTEGER, SAVE:: nInjIn, nOldIn, nWScheme
TYPE(particle), ALLOCATABLE, SAVE:: partTemp(:)
@ -347,6 +348,7 @@ MODULE moduleSolver
!$OMP SECTIONS
!$OMP SECTION
!Reset particles from injection
nn = 0
DO n = 1, nPartInj
IF (partInj(n)%n_in) THEN
@ -358,6 +360,7 @@ MODULE moduleSolver
END DO
!$OMP SECTION
!Reset particles from previous iteration
nn = nInjIn
DO n = 1, nPartOld
IF (partTemp(n)%n_in) THEN
@ -368,6 +371,7 @@ MODULE moduleSolver
END DO
!$OMP SECTION
!Reset particles from weighting scheme
nn = nInjIn + nOldIn
partCurr => partWScheme%head
DO n = 1, nWScheme
@ -381,6 +385,12 @@ MODULE moduleSolver
IF (ASSOCIATED(partWScheme%tail)) NULLIFY(partWScheme%tail)
partWScheme%amount = 0
!$OMP SECTION
!Reset output in nodes
DO e = 1, mesh%numNodes
CALL mesh%nodes(e)%obj%resetOutput()
END DO
!$OMP END SECTIONS
!$OMP SINGLE
@ -622,6 +632,5 @@ MODULE moduleSolver
END SUBROUTINE doOutput
END MODULE moduleSolver