Rename of boundary module as I start with EM

src/modules/mesh/1DCart/moduleMesh1DCart.f90

@@ -4,11 +4,9 @@
 !              z == unused
 MODULE moduleMesh1DCart
   USE moduleMesh
+  use moduleMeshCommon
   IMPLICIT NONE
 
-  REAL(8), PARAMETER:: corSeg(1:3) = (/ -DSQRT(3.D0/5.D0), 0.D0,       DSQRT(3.D0/5.D0) /)
-  REAL(8), PARAMETER:: wSeg(1:3)   = (/        5.D0/9.D0 , 8.D0/9.D0,        5.D0/9.D0  /)
-
   TYPE, PUBLIC, EXTENDS(meshNode):: meshNode1DCart
     !Element coordinates
     REAL(8):: x = 0.D0
@@ -166,18 +164,6 @@ MODULE moduleMesh1DCart
 
     END FUNCTION randPosEdge
 
-    !Compute element functions at point Xi
-    pure function fPsiPoint(Xi, nNodes) RESULT(fPsi)
-      implicit none
-
-      real(8), intent(in):: Xi(1:3)
-      integer, intent(in):: nNodes
-      real(8):: fPsi(1:nNodes)
-
-      fPsi = (/ 1.D0 /)
-
-    end function fPsiPoint
-
     !CELL FUNCTIONS
     !SEGMENT FUNCTIONS
     !Init element
@@ -242,35 +228,6 @@ MODULE moduleMesh1DCart
 
     END FUNCTION randPos1DCartSegm
 
