@@ -110,22 +110,25 @@ module test_linalg_least_squares
110
110
! Dimension of the problem.
111
111
integer(ilp), parameter :: n = 42
112
112
! Data for the least-squares problem.
113
- complex(dp) :: A(n+1, n), b(n+1), x_true(n), x_lstsq(n+1 )
113
+ complex(dp) :: A(n+1, n), b(n+1), x_true(n), x_lstsq(n)
114
114
! Internal variables.
115
115
real(dp), allocatable :: tmp(:, :, :), tmp_vec(:, :)
116
116
! Error handler
117
117
type(linalg_state_type) :: state
118
118
119
119
! Zero-out data.
120
- A = 0.0_dp ; b = 0.0_dp ; x_true = 0.0_dp ; x_lstsq = 0.0_dp
121
- allocate(tmp(n+1, n, 2)) ; tmp = 0.0_dp
122
- allocate(tmp_vec(n, 2)) ; tmp_vec = 0.0_dp
120
+ A = 0.0_dp
121
+ b = 0.0_dp
122
+ x_lstsq = 0.0_dp
123
+ allocate(tmp(n+1, n, 2), tmp_vec(n, 2), source=0.0_dp)
123
124
124
125
! Generate a random complex least-squares problem of size (n+1, n).
125
- call random_number(tmp) ; call random_number(tmp_vec)
126
- A = cmplx(tmp(:, :, 1), tmp(:, :, 2), kind=dp)
126
+ call random_number(tmp)
127
+ call random_number(tmp_vec)
128
+
129
+ A = cmplx(tmp(:, :, 1), tmp(:, :, 2), kind=dp)
127
130
x_true = cmplx(tmp_vec(:, 1), tmp_vec(:, 2), kind=dp)
128
- b = matmul(A, x_true)
131
+ b = matmul(A, x_true)
129
132
130
133
! Solve the lstsq problem.
131
134
call solve_lstsq(A, b, x_lstsq, err=state)
@@ -134,6 +137,10 @@ module test_linalg_least_squares
134
137
call check(error,state%ok(),'issue 823 returned '//state%print())
135
138
if (allocated(error)) return
136
139
140
+ ! Check that least squares are verified
141
+ call check(error,all(abs(x_true-x_lstsq)<sqrt(epsilon(0.0_dp))),'issue 823 results')
142
+ if (allocated(error)) return
143
+
137
144
end subroutine test_issue_823
138
145
139
146
end module test_linalg_least_squares
0 commit comments