Index: /trunk/tools/module_scientific.f90
===================================================================
--- /trunk/tools/module_scientific.f90	(revision 2662)
+++ /trunk/tools/module_scientific.f90	(revision 2663)
@@ -27,4 +27,5 @@
 ! fill3DR_2Dvec: Subroutine to fill a 3D float matrix from a series of indices from a 2D matrix
 ! gradient2D_1o: Subroutine to compute the first order gradient of a 2D field
+! gradient3D_1o: Subroutine to compute the first order gradient of a 3D field
 ! gradient2D_c1o: Subroutine to compute the first order centered gradient of a 2D field
 ! grid_within_polygon: Subroutine to determine which grid cells from a matrix lay inside a polygon
@@ -110,4 +111,11 @@
 ! write_overlap_polys_ascii: Subroutine to write to an ASCII file the associated polygons at a given time-step
 ! write_overlap_tracks_ascii: Subroutine to write to an ASCII the polygons associated to a trajectory at a given time step
+
+!! Vectorial calculus
+! deformation3D_1o: Subroutine to compute the first order deformation of a 3D field
+! gradient3D_1o: Subroutine to compute the first order gradient of a 3D field
+! matmodule3D: Subroutine to compute the module of a 3D matrix with 3 components
+! tilting3D_1o: Subroutine to compute the first order tilting of a 3D field
+! vecmodule3D: Function to compute the module of a 3D vector
 
 !!! *Functions/Subroutines to sort values adpated. The method used is usually referred to as "selection" method.
@@ -7641,4 +7649,39 @@
   END SUBROUTINE gradient2D_1o
 
+  SUBROUTINE gradient3D_1o(dx, dy, dz, var, dsx, dsy, dsz, grad)
+  ! Subroutine to compute the first order gradient of a 3D field
+  !   FROM: From: https://en.wikipedia.org/wiki/Finite_difference_coefficient
+
+    IMPLICIT NONE
+
+    INTEGER, INTENT(in)                                  :: dx, dy, dz
+    REAL(r_k), DIMENSION(dx,dy,dz), INTENT(in)           :: var, dsz
+    REAL(r_k), DIMENSION(dx,dy), INTENT(in)              :: dsx, dsy
+    REAL(r_k), DIMENSION(dx,dy,dz,3), INTENT(out)        :: grad
+
+! Local 
+    INTEGER                                              :: i, j, k
+
+!!!!!!! Variables
+! dx, dy, dz: shape of the 3D field
+! var: variable
+! dsx, dsy, dsz: matrices of distances betweeen grid points along each axis
+! grad: gradient
+
+    fname = 'gradient3D_1o'
+
+    grad = fillval64
+
+    DO i=1, dx-1
+      DO j=1, dy-1
+        DO k=1, dz-1
+          grad(i,j,k,:) = (/ (var(i+1,j,k)-var(i,j,k))/dsx(i,j), (var(i,j+1,k)-var(i,j,k))/dsy(i,j),  &
+             (var(i,j,k+1)-var(i,j,k))/dsz(i,j,k) /)
+        END DO
+      END DO
+    END DO
+
+  END SUBROUTINE gradient3D_1o
+
   SUBROUTINE gradient2D_c1o(dx, dy, var, dsx, dsy, grad)
   ! Subroutine to compute the first order centerd gradient of a 2D field
@@ -7833,5 +7876,4 @@
   END SUBROUTINE lap2D_1o
 
-
   SUBROUTINE lap2D_c1o(dx, dy, dsx, dsy, var, lap)
   ! Subroutine to compute the first order centered laplacian of a 2D vectorial field
@@ -7867,3 +7909,89 @@
   END SUBROUTINE lap2D_c1o
 
