Skip to content

Commit

Permalink
Merge branch 'fortran-lang:master' into sparse
Browse files Browse the repository at this point in the history
  • Loading branch information
jalvesz authored May 28, 2024
2 parents 1d9dabc + dc707f8 commit 8278f38
Show file tree
Hide file tree
Showing 2 changed files with 45 additions and 3 deletions.
4 changes: 2 additions & 2 deletions src/stdlib_linalg_least_squares.fypp
Original file line number Diff line number Diff line change
Expand Up @@ -305,7 +305,7 @@ submodule (stdlib_linalg) stdlib_linalg_least_squares
else
allocate(cwork(lcwork))
endif
ncs = size(iwork,kind=ilp)
ncs = size(cwork,kind=ilp)
#:endif

if (nrs<lrwork .or. nis<liwork#{if rt.startswith('c')}# .or. ncs<lcwork#{endif}# &
Expand All @@ -322,7 +322,7 @@ submodule (stdlib_linalg) stdlib_linalg_least_squares
! Solve system using singular value decomposition
call gelsd(m,n,nrhs,amat,lda,xmat,ldb,singular,rcond,arank, &
#:if rt.startswith('complex')
cwork,nrs,rwork,iwork,info)
cwork,ncs,rwork,iwork,info)
#:else
rwork,nrs,iwork,info)
#:endif
Expand Down
44 changes: 43 additions & 1 deletion test/linalg/test_linalg_lstsq.fypp
Original file line number Diff line number Diff line change
Expand Up @@ -4,7 +4,7 @@
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, only: lstsq,solve_lstsq
use stdlib_linalg_state, only: linalg_state_type

implicit none (type,external)
Expand All @@ -20,6 +20,8 @@ module test_linalg_least_squares
type(unittest_type), allocatable, intent(out) :: tests(:)

allocate(tests(0))

tests = [tests,new_unittest("issue_823",test_issue_823)]

#:for rk,rt,ri in REAL_KINDS_TYPES
#:if rk!="xdp"
Expand Down Expand Up @@ -100,6 +102,46 @@ module test_linalg_least_squares

#:endif
#:endfor

! Test issue #823
subroutine test_issue_823(error)
type(error_type), allocatable, intent(out) :: error

! Dimension of the problem.
integer(ilp), parameter :: n = 42
! Data for the least-squares problem.
complex(dp) :: A(n+1, n), b(n+1), x_true(n), x_lstsq(n)
! Internal variables.
real(dp), allocatable :: tmp(:, :, :), tmp_vec(:, :)
! Error handler
type(linalg_state_type) :: state

! Zero-out data.
A = 0.0_dp
b = 0.0_dp
x_lstsq = 0.0_dp
allocate(tmp(n+1, n, 2), tmp_vec(n, 2), source=0.0_dp)

! Generate a random complex least-squares problem of size (n+1, n).
call random_number(tmp)
call random_number(tmp_vec)

A = cmplx(tmp(:, :, 1), tmp(:, :, 2), kind=dp)
x_true = cmplx(tmp_vec(:, 1), tmp_vec(:, 2), kind=dp)
b = matmul(A, x_true)

! Solve the lstsq problem.
call solve_lstsq(A, b, x_lstsq, err=state)

! Check that no segfault occurred
call check(error,state%ok(),'issue 823 returned '//state%print())
if (allocated(error)) return

! Check that least squares are verified
call check(error,all(abs(x_true-x_lstsq)<sqrt(epsilon(0.0_dp))),'issue 823 results')
if (allocated(error)) return

end subroutine test_issue_823

end module test_linalg_least_squares

Expand Down

0 comments on commit 8278f38

Please sign in to comment.