diff --git a/fvn_sparse/fvn_sparse.f90 b/fvn_sparse/fvn_sparse.f90
index 7fdb5ef..1fc1ad3 100644
--- a/fvn_sparse/fvn_sparse.f90
+++ b/fvn_sparse/fvn_sparse.f90
@@ -7,6 +7,9 @@ interface fvn_sparse_solve
     module procedure fvn_zl_sparse_solve,fvn_zi_sparse_solve,fvn_dl_sparse_solve,fvn_di_sparse_solve
 end interface fvn_sparse_solve
 
+interface fvn_sparse_det
+    module procedure fvn_zl_sparse_det,fvn_zi_sparse_det,fvn_dl_sparse_det,fvn_di_sparse_det
+end interface fvn_sparse_det
 contains
 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 !
@@ -30,7 +33,7 @@ contains
 ! * = zl : double complex + integer(kind=dp_kind)
 ! * = zi : double complex + integer(kind=sp_kind)
 !
-subroutine fvn_zl_sparse_solve(n,nz,T,Ti,Tj,B,x,status)
+subroutine fvn_zl_sparse_solve(n,nz,T,Ti,Tj,B,x,status,det)
 implicit none
 integer(kind=dp_kind), intent(in) :: n,nz
 complex(kind=dp_kind),dimension(nz),intent(in) :: T
@@ -38,6 +41,7 @@ integer(kind=dp_kind),dimension(nz),intent(in)  :: Ti,Tj
 complex(kind=dp_kind),dimension(n),intent(in) :: B
 complex(kind=dp_kind),dimension(n),intent(out) :: x
 integer(kind=dp_kind), intent(out) :: status
+complex(kind=dp_kind), optional, intent(out) :: det
 
 integer(kind=dp_kind),dimension(:),allocatable :: wTi,wTj
 real(kind=dp_kind),dimension(:),allocatable :: Tx,Tz
@@ -48,6 +52,7 @@ real(kind=dp_kind),dimension(:),allocatable :: xx,xz,bx,bz
 real(kind=dp_kind),dimension(90) :: info
 real(kind=dp_kind),dimension(20) :: control
 integer(kind=dp_kind) :: sys
+real(kind=dp_kind) :: Mx,Mz,Ex
 
 
 status=0
@@ -92,6 +97,19 @@ endif
 ! free the C symbolic pointer
 call umfpack_zl_free_symbolic (symbolic)
 
+! if parameter det is present, the determinant of the matrix is calculated
+if (present(det) ) then
+  call umfpack_zl_get_determinant(Mx,Mz,Ex,numeric,info,status)
+  det=cmplx(Mx,Mz,kind=dp_kind)*10**Ex
+  ! info(1) should be zero
+  if (info(1) /= 0) then
+      write(*,*) "Problem during sparse determinant, returned code : ",info(1)
+      status=info(1)
+  endif
+endif
+
+
+
 allocate(bx(n),bz(n),xx(n),xz(n))
 bx=dble(B)
 bz=aimag(B)
@@ -122,7 +140,7 @@ end subroutine
 
 
 
-subroutine fvn_zi_sparse_solve(n,nz,T,Ti,Tj,B,x,status)
+subroutine fvn_zi_sparse_solve(n,nz,T,Ti,Tj,B,x,status,det)
 implicit none
 integer(kind=sp_kind), intent(in) :: n,nz
 complex(kind=dp_kind),dimension(nz),intent(in) :: T
@@ -130,6 +148,7 @@ integer(kind=sp_kind),dimension(nz),intent(in)  :: Ti,Tj
 complex(kind=dp_kind),dimension(n),intent(in) :: B
 complex(kind=dp_kind),dimension(n),intent(out) :: x
 integer(kind=sp_kind), intent(out) :: status
+complex(kind=dp_kind), optional, intent(out) :: det
 
 integer(kind=sp_kind),dimension(:),allocatable :: wTi,wTj
 real(kind=dp_kind),dimension(:),allocatable :: Tx,Tz
@@ -144,6 +163,7 @@ real(kind=dp_kind),dimension(:),allocatable :: xx,xz,bx,bz
 real(kind=dp_kind),dimension(90) :: info
 real(kind=dp_kind),dimension(20) :: control
 integer(kind=sp_kind) :: sys
+real(kind=dp_kind) :: Mx,Mz,Ex
 
 status=0
 ! we use a working copy of Ti and Tj to perform 1-based to 0-based translation
@@ -186,6 +206,20 @@ endif
 ! free the C symbolic pointer
 call umfpack_zi_free_symbolic (symbolic)
 
