Some progress

Fixed an issue with random integer numbers.
Cylindrical coordinates are not perfect yet:
  - Box (cylinder) with initial constant density loses particles at r =
    0
  - Injection density still low in r = 0

src/modules/common/moduleRandom.f90

@@ -43,7 +43,7 @@ MODULE moduleRandom
       rnd = 0
       CALL RANDOM_NUMBER(rnd01)
 
-      rnd = NINT(REAL(b - a) * rnd01) + a
+      rnd = a + FLOOR((b+1-a)*rnd01)
 
     END FUNCTION randomIntAB
 

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

@@ -226,11 +226,6 @@ MODULE moduleMesh2DCyl
       REAL(8):: p1(1:2), p2(1:2)
 
       rnd = random()
-      ! IF (self%r(1) == 0.D0  .OR. &
-      !     self%r(2) == 0.D0) THEN
-      !   rnd = rnd**(1.D0/3.D0)
-      !
-      ! END IF
 
       p1 = (/self%z(1), self%r(1) /)
       p2 = (/self%z(2), self%r(2) /)

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

@@ -296,6 +296,20 @@ MODULE moduleMeshInputGmsh2
 
       CLOSE(10)
 
+      ! Adjust node volume at axis
+      SELECT CASE(self%geometry)
+      CASE("Cyl")
+        DO n  = 1, numNodes
+          r = self%nodes(n)%obj%getCoordinates()
+          IF (r(2) == 0.D0) THEN
+            self%nodes(n)%obj%v = self%nodes(n)%obj%v * 2.0D0!2.0D0/3.0D0
+
+          END IF
+
+        END DO
+
+      END SELECT
+
       !Call mesh connectivity
       CALL self%connectMesh
 

src/modules/mesh/inout/vtu/moduleMeshInputVTU.f90

@@ -496,6 +496,20 @@ MODULE moduleMeshInputVTU
 
       END DO
 
+      ! Adjust node volume at axis
+      SELECT CASE(self%geometry)
+      CASE("Cyl")
+        DO n  = 1, numNodes
+          r = self%nodes(n)%obj%getCoordinates()
+          IF (r(2) == 0.D0) THEN
+            self%nodes(n)%obj%v = self%nodes(n)%obj%v * 2.0D0!2.0D0/3.0D0
+
+          END IF
+
+        END DO
+
+      END SELECT
+        
       !Call mesh connectivity
       CALL self%connectMesh
 

src/modules/moduleInject.f90

@@ -313,12 +313,12 @@ MODULE moduleInject
 
       !$OMP DO
       DO n = nMin, nMax
-        randomX    =  random(1, self%nEdges)
+        randomX    =  randomWeighted(self%cumWeight, self%sumWeight)
         randomEdge => mesh%edges(self%edges(randomX))%obj
         !Random position in edge
         partInj(n)%r = randomEdge%randPos()
         !Assign weight to particle.
-        partInj(n)%weight = self%species%weight * sqrt(partInj(n)%r(2) / (0.1D0/L_ref))
+        partInj(n)%weight = self%species%weight
         !Volume associated to the edge:
         IF (ASSOCIATED(randomEdge%e1)) THEN
           partInj(n)%cell = randomEdge%e1%n