-    !Compute element functions at point Xi
-    PURE FUNCTION fPsiSegm(Xi, nNodes) RESULT(fPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: fPsi(1:nNodes)
-
-      fPsi = (/ 1.D0 - Xi(1), &
-                1.D0 + Xi(1) /)
-
-      fPsi    = fPsi * 0.50D0
-
-    END FUNCTION fPsiSegm
-
-    !Derivative element function at coordinates Xi
-    PURE FUNCTION dPsiSegm(Xi, nNodes) RESULT(dPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: dPsi(1:3,1:nNodes)
-
-      dPsi = 0.D0
-
-      dPsi(1, 1:2) = (/ -5.D-1, 5.D-1 /)
-
-    END FUNCTION dPsiSegm
-
     !Partial derivative in global coordinates
     PURE FUNCTION partialDerSegm(self, nNodes, dPsi) RESULT(pDer)
       IMPLICIT NONE

src/modules/mesh/1DRad/moduleMesh1DRad.f90

@@ -4,11 +4,9 @@
 !              z == unused
 MODULE moduleMesh1DRad
   USE moduleMesh
+  use moduleMeshCommon
   IMPLICIT NONE
 
-  REAL(8), PARAMETER:: corSeg(1:3) = (/ -DSQRT(3.D0/5.D0), 0.D0,       DSQRT(3.D0/5.D0) /)
-  REAL(8), PARAMETER:: wSeg(1:3)   = (/        5.D0/9.D0 , 8.D0/9.D0,        5.D0/9.D0  /)
-
   TYPE, PUBLIC, EXTENDS(meshNode):: meshNode1DRad
     !Element coordinates
     REAL(8):: r = 0.D0
@@ -166,18 +164,6 @@ MODULE moduleMesh1DRad
 
     END FUNCTION randPos1DRad
 
-    !Compute element functions at point Xi
-    pure function fPsiPoint(Xi, nNodes) RESULT(fPsi)
-      implicit none
-
-      real(8), intent(in):: Xi(1:3)
-      integer, intent(in):: nNodes
-      real(8):: fPsi(1:nNodes)
-
-      fPsi = (/ 1.D0 /)
-
-    end function fPsiPoint
-
     !CELL FUNCTIONS
     !SEGMENT FUNCTIONS
     !Init element
@@ -242,35 +228,6 @@ MODULE moduleMesh1DRad
 
     END FUNCTION randPos1DRadSegm
 
-    !Compute element functions at point Xi
-    PURE FUNCTION fPsiSegm(Xi, nNodes) RESULT(fPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: fPsi(1:nNodes)
-
-      fPsi = (/ 1.D0 - Xi(1), &
-                1.D0 + Xi(1) /)
-
-      fPsi    = fPsi * 0.50D0
-
-    END FUNCTION fPsiSegm
-
-    !Derivative element function at coordinates Xi
-    PURE FUNCTION dPsiSegm(Xi, nNodes) RESULT(dPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: dPsi(1:3,1:nNodes)
-
-      dPsi = 0.D0
-
-      dPsi(1, 1:2) = (/ -5.D-1, 5.D-1 /)
-
-    END FUNCTION dPsiSegm
-
     !Partial derivative in global coordinates
     PURE FUNCTION partialDerSegm(self, nNodes, dPsi) RESULT(pDer)
       IMPLICIT NONE

src/modules/mesh/2DCart/moduleMesh2DCart.f90

@@ -4,16 +4,9 @@
 !                  z == unused
 MODULE moduleMesh2DCart
   USE moduleMesh
+  use moduleMeshCommon
   IMPLICIT NONE
 
-  !Values for Gauss integral
-  REAL(8), PARAMETER:: corQuad(1:3) = (/ -DSQRT(3.D0/5.D0), 0.D0,      DSQRT(3.D0/5.D0) /)
-  REAL(8), PARAMETER:: wQuad(1:3)   = (/        5.D0/9.D0,  8.D0/9.D0,       5.D0/9.D0  /)
-
-  REAL(8), PARAMETER:: Xi1Tria(1:4) = (/   1.D0/3.D0,   1.D0/5.D0,   3.D0/5.D0,   1.D0/5.D0  /)
-  REAL(8), PARAMETER:: Xi2Tria(1:4) = (/   1.D0/3.D0,   1.D0/5.D0,   1.D0/5.D0,   3.D0/5.D0  /)
-  REAL(8), PARAMETER:: wTria(1:4)   = (/ -27.D0/96.D0, 25.D0/96.D0, 25.D0/96.D0, 25.D0/96.D0 /)
-
   TYPE, PUBLIC, EXTENDS(meshNode):: meshNode2DCart
     !Element coordinates
     REAL(8):: x = 0.D0, y = 0.D0
@@ -228,21 +221,6 @@ MODULE moduleMesh2DCart
 
     END FUNCTION randPosEdge
 
-    !Compute element functions at point Xi
-    pure function fPsiSegm(Xi, nNodes) RESULT(fPsi)
-      implicit none
-
-      real(8), intent(in):: Xi(1:3)
-      integer, intent(in):: nNodes
-      real(8):: fPsi(1:nNodes)
-
-      fPsi = (/ 1.D0 - Xi(1), &
-                1.D0 + Xi(1) /)
-
-      fPsi    = fPsi * 0.50D0
-
-    end function fPsiSegm
-
     !VOLUME FUNCTIONS
     !QUAD FUNCTIONS
     !Init element
@@ -319,49 +297,6 @@ MODULE moduleMesh2DCart
 
     END FUNCTION randPosCellQuad
 
-    !Compute element functions in point Xi
-    PURE FUNCTION fPsiQuad(Xi, nNodes) RESULT(fPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: fPsi(1:nNodes)
-
-      fPsi = 0.D0
-
-      fPsi = (/ (1.D0 - Xi(1)) * (1.D0 - Xi(2)), &
-                (1.D0 + Xi(1)) * (1.D0 - Xi(2)), &
-                (1.D0 + Xi(1)) * (1.D0 + Xi(2)), &
-                (1.D0 - Xi(1)) * (1.D0 + Xi(2)) /)
-
-      fPsi    = fPsi * 0.25D0
-
-    END FUNCTION fPsiQuad
-
-    !Derivative element function at coordinates Xi
-    PURE FUNCTION dPsiQuad(Xi, nNodes) RESULT(dPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: dPsi(1:3,1:nNodes)
-
-      dPsi = 0.D0
-
-      dPsi(1, 1:4) = (/ -(1.D0 - Xi(2)), &
-                         (1.D0 - Xi(2)), &
-                         (1.D0 + Xi(2)), &
-                        -(1.D0 + Xi(2)) /)
-
-      dPsi(2, 1:4) = (/ -(1.D0 - Xi(1)), &
-                        -(1.D0 + Xi(1)), &
-                         (1.D0 + Xi(1)), &
-                         (1.D0 - Xi(1)) /)
-
-      dPsi = dPsi * 0.25D0
-
-    END FUNCTION dPsiQuad
-
     !Partial derivative in global coordinates
     PURE FUNCTION partialDerQuad(self, nNodes, dPsi) RESULT(pDer)
       IMPLICIT NONE
@@ -665,35 +600,6 @@ MODULE moduleMesh2DCart
 
     END FUNCTION randPosCellTria
 
-    !Compute element functions in point Xi
-    PURE FUNCTION fPsiTria(Xi, nNodes) RESULT(fPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: fPsi(1:nNodes)
-
-      fPsi(1) = 1.D0 - Xi(1) - Xi(2)
-      fPsi(2) =        Xi(1)
-      fPsi(3) =                Xi(2)
-
-    END FUNCTION fPsiTria
-
-    !Compute element derivative functions in point Xi
-    PURE FUNCTION dPsiTria(Xi, nNodes) RESULT(dPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: dPsi(1:3,1:nNodes)
-
-      dPsi = 0.D0
-
-      dPsi(1,1:3) = (/ -1.D0, 1.D0, 0.D0 /)
-      dPsi(2,1:3) = (/ -1.D0, 0.D0, 1.D0 /)
-
-    END FUNCTION dPsiTria
-
     !Compute the derivatives in global coordinates
     PURE FUNCTION partialDerTria(self, nNodes, dPsi) RESULT(pDer)
       IMPLICIT NONE

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

@@ -4,16 +4,9 @@
 !               z == theta (unused)
 MODULE moduleMesh2DCyl
   USE moduleMesh
+  use moduleMeshCommon
   IMPLICIT NONE
 
-  !Values for Gauss integral
-  REAL(8), PARAMETER:: corQuad(1:3) = (/ -DSQRT(3.D0/5.D0), 0.D0,      DSQRT(3.D0/5.D0) /)
-  REAL(8), PARAMETER:: wQuad(1:3)   = (/        5.D0/9.D0,  8.D0/9.D0,       5.D0/9.D0  /)
-
-  REAL(8), PARAMETER:: Xi1Tria(1:4) = (/   1.D0/3.D0,   1.D0/5.D0,   3.D0/5.D0,   1.D0/5.D0  /)
-  REAL(8), PARAMETER:: Xi2Tria(1:4) = (/   1.D0/3.D0,   1.D0/5.D0,   1.D0/5.D0,   3.D0/5.D0  /)
-  REAL(8), PARAMETER:: wTria(1:4)   = (/ -27.D0/96.D0, 25.D0/96.D0, 25.D0/96.D0, 25.D0/96.D0 /)
-
   TYPE, PUBLIC, EXTENDS(meshNode):: meshNode2DCyl
     !Element coordinates
     REAL(8):: r = 0.D0, z = 0.D0
@@ -237,21 +230,6 @@ MODULE moduleMesh2DCyl
 
     END FUNCTION randPosEdge
 
-    !Compute element functions at point Xi
-    pure function fPsiSegm(Xi, nNodes) RESULT(fPsi)
-      implicit none
-
-      real(8), intent(in):: Xi(1:3)
-      integer, intent(in):: nNodes
-      real(8):: fPsi(1:nNodes)
-
-      fPsi = (/ 1.D0 - Xi(1), &
-                1.D0 + Xi(1) /)
-
-      fPsi    = fPsi * 0.50D0
-
-    end function fPsiSegm
-
     !CELL FUNCTIONS
     !QUAD FUNCTIONS
     !Init element
@@ -327,49 +305,6 @@ MODULE moduleMesh2DCyl
 
     END FUNCTION randPosCellQuad
 
-    !Computes element functions in point Xi
-    PURE FUNCTION fPsiQuad(Xi, nNodes) RESULT(fPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: fPsi(1:nNodes)
-
-      fPsi = 0.D0
-
-      fPsi = (/ (1.D0 - Xi(1)) * (1.D0 - Xi(2)), &
-                (1.D0 + Xi(1)) * (1.D0 - Xi(2)), &
-                (1.D0 + Xi(1)) * (1.D0 + Xi(2)), &
-                (1.D0 - Xi(1)) * (1.D0 + Xi(2)) /)
-
-      fPsi    = fPsi * 0.25D0
-
-    END FUNCTION fPsiQuad
-
-    !Derivative element function at coordinates Xi
-    PURE FUNCTION dPsiQuad(Xi, nNodes) RESULT(dPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: dPsi(1:3,1:nNodes)
-
-      dPsi = 0.D0
-
-      dPsi(1, 1:4) = (/ -(1.D0 - Xi(2)), &
-                         (1.D0 - Xi(2)), &
-                         (1.D0 + Xi(2)), &
-                        -(1.D0 + Xi(2)) /)
-
-      dPsi(2, 1:4) = (/ -(1.D0 - Xi(1)), &
-                        -(1.D0 + Xi(1)), &
-                         (1.D0 + Xi(1)), &
-                         (1.D0 - Xi(1)) /)
-
-      dPsi = dPsi * 0.25D0
-
-    END FUNCTION dPsiQuad
-
     !Partial derivative in global coordinates
     PURE FUNCTION partialDerQuad(self, nNodes, dPsi) RESULT(pDer)
       IMPLICIT NONE
@@ -692,35 +627,6 @@ MODULE moduleMesh2DCyl
 
     END FUNCTION randPosCellTria
 
-    !Compute element functions in point Xi
-    PURE FUNCTION fPsiTria(Xi, nNodes) RESULT(fPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: fPsi(1:nNodes)
-
-      fPsi(1) = 1.D0 - Xi(1) - Xi(2)
-      fPsi(2) =        Xi(1)
-      fPsi(3) =                Xi(2)
-
-    END FUNCTION fPsiTria
-
-    !Compute element derivative functions in point Xi
-    PURE FUNCTION dPsiTria(Xi, nNodes) RESULT(dPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: dPsi(1:3,1:nNodes)
-
-      dPsi = 0.D0
-
-      dPsi(1,1:3) = (/ -1.D0, 1.D0, 0.D0 /)
-      dPsi(2,1:3) = (/ -1.D0, 0.D0, 1.D0 /)
-
-    END FUNCTION dPsiTria
-
     !Compute the derivatives in global coordinates
     PURE FUNCTION partialDerTria(self, nNodes, dPsi) RESULT(pDer)
       IMPLICIT NONE

src/modules/mesh/3DCart/moduleMesh3DCart.f90

@@ -4,6 +4,7 @@
 !                  z == z
 MODULE moduleMesh3DCart
   USE moduleMesh
+  use moduleMeshCommon
   IMPLICIT NONE
 
   TYPE, PUBLIC, EXTENDS(meshNode):: meshNode3DCart
@@ -204,20 +205,6 @@ MODULE moduleMesh3DCart
 
     END FUNCTION randPosEdgeTria
 
-    !Shape functions for triangular surface
-    PURE FUNCTION fPsiTria(Xi, nNodes) RESULT(fPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: fPsi(1:nNodes)
-
-      fPsi(1) = 1.D0 - Xi(1) - Xi(2)
-      fPsi(2) =        Xi(1)
-      fPsi(3) =                Xi(2)
-
-    END FUNCTION fPsiTria
-
     !VOLUME FUNCTIONS
     !TETRA FUNCTIONS
     !Init element
@@ -295,37 +282,6 @@ MODULE moduleMesh3DCart
 
     END FUNCTION randPosCellTetra
 
-    !Compute element functions in point Xi
-    PURE FUNCTION fPsiTetra(Xi, nNodes) RESULT(fPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: fPsi(1:nNodes)
-
-      fPsi(1) = 1.D0 - Xi(1) - Xi(2) - Xi(3)
-      fPsi(2) =        Xi(1)
-      fPsi(3) =                Xi(2)
-      fPsi(4) =                        Xi(3)
-
-    END FUNCTION fPsiTetra
-
-    !Compute element derivative functions in point Xi
-    PURE FUNCTION dPsiTetra(Xi, nNodes) RESULT(dPsi)
-      IMPLICIT NONE
-
-      REAL(8), INTENT(in):: Xi(1:3)
-      INTEGER, INTENT(in):: nNodes
-      REAL(8):: dPsi(1:3, 1:nNodes)
-
-      dPsi = 0.D0
-
-      dPsi(1,1:4) = (/ -1.D0, 1.D0, 0.D0, 0.D0 /)
-      dPsi(2,1:4) = (/ -1.D0, 0.D0, 1.D0, 0.D0 /)
-      dPsi(3,1:4) = (/ -1.D0, 0.D0, 0.D0, 1.D0 /)
-
-    END FUNCTION dPsiTetra
-
     !Compute the derivatives in global coordinates
     PURE FUNCTION partialDerTetra(self, nNodes, dPsi) RESULT(pDer)
       IMPLICIT NONE

src/modules/mesh/makefile

@@ -18,7 +18,7 @@ all: moduleMesh.o inout.o 3DCart.o 2DCyl.o 2DCart.o 1DRad.o 1DCart.o 0D.o
 0D.o: moduleMesh.o
 	$(MAKE) -C 0D all
 
-moduleMesh.o: moduleMesh.f90
+moduleMesh.o: moduleMeshCommon.o moduleMesh.f90
 	$(FC) $(FCFLAGS) -c $(subst .o,.f90,$@) -o $(OBJDIR)/$@
 	$(FC) $(FCFLAGS) -c moduleMesh@elements.f90 -o $(OBJDIR)/$@
 	$(FC) $(FCFLAGS) -c moduleMesh@boundary.f90 -o $(OBJDIR)/$@

src/modules/mesh/moduleMesh@boundaryParticle.f90

@@ -1,5 +1,5 @@
 !moduleMeshBoundary: Boundary functions for the mesh edges
-submodule(moduleMesh) boundary
+submodule(moduleMesh) boundaryParticle
   contains
     module function boundaryParticleName_to_Index(boundaryName) result(bp)
       use moduleErrors
@@ -502,4 +502,4 @@ submodule(moduleMesh) boundary
 
     end function physicalSurface_to_id
 
-end submodule boundary
+end submodule boundaryParticle

src/modules/mesh/moduleMesh@common.f90

@@ -1,3 +0,0 @@
-module moduleMeshCommon
-
-end module moduleMeshCommon

src/modules/mesh/moduleMeshCommon.f90

@@ -0,0 +1,164 @@
+! Common subroutines and parameters for mesh elements
+module moduleMeshCommon
+  ! values for integration gauss integral
+  ! Segment
+  real(8), parameter:: corSeg(1:3)  = (/ -dsqrt(3.D0/5.D0),  0.D0,      dsqrt(3.D0/5.D0) /)
+  real(8), parameter:: wSeg(1:3)    = (/        5.D0/9.D0,   8.D0/9.D0,       5.D0/9.D0  /)
+
+  ! Quad
+  real(8), parameter:: corQuad(1:3) = (/ -dsqrt(3.D0/5.D0),  0.D0,      dsqrt(3.D0/5.D0) /)
+  real(8), parameter:: wQuad(1:3)   = (/        5.D0/9.D0,   8.D0/9.D0,       5.D0/9.D0  /)
+
+  ! tria
+  real(8), parameter:: Xi1Tria(1:4) = (/        1.D0/3.D0,   1.D0/5.D0,       3.D0/5.D0,   1.D0/5.D0  /)
+  real(8), parameter:: Xi2Tria(1:4) = (/        1.D0/3.D0,   1.D0/5.D0,       1.D0/5.D0,   3.D0/5.D0  /)
+  real(8), parameter:: wTria(1:4)   = (/      -27.D0/96.D0, 25.D0/96.D0,     25.D0/96.D0, 25.D0/96.D0 /)
+
+  contains
+    ! ELEMENT FUNCTIONS
+    ! Point
+    pure function fPsiPoint(Xi, nNodes) RESULT(fPsi)
+      implicit none
+
+      real(8), intent(in):: Xi(1:3)
+      integer, intent(in):: nNodes
+      real(8):: fPsi(1:nNodes)
+
+      fPsi = (/ 1.D0 /)
+
+    end function fPsiPoint
+
+    ! Segment
+    pure function fPsiSegm(Xi, nNodes) result(fPsi)
+      implicit none
+
+      real(8), intent(in):: Xi(1:3)
+      integer, intent(in):: nNodes
+      real(8):: fPsi(1:nNodes)
+
+      fPsi = (/ 1.D0 - Xi(1), &
+                1.D0 + Xi(1) /)
+
+      fPsi    = fPsi * 0.50D0
+
+    end function fPsiSegm
+
+    ! Quad
+    pure function fPsiQuad(Xi, nNodes) result(fPsi)
+      implicit none
+
+      real(8), intent(in):: Xi(1:3)
+      integer, intent(in):: nNodes
+      real(8):: fPsi(1:nNodes)
+
+      fPsi = 0.D0
+
+      fPsi = (/ (1.D0 - Xi(1)) * (1.D0 - Xi(2)), &
+                (1.D0 + Xi(1)) * (1.D0 - Xi(2)), &
+                (1.D0 + Xi(1)) * (1.D0 + Xi(2)), &
+                (1.D0 - Xi(1)) * (1.D0 + Xi(2)) /)
+
+      fPsi    = fPsi * 0.25D0
+
+    end function fPsiQuad
+
+    ! Tria
+    pure function fPsiTria(Xi, nNodes) result(fPsi)
+      IMPLICIT NONE
+
+      real(8), intent(in):: Xi(1:3)
+      integer, intent(in):: nNodes
+      real(8):: fPsi(1:nNodes)
+
+      fPsi(1) = 1.D0 - Xi(1) - Xi(2)
+      fPsi(2) =        Xi(1)
+      fPsi(3) =                Xi(2)
+
+    end function fPsiTria
+
+    ! Tetra
+    pure function fPsiTetra(Xi, nNodes) result(fPsi)
+      implicit none
+
+      real(8), intent(in):: Xi(1:3)
+      integer, intent(in):: nNodes
+      real(8):: fPsi(1:nNodes)
+
+      fPsi(1) = 1.D0 - Xi(1) - Xi(2) - Xi(3)
+      fPsi(2) =        Xi(1)
+      fPsi(3) =                Xi(2)
+      fPsi(4) =                        Xi(3)
+
+    end function fPsiTetra
+
+    ! DERIVATIVE FOR ELEMENT FUNCTIONS
+    ! Segment
+    pure function dPsiSegm(Xi, nNodes) result(dPsi)
+      implicit none
+
+      real(8), intent(in):: Xi(1:3)
+      integer, intent(in):: nNodes
+      real(8):: dPsi(1:3,1:nNodes)
+
+      dPsi = 0.D0
+
+      dPsi(1, 1:2) = (/ -5.D-1, 5.D-1 /)
+
+    end function dPsiSegm
+
+    ! Quad
+    pure function dPsiQuad(Xi, nNodes) result(dPsi)
+      implicit none
+
+      real(8), intent(in):: Xi(1:3)
+      integer, intent(in):: nNodes
+      real(8):: dPsi(1:3,1:nNodes)
+
+      dPsi = 0.D0
+
+      dPsi(1, 1:4) = (/ -(1.D0 - Xi(2)), &
+                         (1.D0 - Xi(2)), &
+                         (1.D0 + Xi(2)), &
+                        -(1.D0 + Xi(2)) /)
+
+      dPsi(2, 1:4) = (/ -(1.D0 - Xi(1)), &
+                        -(1.D0 + Xi(1)), &
+                         (1.D0 + Xi(1)), &
+                         (1.D0 - Xi(1)) /)
+
+      dPsi = dPsi * 0.25D0
+
+    end function dPsiQuad
+
+    ! Tria
+    pure function dPsiTria(Xi, nNodes) result(dPsi)
+      implicit none
+
+      real(8), intent(in):: Xi(1:3)
+      integer, intent(in):: nNodes
+      real(8):: dPsi(1:3,1:nNodes)
+
+      dPsi = 0.D0
+
+      dPsi(1,1:3) = (/ -1.D0, 1.D0, 0.D0 /)
+      dPsi(2,1:3) = (/ -1.D0, 0.D0, 1.D0 /)
+
+    end function dPsiTria
+
+    !Compute element derivative functions in point Xi
+    pure function dPsiTetra(Xi, nNodes) result(dPsi)
+      implicit none
+
+      real(8), intent(in):: Xi(1:3)
+      integer, intent(in):: nNodes
+      real(8):: dPsi(1:3, 1:nNodes)
+
+      dPsi = 0.D0
+
+      dPsi(1,1:4) = (/ -1.D0, 1.D0, 0.D0, 0.D0 /)
+      dPsi(2,1:4) = (/ -1.D0, 0.D0, 1.D0, 0.D0 /)
+      dPsi(3,1:4) = (/ -1.D0, 0.D0, 0.D0, 1.D0 /)
+
+    end function dPsiTetra
+
+end module moduleMeshCommon