+! if parameter det is present, the determinant of the matrix is calculated
+if (present(det) ) then
+  call umfpack_zi_get_determinant(Mx,Mz,Ex,numeric,info,status)
+  det=cmplx(Mx,Mz,kind=dp_kind)*10**Ex
+  ! info(1) should be zero
+  if (info(1) /= 0) then
+      write(*,*) "Problem during sparse determinant, returned code : ",info(1)
+      status=info(1)
+  endif
+endif
+
+
+
+
 allocate(bx(n),bz(n),xx(n),xz(n))
 bx=dble(B)
 bz=aimag(B)
@@ -216,7 +250,7 @@ end subroutine
 
 
 
-subroutine fvn_dl_sparse_solve(n,nz,T,Ti,Tj,B,x,status)
+subroutine fvn_dl_sparse_solve(n,nz,T,Ti,Tj,B,x,status,det)
 implicit none
 integer(kind=dp_kind), intent(in) :: n,nz
 real(kind=dp_kind),dimension(nz),intent(in) :: T
@@ -224,6 +258,7 @@ integer(kind=dp_kind),dimension(nz),intent(in)  :: Ti,Tj
 real(kind=dp_kind),dimension(n),intent(in) :: B
 real(kind=dp_kind),dimension(n),intent(out) :: x
 integer(kind=dp_kind), intent(out) :: status
+real(kind=dp_kind), optional, intent(out) :: det
 
 integer(kind=dp_kind),dimension(:),allocatable :: wTi,wTj
 real(kind=dp_kind),dimension(:),allocatable :: A
@@ -233,6 +268,7 @@ integer(kind=dp_kind) :: symbolic,numeric
 real(kind=dp_kind),dimension(90) :: info
 real(kind=dp_kind),dimension(20) :: control
 integer(kind=dp_kind) :: sys
+real(kind=dp_kind) :: Mx,Ex
 
 status=0
 ! we use a working copy of Ti and Tj to perform 1-based to 0-based translation
@@ -271,6 +307,18 @@ endif
 ! free the C symbolic pointer
 call umfpack_dl_free_symbolic (symbolic)
 
+! if parameter det is present, the determinant of the matrix is calculated
+if (present(det) ) then
+  call umfpack_dl_get_determinant(Mx,Ex,numeric,info,status)
+  det=Mx*10**Ex
+  ! info(1) should be zero
+  if (info(1) /= 0) then
+      write(*,*) "Problem during sparse determinant, returned code : ",info(1)
+      status=info(1)
+  endif
+endif
+
+
 sys=0
 ! sys may be used to define type of solving -> see umfpack.h
 
@@ -294,7 +342,7 @@ end subroutine
 
 
 
-subroutine fvn_di_sparse_solve(n,nz,T,Ti,Tj,B,x,status)
+subroutine fvn_di_sparse_solve(n,nz,T,Ti,Tj,B,x,status,det)
 implicit none
 integer(kind=sp_kind), intent(in) :: n,nz
 real(kind=dp_kind),dimension(nz),intent(in) :: T
@@ -302,6 +350,7 @@ integer(kind=sp_kind),dimension(nz),intent(in)  :: Ti,Tj
 real(kind=dp_kind),dimension(n),intent(in) :: B
 real(kind=dp_kind),dimension(n),intent(out) :: x
 integer(kind=sp_kind), intent(out) :: status
+real(kind=dp_kind), optional, intent(out) :: det
 
 integer(kind=sp_kind),dimension(:),allocatable :: wTi,wTj
 real(kind=dp_kind),dimension(:),allocatable :: A
@@ -314,6 +363,7 @@ integer(kind=sp_kind),dimension(2) :: symbolic,numeric
 real(kind=dp_kind),dimension(90) :: info
 real(kind=dp_kind),dimension(20) :: control
 integer(kind=sp_kind) :: sys
+real(kind=dp_kind) :: Mx,Ex
 
 status=0
 ! we use a working copy of Ti and Tj to perform 1-based to 0-based translation
@@ -352,9 +402,19 @@ endif
 ! free the C symbolic pointer
 call umfpack_di_free_symbolic (symbolic)
 
+! if parameter det is present, the determinant of the matrix is calculated
+if (present(det) ) then
+  call umfpack_di_get_determinant(Mx,Ex,numeric,info,status)
+  det=Mx*10**Ex
+  ! info(1) should be zero
+  if (info(1) /= 0) then
+      write(*,*) "Problem during sparse determinant, returned code : ",info(1)
+      status=info(1)
+  endif
+endif
+
 sys=0
 ! sys may be used to define type of solving -> see umfpack.h
