diff --git a/src/modules/mesh/2DCyl/moduleMesh2DCyl.f90 b/src/modules/mesh/2DCyl/moduleMesh2DCyl.f90
index 2d300b9..ccf0c0f 100644
--- a/src/modules/mesh/2DCyl/moduleMesh2DCyl.f90
+++ b/src/modules/mesh/2DCyl/moduleMesh2DCyl.f90
@@ -247,7 +247,7 @@ MODULE moduleMesh2DCyl
     !QUAD FUNCTIONS
     !Init element
     SUBROUTINE initCellQuad2DCyl(self, n, p, nodes)
-      USE moduleRefParam
+      USE moduleRefParam, ONLY: L_ref
       IMPLICIT NONE
 
       CLASS(meshCell2DCylQuad), INTENT(out):: self
@@ -594,17 +594,17 @@ MODULE moduleMesh2DCyl
       r = DOT_PRODUCT(fPsi,self%r)
       !Computes total volume of the cell
       self%volume  = r*detJ*PI8 !2*pi * 4 (weight of 1 point 2D-Gaussian integral)
-      !Computes volume per node
-      Xi = (/-5.D-1, -5.D-1, 0.D0/)
+      !Computes volume per node. Change the radius point to calculate the area to improve accuracy near the axis.
+      Xi = (/-5.D-1, -0.25D0, 0.D0/)
       r = self%gatherF(Xi, 4, self%r)
       self%n1%v = self%n1%v + fPsi(1)*r*detJ*PI8
-      Xi = (/ 5.D-1, -5.D-1, 0.D0/)
+      Xi = (/ 5.D-1, -0.25D0, 0.D0/)
       r = self%gatherF(Xi, 4, self%r)
       self%n2%v = self%n2%v + fPsi(2)*r*detJ*PI8
-      Xi = (/ 5.D-1,  5.D-1, 0.D0/)
+      Xi = (/ 5.D-1,  0.75D0, 0.D0/)
       r = self%gatherF(Xi, 4, self%r)
       self%n3%v = self%n3%v + fPsi(3)*r*detJ*PI8
-      Xi = (/-5.D-1,  5.D-1, 0.D0/)
+      Xi = (/-5.D-1,  0.75D0, 0.D0/)
       r = self%gatherF(Xi, 4, self%r)
       self%n4%v = self%n4%v + fPsi(4)*r*detJ*PI8