Merge branch 'issue/volumeNode' into 'feature/average'

Fix for node volume

See merge request JorgeGonz/fpakc!28

src/modules/mesh/2DCyl/moduleMesh2DCyl.f90

@@ -249,7 +249,7 @@ MODULE moduleMesh2DCyl
       dr  = self%r(2) - self%r(1)
       dz  = self%z(2) - self%z(1)
       IF (dr /= 0.D0) THEN
-        r(2) = dr*DSQRT(rnd) + self%r(1)
+        r(2) = dr * DSQRT(rnd) + self%r(1)
         r(1) = dz * (r(2) - self%r(1))/dr + self%z(1)
 
       ELSE
@@ -320,7 +320,18 @@ MODULE moduleMesh2DCyl
       fPsi = self%fPsi(xi)
       r = DOT_PRODUCT(fPsi,self%r)
       self%volume  =      r*detJ
-      self%arNodes = fPsi*r*detJ
+      xi = (/-5.D-1, -5.D-1, 0.D0/)
+      r = DOT_PRODUCT(self%fPsi(xi),self%r)
+      self%arNodes(1) = fPsi(1)*r*detJ
+      xi = (/ 5.D-1, -5.D-1, 0.D0/)
+      r = DOT_PRODUCT(self%fPsi(xi),self%r)
+      self%arNodes(2) = fPsi(2)*r*detJ
+      xi = (/ 5.D-1,  5.D-1, 0.D0/)
+      r = DOT_PRODUCT(self%fPsi(xi),self%r)
+      self%arNodes(3) = fPsi(3)*r*detJ
+      xi = (/-5.D-1,  5.D-1, 0.D0/)
+      r = DOT_PRODUCT(self%fPsi(xi),self%r)
+      self%arNodes(4) = fPsi(4)*r*detJ
 
     END SUBROUTINE areaQuad
 

src/modules/mesh/inout/gmsh2/moduleMeshOutputGmsh2.f90

@@ -151,8 +151,9 @@ MODULE moduleMeshOutputGmsh2
       INTEGER:: n
       REAL(8):: time
       CHARACTER(:), ALLOCATABLE:: fileName
-      CHARACTER (LEN=iterationDigits):: tstring
+      CHARACTER (LEN=iterationDigits):: tString
       CHARACTER(:), ALLOCATABLE:: title
+      CHARACTER (LEN=2):: cString
 
       SELECT TYPE(self)
       TYPE IS(meshParticles)
@@ -168,9 +169,9 @@ MODULE moduleMeshOutputGmsh2
 
       IF (collOutput) THEN
         time = DBLE(t)*tauMin*ti_ref
-        WRITE(tstring, iterationFormat) t
+        WRITE(tString, iterationFormat) t
 
-        fileName='OUTPUT_' // tstring// '_Collisions.msh'
+        fileName='OUTPUT_' // tString// '_Collisions.msh'
         WRITE(*, "(6X,A15,A)") "Creating file: ", fileName
         OPEN (60, file = path // folder // '/' //  fileName)
 
@@ -178,7 +179,8 @@ MODULE moduleMeshOutputGmsh2
 
         DO k = 1, nCollPairs
           DO c = 1, interactionMatrix(k)%amount
-            WRITE(title, "(5A,I2)") '"Pair ', interactionMatrix(k)%sp_i%name, '-', interactionMatrix(k)%sp_j%name, ' collision ', c
+            WRITE(cString, "(I2)") c
+            title = '"Pair ' // interactionMatrix(k)%sp_i%name // '-' // interactionMatrix(k)%sp_j%name // ' collision ' // cString
             CALL writeGmsh2HeaderElementData(60, title, t, time, 1, self%numVols)
             DO n=1, self%numVols
               WRITE(60, "(I6,I10)") n + numEdges, self%vols(n)%obj%tallyColl(k)%tally(c)

src/modules/moduleSolver.f90

@@ -660,24 +660,18 @@ MODULE moduleSolver
       TYPE(particle), INTENT(inout):: part
       CLASS(meshVol), POINTER, INTENT(in):: volOld
       CLASS(meshVol), POINTER, INTENT(inout):: volNew
-      REAL(8):: fractionVolume, fractionWeight
-      INTEGER:: nSplit
+      REAL(8):: fractionVolume, pSplit
 
-      !If particle has change cell, call Weighting scheme
-      IF (volOld%n /= volNew%n) THEN
+      !If particle changes volume to smaller cell
+      IF (volOld%volume > volNew%volume) THEN
         fractionVolume = volOld%volume/volNew%volume
 
-        part%weight = part%weight * fractionVolume
+        !Calculate probability of splitting particle
+        pSplit = 1.D0 - DEXP(-fractionVolume)
 
-        fractionWeight = part%weight / part%species%weight
-
-        IF (fractionWeight >= 2.D0) THEN
-          nSplit = FLOOR(fractionWeight)
-          CALL splitParticle(part, nSplit, volNew)
-
-        ELSEIF (part%weight < 1.D0) THEN
-          !Particle has lost statistical meaning and will be terminated
-          part%n_in = .FALSE.
+        IF (random() < pSplit THEN
+          !Split particle in two
+          CALL splitParticle(part, 2, volNew)
 
         END IF