-
 ! Solving
 call umfpack_di_solve (sys, Ap, Ai, A, x, B, numeric, control, info)
 ! info(1) should be zero
@@ -371,4 +431,313 @@ deallocate(wTi,wTj)
 end subroutine
 
 
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+!
+! SPARSE DETERMINANT
+!
+!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+subroutine fvn_zl_sparse_det(n,nz,T,Ti,Tj,det,status)
+implicit none
+integer(kind=dp_kind), intent(in) :: n,nz
+complex(kind=dp_kind),dimension(nz),intent(in) :: T
+integer(kind=dp_kind),dimension(nz),intent(in)  :: Ti,Tj
+complex(kind=dp_kind),intent(out) :: det
+integer(kind=dp_kind), intent(out) :: status
+
+integer(kind=dp_kind),dimension(:),allocatable :: wTi,wTj
+real(kind=dp_kind),dimension(:),allocatable :: Tx,Tz
+real(kind=dp_kind),dimension(:),allocatable :: Ax,Az
+integer(kind=dp_kind),dimension(:),allocatable :: Ap,Ai
+integer(kind=dp_kind) :: symbolic,numeric
+real(kind=dp_kind),dimension(90) :: info
+real(kind=dp_kind),dimension(20) :: control
+real(kind=dp_kind) :: Mx,Mz,Ex
+
+
+status=0
+
+! we use a working copy of Ti and Tj to perform 1-based to 0-based translation
+! Tx and Tz are the real and imaginary parts of T
+allocate(wTi(nz),wTj(nz))
+allocate(Tx(nz),Tz(nz))
+Tx=dble(T)
+Tz=aimag(T)
+wTi=Ti-1
+wTj=Tj-1
+allocate(Ax(nz),Az(nz))
+allocate(Ap(n+1),Ai(nz))
+
+! perform the triplet to compressed column form -> Ap,Ai,Ax,Az
+call umfpack_zl_triplet_to_col(n,n,nz,wTi,wTj,Tx,Tz,Ap,Ai,Ax,Az,status)
+! if status is not zero a problem has occured
+if (status /= 0) then
+    write(*,*) "Problem during umfpack_zl_triplet_to_col"
+endif
+
+! Define defaults control values
+call umfpack_zl_defaults(control)
+
+! Symbolic analysis
+call umfpack_zl_symbolic(n,n,Ap,Ai,Ax,Az,symbolic, control, info)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during symbolic analysis"
+    status=info(1)
+endif
+
+! Numerical factorization
+call umfpack_zl_numeric (Ap, Ai, Ax, Az, symbolic, numeric, control, info)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during numerical factorization"
+    status=info(1)
+endif
+
+! free the C symbolic pointer
+call umfpack_zl_free_symbolic (symbolic)
+
+call umfpack_zl_get_determinant(Mx,Mz,Ex,numeric,info,status)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during sparse determinant, returned code : ",info(1)
+    status=info(1)
+endif
+det=cmplx(Mx,Mz,kind=dp_kind)*10**Ex
+
+! free the C numeric pointer
+call umfpack_zl_free_numeric (numeric)
+
+deallocate(Ax,Az)
+deallocate(Tx,Tz)
+deallocate(wTi,wTj)
+end subroutine
+
+
+subroutine fvn_zi_sparse_det(n,nz,T,Ti,Tj,det,status)
+implicit none
+integer(kind=sp_kind), intent(in) :: n,nz
+complex(kind=dp_kind),dimension(nz),intent(in) :: T
+integer(kind=sp_kind),dimension(nz),intent(in)  :: Ti,Tj
+integer(kind=sp_kind), intent(out) :: status
+complex(kind=dp_kind), intent(out) :: det
+
+integer(kind=sp_kind),dimension(:),allocatable :: wTi,wTj
+real(kind=dp_kind),dimension(:),allocatable :: Tx,Tz
+real(kind=dp_kind),dimension(:),allocatable :: Ax,Az
+integer(kind=sp_kind),dimension(:),allocatable :: Ap,Ai
+!integer(kind=dp_kind) :: symbolic,numeric
+integer(kind=sp_kind),dimension(2) :: symbolic,numeric
+! As symbolic and numeric are used to store a C pointer, it is necessary to
+! still use an integer(kind=dp_kind) for 64bits machines
+! An other possibility : integer(kind=sp_kind),dimension(2) :: symbolic,numeric
+real(kind=dp_kind),dimension(90) :: info
+real(kind=dp_kind),dimension(20) :: control
+real(kind=dp_kind) :: Mx,Mz,Ex
+
+status=0
+! we use a working copy of Ti and Tj to perform 1-based to 0-based translation
+! Tx and Tz are the real and imaginary parts of T
+allocate(wTi(nz),wTj(nz))
+allocate(Tx(nz),Tz(nz))
+Tx=dble(T)
+Tz=aimag(T)
+wTi=Ti-1
+wTj=Tj-1
+allocate(Ax(nz),Az(nz))
+allocate(Ap(n+1),Ai(nz))
+
+! perform the triplet to compressed column form -> Ap,Ai,Ax,Az
+call umfpack_zi_triplet_to_col(n,n,nz,wTi,wTj,Tx,Tz,Ap,Ai,Ax,Az,status)
+! if status is not zero a problem has occured
+if (status /= 0) then
+    write(*,*) "Problem during umfpack_zl_triplet_to_col"
+endif
+
+! Define defaults control values
+call umfpack_zi_defaults(control)
+
+! Symbolic analysis
+call umfpack_zi_symbolic(n,n,Ap,Ai,Ax,Az,symbolic, control, info)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during symbolic analysis"
+    status=info(1)
+endif
+
+! Numerical factorization
+call umfpack_zi_numeric (Ap, Ai, Ax, Az, symbolic, numeric, control, info)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during numerical factorization"
+    status=info(1)
+endif
+
+! free the C symbolic pointer
+call umfpack_zi_free_symbolic (symbolic)
+
+call umfpack_zi_get_determinant(Mx,Mz,Ex,numeric,info,status)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during sparse determinant, returned code : ",info(1)
+    status=info(1)
+endif
+det=cmplx(Mx,Mz,kind=dp_kind)*10**Ex
+
+! free the C numeric pointer
+call umfpack_zi_free_numeric (numeric)
+
+deallocate(Ax,Az)
+deallocate(Tx,Tz)
+deallocate(wTi,wTj)
+end subroutine
+
+subroutine fvn_dl_sparse_det(n,nz,T,Ti,Tj,det,status)
+implicit none
+integer(kind=dp_kind), intent(in) :: n,nz
+real(kind=dp_kind),dimension(nz),intent(in) :: T
+integer(kind=dp_kind),dimension(nz),intent(in)  :: Ti,Tj
+integer(kind=dp_kind), intent(out) :: status
+real(kind=dp_kind), intent(out) :: det
+
+integer(kind=dp_kind),dimension(:),allocatable :: wTi,wTj
+real(kind=dp_kind),dimension(:),allocatable :: A
+integer(kind=dp_kind),dimension(:),allocatable :: Ap,Ai
+!integer(kind=dp_kind) :: symbolic,numeric
+integer(kind=dp_kind) :: symbolic,numeric
+real(kind=dp_kind),dimension(90) :: info
+real(kind=dp_kind),dimension(20) :: control
+real(kind=dp_kind) :: Mx,Ex
+
+status=0
+! we use a working copy of Ti and Tj to perform 1-based to 0-based translation
+allocate(wTi(nz),wTj(nz))
+wTi=Ti-1
+wTj=Tj-1
+allocate(A(nz))
+allocate(Ap(n+1),Ai(nz))
+
+! perform the triplet to compressed column form -> Ap,Ai,Ax,Az
+call umfpack_dl_triplet_to_col(n,n,nz,wTi,wTj,T,Ap,Ai,A,status)
+! if status is not zero a problem has occured
+if (status /= 0) then
+    write(*,*) "Problem during umfpack_dl_triplet_to_col"
+endif
+
+! Define defaults control values
+call umfpack_dl_defaults(control)
+
+! Symbolic analysis
+call umfpack_dl_symbolic(n,n,Ap,Ai,A,symbolic, control, info)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during symbolic analysis"
+    status=info(1)
+endif
+
+! Numerical factorization
+call umfpack_dl_numeric (Ap, Ai, A, symbolic, numeric, control, info)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during numerical factorization"
+    status=info(1)
+endif
+
+! free the C symbolic pointer
+call umfpack_dl_free_symbolic (symbolic)
+
+call umfpack_dl_get_determinant(Mx,Ex,numeric,info,status)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during sparse determinant, returned code : ",info(1)
+    status=info(1)
+endif
+det=Mx*10**Ex
+
+! free the C numeric pointer
+call umfpack_dl_free_numeric (numeric)
+
+deallocate(A)
+deallocate(wTi,wTj)
+end subroutine
+
+
+subroutine fvn_di_sparse_det(n,nz,T,Ti,Tj,det,status)
+implicit none
+integer(kind=sp_kind), intent(in) :: n,nz
+real(kind=dp_kind),dimension(nz),intent(in) :: T
+integer(kind=sp_kind),dimension(nz),intent(in)  :: Ti,Tj
+integer(kind=sp_kind), intent(out) :: status
+real(kind=dp_kind), intent(out) :: det
+
+integer(kind=sp_kind),dimension(:),allocatable :: wTi,wTj
+real(kind=dp_kind),dimension(:),allocatable :: A
+integer(kind=sp_kind),dimension(:),allocatable :: Ap,Ai
+!integer(kind=dp_kind) :: symbolic,numeric
+integer(kind=sp_kind),dimension(2) :: symbolic,numeric
+! As symbolic and numeric are used to store a C pointer, it is necessary to
+! still use an integer(kind=dp_kind) for 64bits machines
+! An other possibility : integer(kind=sp_kind),dimension(2) :: symbolic,numeric
+real(kind=dp_kind),dimension(90) :: info
+real(kind=dp_kind),dimension(20) :: control
+integer(kind=sp_kind) :: sys
+real(kind=dp_kind) :: Mx,Ex
+
+
+status=0
+! we use a working copy of Ti and Tj to perform 1-based to 0-based translation
+allocate(wTi(nz),wTj(nz))
+wTi=Ti-1
+wTj=Tj-1
+allocate(A(nz))
+allocate(Ap(n+1),Ai(nz))
+
+! perform the triplet to compressed column form -> Ap,Ai,Ax,Az
+call umfpack_di_triplet_to_col(n,n,nz,wTi,wTj,T,Ap,Ai,A,status)
+! if status is not zero a problem has occured
+if (status /= 0) then
+    write(*,*) "Problem during umfpack_di_triplet_to_col"
+endif
+
+! Define defaults control values
+call umfpack_di_defaults(control)
+
+! Symbolic analysis
+call umfpack_di_symbolic(n,n,Ap,Ai,A,symbolic, control, info)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during symbolic analysis"
+    status=info(1)
+endif
+
+! Numerical factorization
+call umfpack_di_numeric (Ap, Ai, A, symbolic, numeric, control, info)
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during numerical factorization"
+    status=info(1)
+endif
+
+
+! free the C symbolic pointer
+call umfpack_di_free_symbolic (symbolic)
+
+call umfpack_di_get_determinant(Mx,Ex,numeric,info,status)
+det=Mx*10**Ex
+! info(1) should be zero
+if (info(1) /= 0) then
+    write(*,*) "Problem during sparse determinant, returned code : ",info(1)
+    status=info(1)
+endif
+
+! free the C numeric pointer
+call umfpack_di_free_numeric (numeric)
+
+deallocate(A)
+deallocate(wTi,wTj)
+end subroutine
+
+
+
+
+
 end module fvn_sparse
