add test

Author: perazz

test/linalg/CMakeLists.txt

@@ -2,10 +2,12 @@ set(
   fppFiles
   "test_linalg.fypp"
   "test_blas_lapack.fypp"
+  "test_linalg_lstsq.fypp"
   "test_linalg_matrix_property_checks.fypp"
 )
 fypp_f90("${fyppFlags}" "${fppFiles}" outFiles)
 
 ADDTEST(linalg)
 ADDTEST(linalg_matrix_property_checks)
+ADDTEST(linalg_lstsq)
 ADDTEST(blas_lapack)

test/linalg/test_linalg_lstsq.fypp

@@ -0,0 +1,96 @@
+#:include "common.fypp"
+#:set RC_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES
+! Test least squares solver
+module test_linalg_least_squares
+    use testdrive, only: error_type, check, new_unittest, unittest_type
+    use stdlib_linalg_constants
+    use stdlib_linalg, only: lstsq
+    use stdlib_linalg_state, only: linalg_state_type
+
+    implicit none (type,external)
+    private
+    
+    public :: test_least_squares 
+
+    contains
+    
+
+    !> Solve sample least squares problems
+    subroutine test_least_squares(tests)
+        !> Collection of tests
+        type(unittest_type), allocatable, intent(out) :: tests(:)
+        
+        allocate(tests(0))
+
+        #:for rk,rt,ri in REAL_KINDS_TYPES
+        #:if rk!="xdp"
+        tests = [tests,new_unittest("lease_squares_${ri}$",test_lstsq_one_${ri}$)]
+        #:endif
+        #:endfor
+
+    end subroutine test_least_squares
+
+    !> Simple polynomial fit
+    #:for rk,rt,ri in REAL_KINDS_TYPES
+    #:if rk!="xdp"
+    subroutine test_lstsq_one_${ri}$(error)
+        type(error_type), allocatable, intent(out) :: error
+
+        type(linalg_state_type) :: state
+        integer(ilp) :: rank
+
+        !> Example scattered data
+        ${rt}$, parameter :: x(*)  = real([1.0, 2.5, 3.5, 4.0, 5.0, 7.0, 8.5], ${rk}$)
+        ${rt}$, parameter :: y(*)  = real([0.3, 1.1, 1.5, 2.0, 3.2, 6.6, 8.6], ${rk}$)
+        ${rt}$, parameter :: ab(*) = real([0.20925829,  0.12013861], ${rk}$)
+
+        ${rt}$ :: M(size(x),2),p(2)
+
+        ! Coefficient matrix for polynomial y = a + b*x**2
+        M(:,1) = x**0
+        M(:,2) = x**2
+
+        ! Find polynomial
+        p = lstsq(M,y,rank=rank,err=state)
+
+        call check(error,state%ok(),state%print())
+        if (allocated(error)) return
+        
+        call check(error, all(abs(p-ab)<1.0e-6_${rk}$), 'data converged')
+        if (allocated(error)) return
+        
+        call check(error, rank==2, 'matrix rank == 2')
+        if (allocated(error)) return
+
+    end subroutine test_lstsq_one_${ri}$
+    
+    #:endif
+    #:endfor
+
+end module test_linalg_least_squares
+
+program test_lstsq
+    use, intrinsic :: iso_fortran_env, only : error_unit
+    use testdrive, only : run_testsuite, new_testsuite, testsuite_type
+    use test_linalg_least_squares, only : test_least_squares
+    implicit none
+    integer :: stat, is
+    type(testsuite_type), allocatable :: testsuites(:)
+    character(len=*), parameter :: fmt = '("#", *(1x, a))'
+
+    stat = 0
+
+    testsuites = [ &
+        new_testsuite("linalg_least_squares", test_least_squares) &
+        ]
+
+    do is = 1, size(testsuites)
+        write(error_unit, fmt) "Testing:", testsuites(is)%name
+        call run_testsuite(testsuites(is)%collect, error_unit, stat)
+    end do
+
+    if (stat > 0) then
+        write(error_unit, '(i0, 1x, a)') stat, "test(s) failed!"
+        error stop
+    end if
+end program test_lstsq