+  REAL(r_k) FUNCTION vecmodule3D(vec)
+  ! Function to compute the module of a 3D vector
+
+    IMPLICIT NONE
+
+    REAL(r_k), DIMENSION(3), INTENT(in)                  :: vec
+
+    fname = 'vecmodule3D'
+
+    vecmodule3D = SQRT(vec(1)*vec(1) + vec(2)*vec(2) + vec(3)*vec(3))
+
+  END FUNCTION vecmodule3D
+
+  SUBROUTINE FUNCTION matmodule3D(d1,d2,d3,mat)
+  ! Subroutine to compute the module of a 3D matrix with 3 components
+
+    IMPLICIT NONE
+
+    REAL(r_k), DIMENSION(d1,d2,d3,3), INTENT(in)         :: mat
+    REAL(r_k), DIMENSION(d1,d2,d3), INTENT(out)          :: matmodule
+
+    fname = 'matmodule3D'
+
+    matmodule3D = SQRT(vec(:,:,:,1)*vec(:,:,:,1) + vec(:,:,:,2)*vec(:,:,:,2) +                        &
+      vec(:,:,:,3)*vec(:,:,:,3))
+
+  END SUBROUTINE matmodule3D
+
+  SUBROUTINE deformation3D_1o(dx, dy, dz, vargrad, uagrad, vagrad, def)
+  ! Subroutine to compute the first order deformation of a 3D field
+  !   FROM: From: https://en.wikipedia.org/wiki/Finite_difference_coefficient
+
+    IMPLICIT NONE
+
+    INTEGER, INTENT(in)                                  :: dx, dy, dz
+    REAL(r_k), DIMENSION(dx,dy,dz,3), INTENT(in)         :: vargrad, uagrad, vagrad
+    REAL(r_k), DIMENSION(dx,dy,dz,3), INTENT(out)        :: def
+
+! Local 
+    INTEGER                                              :: i, j, k
+
+!!!!!!! Variables
+! dx, dy, dz: shape of the 3D field
+! vargrad: gradient of the variable
+! uagrad, vagrad: gradient of x/y components of the wind field
+! def: deformation
+
+    fname = 'deformation3D_1o'
+
+    def = fillval64
+
+    def(:,:,:,1) = - uagrad(:,:,:,1)*vargrad(:,:,:1) - vagrad(:,:,:,1)*vargrad(:,:,:2)
+    def(:,:,:,2) = - uagrad(:,:,:,2)*vargrad(:,:,:1) - vagrad(:,:,:,2)*vargrad(:,:,:2)
+    def(:,:,:,3) = - uagrad(:,:,:,3)*vargrad(:,:,:1) - vagrad(:,:,:,3)*vargrad(:,:,:2)
+
+  END SUBROUTINE deformation3D_1o
+
+  SUBROUTINE tilting3D_1o(dx, dy, dz, vargrad, wagrad, tilt)
+  ! Subroutine to compute the first order tilting of a 3D field
+  !   FROM: From: https://en.wikipedia.org/wiki/Finite_difference_coefficient
+
+    IMPLICIT NONE
+
+    INTEGER, INTENT(in)                                  :: dx, dy, dz
+    REAL(r_k), DIMENSION(dx,dy,dz,3), INTENT(in)         :: vargrad, wagrad
+    REAL(r_k), DIMENSION(dx,dy,dz,3), INTENT(out)        :: tilt
+
+! Local 
+    INTEGER                                              :: i, j, k
+
+!!!!!!! Variables
+! dx, dy, dz: shape of the 3D field
+! vargrad: gradient of the variable
+! wagrad: gradient of z component of the wind field
+! tilt: tilting
+
+    fname = 'tilting3D_1o'
+
+    tilt = fillval64
+
+    tilt(:,:,:,1) = - uagrad(:,:,:,1)*vargrad(:,:,:3)
+    tilt(:,:,:,2) = - vagrad(:,:,:,2)*vargrad(:,:,:3)
+    tilt(:,:,:,3) = - wagrad(:,:,:,3)*vargrad(:,:,:3)
+
+  END SUBROUTINE tilting3D_1o
+
 END MODULE module_scientific