diff --git a/fvn_sparse/umfpack_wrapper.c b/fvn_sparse/umfpack_wrapper.c
index a0250b2..edcada7 100644
--- a/fvn_sparse/umfpack_wrapper.c
+++ b/fvn_sparse/umfpack_wrapper.c
@@ -72,6 +72,11 @@ void umfpack_zl_triplet_to_col_ (UF_long *m, UF_long *n, UF_long *nz, UF_long Ti
 }
 
 
+/* get determinant */
+void umfpack_zl_get_determinant_ ( double *Mx, double *Mz, double *Ex, void **Numeric, double Info [UMFPACK_INFO], UF_long *status )
+{
+    *status = umfpack_zl_get_determinant (Mx,Mz,Ex,*Numeric, Info);
+}
 
 
 /* complex(8) and integer(4) */
@@ -131,6 +136,14 @@ void umfpack_zi_triplet_to_col_ (int *m, int *n, int *nz, int Ti [ ], int Tj [ ]
     *status = umfpack_zi_triplet_to_col (*m, *n, *nz, Ti, Tj, Tx, Tz, Ap, Ai, Ax, Az, (int *) NULL);
 }
 
+
+/* get determinant */
+void umfpack_zi_get_determinant_ ( double *Mx, double *Mz, double *Ex, void **Numeric, double Info [UMFPACK_INFO], int *status )
+{
+    *status = umfpack_zi_get_determinant (Mx,Mz,Ex,*Numeric, Info);
+}
+
+
 /* real(8) and integer(8) */
 
 /* defaults  */
