4
4
module test_linalg_inverse
5
5
use testdrive, only: error_type, check, new_unittest, unittest_type
6
6
use stdlib_linalg_constants
7
- use stdlib_linalg, only: inv,invert,operator(.inv.)
8
- use stdlib_linalg_state, only: linalg_state_type
7
+ use stdlib_linalg, only: inv,invert,operator(.inv.),eye
8
+ use stdlib_linalg_state, only: linalg_state_type,LINALG_ERROR
9
9
10
10
implicit none (type,external)
11
11
private
@@ -21,29 +21,27 @@ module test_linalg_inverse
21
21
22
22
allocate(tests(0))
23
23
24
- #:for rk,rt,ri in REAL_KINDS_TYPES
24
+ #:for rk,rt,ri in RC_KINDS_TYPES
25
25
#:if rk!="xdp"
26
- tests = [tests,new_unittest("inverse_${ri}$_eye_inverse",test_${ri}$_eye_inverse)]
26
+ tests = [tests,new_unittest("${ri}$_eye_inverse",test_${ri}$_eye_inverse), &
27
+ new_unittest("${ri}$_singular_inverse",test_${ri}$_singular_inverse)]
27
28
#:endif
28
29
#:endfor
29
30
30
31
end subroutine test_inverse_matrix
31
32
32
- !> Invert real identity matrix
33
33
#:for rk,rt,ri in REAL_KINDS_TYPES
34
34
#:if rk!="xdp"
35
+ !> Invert real identity matrix
35
36
subroutine test_${ri}$_eye_inverse(error)
36
37
type(error_type), allocatable, intent(out) :: error
37
38
38
39
type(linalg_state_type) :: state
39
40
40
41
integer(ilp), parameter :: n = 25_ilp
41
- integer(ilp) :: i,j
42
42
${rt}$ :: a(n,n),inva(n,n)
43
43
44
- do concurrent (i=1:n,j=1:n)
45
- a(i,j) = merge(1.0_${rk}$,0.0_${rk}$,i==j)
46
- end do
44
+ a = eye(n)
47
45
48
46
!> Inverse function
49
47
inva = inv(a,err=state)
@@ -64,6 +62,27 @@ module test_linalg_inverse
64
62
65
63
end subroutine test_${ri}$_eye_inverse
66
64
65
+ !> Invert singular matrix
66
+ subroutine test_${ri}$_singular_inverse(error)
67
+ type(error_type), allocatable, intent(out) :: error
68
+
69
+ type(linalg_state_type) :: err
70
+
71
+ integer(ilp), parameter :: n = 25_ilp
72
+ ${rt}$ :: a(n,n)
73
+
74
+ a = eye(n)
75
+
76
+ !> Make rank-deficient
77
+ a(12,12) = 0
78
+
79
+ !> Inverse function
80
+ call invert(a,err=err)
81
+ call check(error,err%state==LINALG_ERROR,'singular ${rt}$ inverse returned '//err%print())
82
+ if (allocated(error)) return
83
+
84
+ end subroutine test_${ri}$_singular_inverse
85
+
67
86
#:endif
68
87
#:endfor
69
88
@@ -132,6 +151,24 @@ module test_linalg_inverse
132
151
133
152
end subroutine test_${ci}$_eye_inverse
134
153
154
+ !> Invert singular matrix
155
+ subroutine test_${ci}$_singular_inverse(error)
156
+ type(error_type), allocatable, intent(out) :: error
157
+
158
+ type(linalg_state_type) :: err
159
+
160
+ integer(ilp), parameter :: n = 25_ilp
161
+ ${ct}$ :: a(n,n)
162
+
163
+ a = (0.0_${ck}$,0.0_${ck}$)
164
+
165
+ !> Inverse function
166
+ call invert(a,err=err)
167
+ call check(error,err%state==LINALG_ERROR,'singular ${ct}$ inverse returned '//err%print())
168
+ if (allocated(error)) return
169
+
170
+ end subroutine test_${ci}$_singular_inverse
171
+
135
172
#:endif
136
173
#:endfor
137
174
0 commit comments