First DLRA and LAPACK test

Makefile

@@ -93,7 +93,7 @@ $(OBJDIR)/moments_eq_rhs_mod.o $(OBJDIR)/numerics_mod.o $(OBJDIR)/parallel_mod.o
 $(OBJDIR)/ppexit.o $(OBJDIR)/prec_const_mod.o \
 $(OBJDIR)/processing_mod.o $(OBJDIR)/readinputs.o $(OBJDIR)/restarts_mod.o \
 $(OBJDIR)/solve_EM_fields.o $(OBJDIR)/species_mod.o $(OBJDIR)/stepon.o $(OBJDIR)/tesend.o \
-$(OBJDIR)/time_integration_mod.o $(OBJDIR)/utility_mod.o
+$(OBJDIR)/time_integration_mod.o $(OBJDIR)/utility_mod.o $(OBJDIR)/DLRA_mod.o
 
 # To compile the executable
  compile: $(FOBJ)
@@ -282,7 +282,8 @@ $(OBJDIR)/time_integration_mod.o $(OBJDIR)/utility_mod.o
    $(OBJDIR)/basic_mod.o $(OBJDIR)/nonlinear_mod.o $(OBJDIR)/grid_mod.o \
 	 $(OBJDIR)/array_mod.o $(OBJDIR)/numerics_mod.o $(OBJDIR)/fields_mod.o \
 	 $(OBJDIR)/ghosts_mod.o $(OBJDIR)/moments_eq_rhs_mod.o $(OBJDIR)/solve_EM_fields.o\
-	 $(OBJDIR)/utility_mod.o $(OBJDIR)/model_mod.o $(OBJDIR)/time_integration_mod.o
+	 $(OBJDIR)/utility_mod.o $(OBJDIR)/model_mod.o $(OBJDIR)/time_integration_mod.o \
+	 $(OBJDIR)/DLRA_mod.o
 	$(F90) -c $(F90FLAGS) $(FPPFLAGS) $(EXTMOD) $(EXTINC) src/stepon.F90 -o $@
 
  $(OBJDIR)/tesend.o : src/tesend.F90 \
@@ -296,3 +297,7 @@ $(OBJDIR)/time_integration_mod.o $(OBJDIR)/utility_mod.o
  $(OBJDIR)/utility_mod.o : src/utility_mod.F90  \
    $(OBJDIR)/grid_mod.o $(OBJDIR)/basic_mod.o $(OBJDIR)/prec_const_mod.o
 	$(F90) -c $(F90FLAGS) $(FPPFLAGS) $(EXTMOD) $(EXTINC) src/utility_mod.F90 -o $@
+
+ $(OBJDIR)/DLRA_mod.o : src/DLRA_mod.F90  \
+   $(OBJDIR)/basic_mod.o $(OBJDIR)/prec_const_mod.o
+	$(F90) -c $(F90FLAGS) $(FPPFLAGS) $(EXTMOD) $(EXTINC) src/DLRA_mod.F90 -o $@

local/make.inc

@@ -78,10 +78,10 @@ ifdef FFTW3DIR
 endif
 
 # Add lapack local lib
-ifdef LAPACKDIR
-      LIBS += -llapack -lblas
-      LDIRS += -L$(LAPACKDIR)/lib
-endif
+#ifdef LAPACKDIR
+      LIBS  += -llapack -lblas
+#      LDIRS += -L$(LAPACKDIR)/lib
+#endif
 
 # FM library
 ifdef FMDIR

src/DLRA_mod.F90

 module DLRA
     USE prec_const
-    USE lapack
     implicit none
 
-    PUBLIC :: filter_singular_value_ky_pj
+    PUBLIC :: filter_singular_value_ky_pj, test_SVD
 
     CONTAINS
 
@@ -17,4 +16,79 @@ module DLRA
         ! Singular value decomposition
         ! CALL SVD(array_ky_pj,singular_values)
     END SUBROUTINE
+
+    SUBROUTINE test_SVD
+        ! Program to perform Singular Value Decomposition (SVD)
+        ! using LAPACK library
+
+        ! Specify the dimensions of the input matrix A
+        INTEGER, PARAMETER :: m = 3, n = 2
+        INTEGER, PARAMETER :: lda = m
+
+        ! Declare the input matrix A
+        REAL, DIMENSION(lda,n) :: A
+
+        ! Specify the dimensions of the output matrices
+        INTEGER, PARAMETER :: ldu = m, ldvt = n
+        INTEGER, PARAMETER :: lwork = 5*n
+        REAL, DIMENSION(ldu,m) :: U
+        REAL, DIMENSION(n) :: S
+        REAL, DIMENSION(ldvt,n) :: VT
+        REAL, DIMENSION(lwork) :: work
+        INTEGER :: info, i,j
+
+        ! Define the input matrix A
+        A = RESHAPE((/ 1.0, 2.0, 3.0, 4.0, 5.0, 6.0 /), SHAPE(A))
+
+        WRITE(*,*) 'Input matrix A = '
+
+        DO i = 1, m
+            WRITE(*,'(2X, 3F8.3)') (A(i,j), j=1,n)
+        END DO
+
+        ! Compute the SVD of A using the LAPACK subroutine SGESVD
+        CALL SGESVD('A', 'A', m, n, A, lda, S, U, ldu, VT, ldvt, work, lwork, info)
+
+        ! Print the results
+        WRITE(*,*) 'U = '
+        DO i = 1, m
+            WRITE(*,'(6F8.3)') (U(i,j), j=1,m)
+        END DO
+        WRITE(*,*)
+        WRITE(*,*) 'S = '
+        WRITE(*,'(2F8.3)') (S(i), i=1,n)
+        WRITE(*,*)
+        WRITE(*,*) 'VT = '
+        DO i = 1, n
+            WRITE(*,'(6F8.3)') (VT(i,j), j=1,n)
+        END DO
+
+        ! Reconstruct A from its SVD
+        A = MATMUL(U, MATMUL(diagmat(S,m,n), TRANSPOSE(VT)))
+        ! Print the reconstructed matrix A
+
+        WRITE(*,*) 'Reconstructed matrix A = '
+
+        DO i = 1, m
+            WRITE(*,'(2X, 3F8.3)') (A(i,j), j=1,n)
+        END DO
+
+        print*, "this was a test of the SVD using LAPACK. End run."
+        stop
+    END SUBROUTINE test_SVD
+
+    FUNCTION diagmat(v, m, n) RESULT(A)
+        REAL, DIMENSION(:), INTENT(IN) :: v
+        INTEGER, INTENT(IN) :: m, n
+        REAL, DIMENSION(m,n) :: A
+        INTEGER :: i, j
+      
+        A = 0.0   ! Initialize A to a zero matrix
+      
+        DO i = 1, MIN(m,n)
+          A(i,i) = v(i)   ! Set the diagonal elements of A to the values in v
+        END DO
+      
+      END FUNCTION diagmat
+
 end module DLRA
\ No newline at end of file

src/stepon.F90

@@ -8,6 +8,7 @@ SUBROUTINE stepon
    use mpi,                   ONLY: MPI_COMM_WORLD
    USE time_integration,      ONLY: ntimelevel
    USE prec_const,            ONLY: xp
+   USE DLRA,                  ONLY: test_SVD
    IMPLICIT NONE
 
    INTEGER :: num_step, ierr
@@ -56,6 +57,7 @@ SUBROUTINE stepon
 
       !! TEST SINGULAR VALUE DECOMPOSITION
       ! CALL filter_singular_value_ky_pj(nsv,moments)
+      ! CALL test_SVD
       
       IF( nlend ) EXIT ! exit do loop