@@ -178,6 +191,14 @@ void umfpack_dl_triplet_to_col_ (UF_long *m, UF_long *n, UF_long *nz, UF_long Ti
     *status = umfpack_dl_triplet_to_col (*m, *n, *nz, Ti, Tj, T, Ap, Ai, A, (UF_long *) NULL);
 }
 
+/* get determinant */
+void umfpack_dl_get_determinant_ ( double *Mx, double *Ex, void **Numeric, double Info [UMFPACK_INFO], UF_long *status )
+{
+    *status = umfpack_dl_get_determinant (Mx,Ex,*Numeric, Info);
+}
+
+
+
 /* real(8) and integer(4) */
 
 /* defaults  */
@@ -225,3 +246,9 @@ void umfpack_di_triplet_to_col_ (int *m, int *n, int *nz, int Ti [ ], int Tj [ ]
     *status = umfpack_di_triplet_to_col (*m, *n, *nz, Ti, Tj, T, Ap, Ai, A, (int *) NULL);
 }
 
+/* get determinant */
+void umfpack_di_get_determinant_ ( double *Mx, double *Ex, void **Numeric, double Info [UMFPACK_INFO], int *status )
+{
+    *status = umfpack_di_get_determinant (Mx,Ex,*Numeric, Info);
+}
+
diff --git a/fvn_test/Makefile b/fvn_test/Makefile
index 19f1a31..97ea79d 100644
--- a/fvn_test/Makefile
+++ b/fvn_test/Makefile
@@ -2,10 +2,11 @@
 include $(BTREE)/Make.inc
 
 programs = test_fac$(exext) test_matinv$(exext) test_specfunc$(exext)  \
