@@ -13,7 +13,6 @@ module test_linalg_least_squares
1313 public :: test_least_squares
1414
1515 contains
16-
1716
1817 !> Solve sample least squares problems
1918 subroutine test_least_squares(tests)
@@ -24,15 +23,16 @@ module test_linalg_least_squares
2423
2524 #:for rk,rt,ri in REAL_KINDS_TYPES
2625 #:if rk!="xdp"
27- tests = [tests,new_unittest("lease_squares_${ri}$",test_lstsq_one_${ri}$)]
26+ tests = [tests,new_unittest("least_squares_${ri}$",test_lstsq_one_${ri}$), &
27+ new_unittest("least_squares_randm_${ri}$",test_lstsq_random_${ri}$)]
2828 #:endif
2929 #:endfor
3030
3131 end subroutine test_least_squares
32-
33- !> Simple polynomial fit
32+
3433 #:for rk,rt,ri in REAL_KINDS_TYPES
3534 #:if rk!="xdp"
35+ !> Simple polynomial fit
3636 subroutine test_lstsq_one_${ri}$(error)
3737 type(error_type), allocatable, intent(out) :: error
3838
@@ -64,6 +64,32 @@ module test_linalg_least_squares
6464
6565 end subroutine test_lstsq_one_${ri}$
6666
67+ !> Fit from random array
68+ subroutine test_lstsq_random_${ri}$(error)
69+ type(error_type), allocatable, intent(out) :: error
70+
71+ type(linalg_state_type) :: state
72+ integer(ilp), parameter :: n = 12, m = 3
73+ ${rt}$ :: xsol(m),x(m),y(n),A(n,m)
74+
75+ ! Random coefficient matrix and solution
76+ call random_number(A)
77+ call random_number(xsol)
78+
79+ ! Compute rhs
80+ y = matmul(A,xsol)
81+
82+ ! Find polynomial
83+ x = lstsq(A,y,err=state)
84+
85+ call check(error,state%ok(),state%print())
86+ if (allocated(error)) return
87+
88+ call check(error, all(abs(x-xsol)<1.0e-6_${rk}$), 'data converged')
89+ if (allocated(error)) return
90+
91+ end subroutine test_lstsq_random_${ri}$
92+
6793 #:endif
6894 #:endfor
6995
0 commit comments