Skip to content

Commit 80a1df1

Browse files
perazzperazz
committed
add tests
1 parent 14f771b commit 80a1df1

File tree

2 files changed

+196
-1
lines changed

2 files changed

+196
-1
lines changed

‎test/linalg/CMakeLists.txt

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -10,6 +10,7 @@ set(
1010
"test_linalg_norm.fypp"
1111
"test_linalg_determinant.fypp"
1212
"test_linalg_qr.fypp"
13+
"test_linalg_schur.fypp"
1314
"test_linalg_svd.fypp"
1415
"test_linalg_matrix_property_checks.fypp"
1516
"test_linalg_sparse.fypp"
@@ -26,6 +27,7 @@ ADDTEST(linalg_norm)
2627
ADDTEST(linalg_solve)
2728
ADDTEST(linalg_lstsq)
2829
ADDTEST(linalg_qr)
30+
ADDTEST(linalg_schur)
2931
ADDTEST(linalg_svd)
3032
ADDTEST(blas_lapack)
31-
ADDTEST(linalg_sparse)
33+
ADDTEST(linalg_sparse)

‎test/linalg/test_linalg_schur.fypp

Lines changed: 193 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,193 @@
1+
#:include "common.fypp"
2+
#:set RC_KINDS_TYPES = REAL_KINDS_TYPES + CMPLX_KINDS_TYPES
3+
! Test Schur decomposition
4+
module test_linalg_schur
5+
use testdrive, only: error_type, check, new_unittest, unittest_type
6+
use stdlib_linalg_constants
7+
use stdlib_linalg_state, only: LINALG_VALUE_ERROR,linalg_state_type
8+
use stdlib_linalg, only: schur,schur_space
9+
use ieee_arithmetic, only: ieee_value,ieee_quiet_nan
10+
11+
implicit none (type,external)
12+
13+
public :: test_schur_decomposition
14+
15+
contains
16+
17+
!> schur decomposition tests
18+
subroutine test_schur_decomposition(tests)
19+
!> Collection of tests
20+
type(unittest_type), allocatable, intent(out) :: tests(:)
21+
22+
allocate(tests(0))
23+
24+
#:for rk,rt,ri in RC_KINDS_TYPES
25+
tests = [tests,new_unittest("schur_api_${ri}$",test_schur_api_${ri}$), &
26+
new_unittest("schur_random_${ri}$",test_schur_random_${ri}$)]
27+
#:endfor
28+
29+
end subroutine test_schur_decomposition
30+
31+
!> schur decomposition of a random matrix
32+
#:for rk,rt,ri in RC_KINDS_TYPES
33+
subroutine test_schur_api_${ri}$(error)
34+
type(error_type), allocatable, intent(out) :: error
35+
36+
integer(ilp), parameter :: n = 15_ilp
37+
integer(ilp) :: lwork
38+
real(${rk}$), parameter :: tol = 10*sqrt(epsilon(0.0_${rk}$))
39+
complex(${rk}$) :: eigs(n)
40+
${rt}$, dimension(n,n) :: a,t,z
41+
${rt}$, allocatable :: storage(:)
42+
#:if 'complex' in rt
43+
real(${rk}$) :: rea(n,n),ima(n,n)
44+
#:endif
45+
type(linalg_state_type) :: state
46+
47+
#:if 'complex' in rt
48+
call random_number(rea)
49+
call random_number(ima)
50+
a = cmplx(rea,ima,kind=${rk}$)
51+
#:else
52+
call random_number(a)
53+
#:endif
54+
55+
! Test simple API
56+
call schur(a,t,err=state)
57+
call check(error,state%ok(),state%print())
58+
if (allocated(error)) return
59+
60+
! Test output transformation matrix
61+
call schur(a,t,z,err=state)
62+
call check(error,state%ok(),state%print())
63+
if (allocated(error)) return
64+
65+
! Test output eigenvalues
66+
call schur(a,t,eigvals=eigs,err=state)
67+
call check(error,state%ok(),state%print())
68+
if (allocated(error)) return
69+
70+
! Test storage query
71+
call schur_space(a,lwork,err=state)
72+
call check(error,state%ok(),state%print())
73+
if (allocated(error)) return
74+
75+
! Test with user-defined storage
76+
allocate(storage(lwork))
77+
call schur(a,t,eigvals=eigs,storage=storage,err=state)
78+
call check(error,state%ok(),state%print())
79+
if (allocated(error)) return
80+
81+
end subroutine test_schur_api_${ri}$
82+
83+
subroutine test_schur_random_${ri}$(error)
84+
type(error_type), allocatable, intent(out) :: error
85+
86+
integer(ilp), parameter :: n = 3_ilp
87+
real(${rk}$), parameter :: rtol = 1.0e-4_${rk}$
88+
real(${rk}$), parameter :: eps = sqrt(epsilon(0.0_${rk}$))
89+
integer(ilp) :: lwork
90+
${rt}$, allocatable :: storage(:)
91+
${rt}$, dimension(n,n) :: a,t,z,aorig
92+
#:if 'complex' in rt
93+
real(${rk}$), dimension(n,n) :: a_re,a_im
94+
#:endif
95+
type(linalg_state_type) :: state
96+
97+
#:if 'complex' in rt
98+
call random_number(a_re)
99+
call random_number(a_im)
100+
a = cmplx(a_re,a_im,kind=${rk}$)
101+
#:else
102+
call random_number(a)
103+
#:endif
104+
aorig = a
105+
106+
! 1) Run schur (standard)
107+
call schur(a,t,z,err=state)
108+
109+
! Check return code
110+
call check(error,state%ok(),state%print())
111+
if (allocated(error)) return
112+
113+
! Check solution
114+
call check(error, all(schur_error(a,z,t)<=max(rtol*abs(a),eps)), &
115+
'converged solution (${rt}$)')
116+
if (allocated(error)) return
117+
118+
! 2) Run schur (overwrite A)
119+
call schur(a,t,z,overwrite_a=.true.,err=state)
120+
121+
! Check return code
122+
call check(error,state%ok(),state%print())
123+
if (allocated(error)) return
124+
125+
! Check solution
126+
call check(error, all(schur_error(aorig,z,t)<=max(rtol*abs(aorig),eps)), &
127+
'converged solution (${rt}$ - overwrite A)')
128+
if (allocated(error)) return
129+
130+
! 3) Use working storage
131+
a = aorig
132+
call schur_space(a,lwork,err=state)
133+
134+
! Check return code
135+
call check(error,state%ok(),state%print())
136+
if (allocated(error)) return
137+
138+
allocate(storage(lwork))
139+
call schur(a,t,z,storage=storage,err=state)
140+
141+
! Check return code
142+
call check(error,state%ok(),state%print())
143+
if (allocated(error)) return
144+
145+
! Check solution
146+
call check(error, all(schur_error(a,z,t)<=max(rtol*abs(a),eps)), &
147+
'converged solution (${rt}$ - external storage)')
148+
if (allocated(error)) return
149+
150+
contains
151+
152+
pure function schur_error(a,z,t) result(err)
153+
${rt}$, intent(in), dimension(:,:) :: a,z,t
154+
real(${rk}$), dimension(size(a,1),size(a,2)) :: err
155+
156+
#:if rt.startswith('real')
157+
err = abs(matmul(matmul(z,t),transpose(z)) - a)
158+
#:else
159+
err = abs(matmul(matmul(z,t),conjg(transpose(z))) - a)
160+
#:endif
161+
end function schur_error
162+
163+
end subroutine test_schur_random_${ri}$
164+
165+
#:endfor
166+
167+
end module test_linalg_schur
168+
169+
program test_schur
170+
use, intrinsic :: iso_fortran_env, only : error_unit
171+
use testdrive, only : run_testsuite, new_testsuite, testsuite_type
172+
use test_linalg_schur, only : test_schur_decomposition
173+
implicit none
174+
integer :: stat, is
175+
type(testsuite_type), allocatable :: testsuites(:)
176+
character(len=*), parameter :: fmt = '("#", *(1x, a))'
177+
178+
stat = 0
179+
180+
testsuites = [ &
181+
new_testsuite("linalg_schur", test_schur_decomposition) &
182+
]
183+
184+
do is = 1, size(testsuites)
185+
write(error_unit, fmt) "Testing:", testsuites(is)%name
186+
call run_testsuite(testsuites(is)%collect, error_unit, stat)
187+
end do
188+
189+
if (stat > 0) then
190+
write(error_unit, '(i0, 1x, a)') stat, "test(s) failed!"
191+
error stop
192+
end if
193+
end program test_schur

0 commit comments

Comments
 (0)