-test_det$(exext) test_matcon$(exext) test_matev$(exext) test_sparse$(exext) test_inter1d$(exext) \
+test_det$(exext) test_matcon$(exext) test_matev$(exext) test_inter1d$(exext) \
 test_inter2d$(exext) test_inter3d$(exext) test_akima$(exext) test_lsp$(exext) test_muller$(exext) \
 test_integ$(exext) test_bsyn$(exext) test_bsjn$(exext) test_bskn$(exext) test_bsin$(exext) test_operators$(exext) test_ze1$(exext) \
-test_besri$(exext) test_besrj$(exext) test_bestime$(exext)
+test_besri$(exext) test_besrj$(exext) test_bestime$(exext) \
+test_sparse_zi$(exext) test_sparse_zl$(exext) test_sparse_di$(exext) test_sparse_dl$(exext)
 
 prog:$(programs)
 
diff --git a/fvn_test/test_sparse.f90 b/fvn_test/test_sparse.f90
deleted file mode 100644
index 09014d8..0000000
--- a/fvn_test/test_sparse.f90
+++ /dev/null
@@ -1,35 +0,0 @@
-program test_sparse
-use fvn_sparse
-implicit none
-integer(kind=ip_kind), parameter :: nz=12
-integer(kind=ip_kind), parameter :: n=5
-real(kind=dp_kind),dimension(nz) :: A
-real(kind=dp_kind),dimension(n,n) :: As
-integer(kind=ip_kind),dimension(nz) :: Ti,Tj
-real(kind=dp_kind),dimension(n) :: B,x
-integer(kind=ip_kind) :: status,i
-! Description of the matrix in triplet form
-A = (/ 2.,3.,3.,-1.,4.,4.,-3.,1.,2.,2.,6.,1. /)
-B = (/ 8., 45., -3., 3., 19./)
-Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /)
-Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)
-
-! Reconstruction of the matrix in standard form
-! just needed for printing the matrix here
-As=0.
-do i=1,nz
-    As(Ti(i),Tj(i))=A(i)
-end do
-write(*,*) "Matrix in standard representation :"
-do i=1,5
-    write(*,'(5f8.4)') As(i,:)
-end do
-write(*,*)
-write(*,'("Right hand side :",5f8.4)') B
-
-!specific routine that will be used here
-!call fvn_di_sparse_solve(n,nz,A,Ti,Tj,B,x,status)
-call fvn_di_sparse_solve(n,nz,A,Ti,Tj,B,x,status)
-write(*,'("Solution :",5f8.4)') x
-write(*,'("Product matrix Solution :",5f8.4)') matmul(As,x)
-end program
\ No newline at end of file
diff --git a/fvn_test/test_sparse_di.f90 b/fvn_test/test_sparse_di.f90
new file mode 100644
index 0000000..868446b
--- /dev/null
+++ b/fvn_test/test_sparse_di.f90
@@ -0,0 +1,49 @@
+program test_sparse
+! test sparse routine di, real(8) and integer(4)
+use fvn
+implicit none
+integer(kind=sp_kind), parameter :: nz=12
+integer(kind=sp_kind), parameter :: n=5
+real(kind=dp_kind),dimension(nz) :: A
+real(kind=dp_kind),dimension(n,n) :: As
+integer(kind=sp_kind),dimension(nz) :: Ti,Tj
+real(kind=dp_kind),dimension(n) :: B,x
+integer(kind=sp_kind) :: status,i
+real(kind=dp_kind) :: det
+! Description of the matrix in triplet form
+A = (/ 2.,3.,3.,-1.,4.,4.,-3.,1.,2.,2.,6.,1. /)
+B = (/ 8., 45., -3., 3., 19./)
+Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /)
+Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)
+
+! Reconstruction of the matrix in standard form
+! just needed for printing the matrix here
+As=0.
+do i=1,nz
+    As(Ti(i),Tj(i))=A(i)
+end do
+write(*,*) "Matrix in standard representation :"
+do i=1,5
+    write(*,'(5f8.4)') As(i,:)
+end do
+write(*,*)
+write(*,*) "Standard determinant =", fvn_det(5,As)
+write(*,*)
+write(*,'("Right hand side :",5f8.4)') B
+
+! can use either specific interface, fvn_di_sparse_det
+! either generic one fvn_sparse_det
+call fvn_di_sparse_det(n,nz,A,Ti,Tj,det,status)
+write(*,*)
+write(*,*) "Sparse Det = ",det
+! can use either specific interface fvn_di_sparse_solve
+! either generic one fvn_sparse_solve
+! parameter det is optional
+call fvn_di_sparse_solve(n,nz,A,Ti,Tj,B,x,status,det)
+write(*,*)
+write(*,*) "Sparse Det as solve option = ",det
+write(*,*)
+write(*,'("Solution :",5f8.4)') x
+write(*,*)
+write(*,'("Product matrix Solution :",5f8.4)') matmul(As,x)
+end program
diff --git a/fvn_test/test_sparse_dl.f90 b/fvn_test/test_sparse_dl.f90
new file mode 100644
index 0000000..c683661
--- /dev/null
+++ b/fvn_test/test_sparse_dl.f90
@@ -0,0 +1,49 @@
+program test_sparse
+! test sparse routine dl, real(8) and integer(8)
+use fvn
+implicit none
+integer(kind=dp_kind), parameter :: nz=12
+integer(kind=dp_kind), parameter :: n=5
+real(kind=dp_kind),dimension(nz) :: A
+real(kind=dp_kind),dimension(n,n) :: As
+integer(kind=dp_kind),dimension(nz) :: Ti,Tj
+real(kind=dp_kind),dimension(n) :: B,x
+integer(kind=dp_kind) :: status,i
+real(kind=dp_kind) :: det
+! Description of the matrix in triplet form
+A = (/ 2.,3.,3.,-1.,4.,4.,-3.,1.,2.,2.,6.,1. /)
+B = (/ 8., 45., -3., 3., 19./)
+Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /)
+Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)
+
+! Reconstruction of the matrix in standard form
+! just needed for printing the matrix here
+As=0.
+do i=1,nz
+    As(Ti(i),Tj(i))=A(i)
+end do
+write(*,*) "Matrix in standard representation :"
+do i=1,5
+    write(*,'(5f8.4)') As(i,:)
+end do
+write(*,*)
+write(*,*) "Standard determinant =", fvn_det(5,As)
+write(*,*)
+write(*,'("Right hand side :",5f8.4)') B
+
+! can use either specific interface, fvn_dl_sparse_det
+! either generic one fvn_sparse_det
+call fvn_dl_sparse_det(n,nz,A,Ti,Tj,det,status)
+write(*,*)
+write(*,*) "Sparse Det = ",det
+! can use either specific interface fvn_dl_sparse_solve
+! either generic one fvn_sparse_solve
+! parameter det is optional
+call fvn_dl_sparse_solve(n,nz,A,Ti,Tj,B,x,status,det)
+write(*,*)
+write(*,*) "Sparse Det as solve option = ",det
+write(*,*)
+write(*,'("Solution :",5f8.4)') x
+write(*,*)
+write(*,'("Product matrix Solution :",5f8.4)') matmul(As,x)
+end program
diff --git a/fvn_test/test_sparse_zi.f90 b/fvn_test/test_sparse_zi.f90
new file mode 100644
index 0000000..1d5ac67
--- /dev/null
+++ b/fvn_test/test_sparse_zi.f90
@@ -0,0 +1,55 @@
+program test_sparse
+use fvn
+implicit none
+integer(kind=sp_kind), parameter :: nz=12
+integer(kind=sp_kind), parameter :: n=5
+complex(kind=dp_kind),dimension(nz) :: A
+complex(kind=dp_kind),dimension(n,n) :: As
+integer(kind=sp_kind),dimension(nz) :: Ti,Tj
+complex(kind=dp_kind),dimension(n) :: B,x
+integer(kind=sp_kind) :: status,i
+complex(kind=dp_kind) :: det
+character(len=80) :: fmcmplx
+
+fmcmplx='(5("(",f8.5,",",f8.5,")  "))'
+
+! Description of the matrix in triplet form
+A = (/ (2.,-1.),(3.,2.),(3.,1.),(-1.,5.),(4.,-7.),(4.,0.),(-3.,-4.),(1.,3.),(2.,0.),(2.,-2.),(6.,4.),(1.,0.) /)
+B = (/ (8.,3.), (45.,1.), (-3.,-2.), (3.,0.), (19.,2.) /)
+Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /)
+Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)
+
+! Reconstruction of the matrix in standard form
+As=0.
+do i=1,nz
+    As(Ti(i),Tj(i))=A(i)
+end do
+
+write(*,*) "Matrix in standard representation :"
+do i=1,5
+    write(*,fmcmplx) As(i,:)
+end do
+write(*,*)
+write(*,*) "Standard determinant : ",fvn_det(5,As)
+write(*,*)
+write(*,*) "Right hand side :"
+write(*,fmcmplx) B
+
+! can use either specific interface, fvn_zi_sparse_det
+! either generic one fvn_sparse_det
+call fvn_zi_sparse_det(n,nz,A,Ti,Tj,det,status)
+write(*,*)
+write(*,*) "Sparse Det = ",det
+! can use either specific interface fvn_zi_sparse_solve
+! either generic one fvn_sparse_solve
+! parameter det is optional
+call fvn_zi_sparse_solve(n,nz,A,Ti,Tj,B,x,status,det)
+write(*,*)
+write(*,*) "Sparse Det as solve option= ",det
+write(*,*)
+write(*,*) "Solution :"
+write(*,fmcmplx) x
+write(*,*)
+write(*,*) "Product matrix Solution :"
+write(*,fmcmplx) matmul(As,x)
+end program
diff --git a/fvn_test/test_sparse_zl.f90 b/fvn_test/test_sparse_zl.f90
new file mode 100644
index 0000000..24e7867
--- /dev/null
+++ b/fvn_test/test_sparse_zl.f90
@@ -0,0 +1,55 @@
+program test_sparse
+! test sparse routine zl, complex(8) and integer(8)
+use fvn
+implicit none
+integer(kind=dp_kind), parameter :: nz=12
+integer(kind=dp_kind), parameter :: n=5
+complex(kind=dp_kind),dimension(nz) :: A
+complex(kind=dp_kind),dimension(n,n) :: As
+integer(kind=dp_kind),dimension(nz) :: Ti,Tj
+complex(kind=dp_kind),dimension(n) :: B,x
+integer(kind=dp_kind) :: status,i
+complex(kind=dp_kind) :: det
+character(len=80) :: fmcmplx
+
+fmcmplx='(5("(",f8.5,",",f8.5,")  "))'
+
+! Description of the matrix in triplet form
+A = (/ (2.,-1.),(3.,2.),(3.,1.),(-1.,5.),(4.,-7.),(4.,0.),(-3.,-4.),(1.,3.),(2.,0.),(2.,-2.),(6.,4.),(1.,0.) /)
+B = (/ (8.,3.), (45.,1.), (-3.,-2.), (3.,0.), (19.,2.) /)
+Ti = (/ 1,2,1,3,5,2,3,4,5,3,2,5 /)
+Tj = (/ 1,1,2,2,2,3,3,3,3,4,5,5 /)
+
+! Reconstruction of the matrix in standard form
+As=0.
+do i=1,nz
+    As(Ti(i),Tj(i))=A(i)
+end do
+
+write(*,*) "Matrix in standard representation :"
+do i=1,5
+    write(*,fmcmplx) As(i,:)
+end do
+write(*,*)
+write(*,*) "Standard determinant : ",fvn_det(5,As)
+write(*,*)
+write(*,*) "Right hand side :"
+write(*,fmcmplx) B
+! can use either specific interface, fvn_zl_sparse_det
+! either generic one fvn_sparse_det
+call fvn_zl_sparse_det(n,nz,A,Ti,Tj,det,status)
+write(*,*)
+write(*,*) "Sparse Det = ",det
+! can use either specific interface fvn_zl_sparse_solve
+! either generic one fvn_sparse_solve
+! parameter det is optional
+call fvn_zl_sparse_solve(n,nz,A,Ti,Tj,B,x,status,det)
+write(*,*)
+write(*,*) "Sparse Det as solve option= ",det
+write(*,*)
+write(*,*) "Solution :"
+write(*,fmcmplx) x
+write(*,*)
+write(*,*) "Product matrix Solution :"
+write(*,fmcmplx) matmul(As,x)
+end program