chol: function interface, implement and test

Author: perazz

src/stdlib_linalg.fypp

@@ -16,6 +16,7 @@ module stdlib_linalg
   implicit none
   private
 
+  public :: chol
   public :: cholesky
   public :: det
   public :: operator(.det.)
@@ -43,6 +44,23 @@ module stdlib_linalg
   ! Export linalg error handling
   public :: linalg_state_type, linalg_error_handling
 
+  interface chol
+     #:for rk,rt,ri in RC_KINDS_TYPES
+     #:if rk!="xdp"
+     pure module function stdlib_linalg_${ri}$_cholesky_fun(a,lower,other_zeroed) result(c)
+         !> Input matrix a[m,n]
+         ${rt}$, intent(in) :: a(:,:)
+         !> [optional] is the lower or upper triangular factor required? Default = lower
+         logical(lk), optional, intent(in) :: lower
+         !> [optional] should the unused half of the return matrix be zeroed out? Default: yes
+         logical(lk), optional, intent(in) :: other_zeroed
+         !> Output matrix with Cholesky factors c[n,n]
+         ${rt}$, allocatable :: c(:,:)                  
+     end function stdlib_linalg_${ri}$_cholesky_fun
+     #:endif
+     #:endfor    
+  end interface chol  
+
   interface cholesky
      #:for rk,rt,ri in RC_KINDS_TYPES
      #:if rk!="xdp"

src/stdlib_linalg_cholesky.fypp

@@ -156,6 +156,24 @@ submodule (stdlib_linalg) stdlib_linalg_cholesky
 
      end subroutine stdlib_linalg_${ri}$_cholesky
 
+     ! Compute the Cholesky factorization of a symmetric / Hermitian matrix, A = L*L^T = U^T*U. 
+     ! Function interface
+     pure module function stdlib_linalg_${ri}$_cholesky_fun(a,lower,other_zeroed) result(c)
+         !> Input matrix a[n,n]
+         ${rt}$, intent(in) :: a(:,:)
+         !> Output matrix with Cholesky factors c[n,n]
+         ${rt}$, allocatable :: c(:,:)         
+         !> [optional] is the lower or upper triangular factor required? Default = lower
+         logical(lk), optional, intent(in) :: lower
+         !> [optional] should the unused half of the return matrix be zeroed out? Default: yes
+         logical(lk), optional, intent(in) :: other_zeroed
+         
+         allocate(c,source=a)
+         
+         call stdlib_linalg_${ri}$_cholesky_inplace(c,lower,other_zeroed)
+
+     end function stdlib_linalg_${ri}$_cholesky_fun
+
      #:endif
      #:endfor
 

test/linalg/test_linalg_cholesky.fypp

@@ -4,13 +4,13 @@
 module test_linalg_cholesky
     use testdrive, only: error_type, check, new_unittest, unittest_type
     use stdlib_linalg_constants
-    use stdlib_linalg, only: cholesky
+    use stdlib_linalg, only: cholesky,chol
     use stdlib_linalg_state, only: linalg_state_type
 
     implicit none (type,external)
     private
 
-    public :: test_cholesky_factorization
+    public :: test_cholesky_factorization 
 
     contains
 
@@ -53,12 +53,17 @@ module test_linalg_cholesky
         ! 1) Cholesky factorization with full matrices
         call cholesky(a, l, other_zeroed=.true., err=state)
 
-        call check(error,state%ok(),state%print())
+        call check(error,state%ok(),'cholesky (subr) :: '//state%print())
         if (allocated(error)) return
 
-        call check(error, all(abs(a-matmul(l,transpose(l)))<tol), 'data converged')
+        call check(error, all(abs(a-matmul(l,transpose(l)))<tol), 'cholesky (subr) :: data converged')
         if (allocated(error)) return        
 
+        ! 2) Function interface 
+        l = chol(a, other_zeroed=.true.)
+        
+        call check(error, all(abs(a-matmul(l,transpose(l)))<tol), 'cholesky (function) :: data converged')
+        if (allocated(error)) return           
 
     end subroutine test_cholesky_${ri}$