@@ -10,7 +10,8 @@ module m_eigen_solver
10
10
11
11
implicit none
12
12
13
- private ; public :: cg, cbal, corth, comqr2, csroot, cdiv, pythag
13
+ private ;
14
+ public :: cg, cbal, corth, comqr2, csroot, cdiv, pythag
14
15
15
16
contains
16
17
@@ -51,9 +52,14 @@ subroutine cg(nm, nl, ar, ai, wr, wi, zr, zi, fv1, fv2, fv3, ierr)
51
52
! this version dated august 1983.
52
53
!
53
54
! ------------------------------------------------------------------
54
- integer nm, nl, is1, is2, ierr
55
- real (kind (0d0 )), dimension (nm, nl) :: ar, ai, zr, zi
56
- real (kind (0d0 )), dimension (nl) :: wr, wi, fv1, fv2, fv3
55
+ integer , intent (in ) :: nm, nl
56
+ real (kind (0d0 )), dimension (nm:nl), intent (inout ) :: ar, ai
57
+ real (kind (0d0 )), dimension (nl), intent (out ) :: wr, wi
58
+ real (kind (0d0 )), dimension (nm, nl), intent (out ) :: zr, zi
59
+ real (kind (0d0 )), dimension (nl), intent (out ) :: fv1, fv2, fv3
60
+ integer , intent (out ) :: ierr
61
+
62
+ integer :: is1, is2
57
63
58
64
if (nl <= nm) go to 10
59
65
ierr = 10 * nl
@@ -125,11 +131,14 @@ subroutine cbal(nm, nl, ar, ai, low, igh, scale)
125
131
! this version dated august 1983.
126
132
!
127
133
! ------------------------------------------------------------------
128
- integer i, j, k, l, ml, nl, jj, nm, igh, low, iexc
129
- real (kind (0d0 )), dimension (nm, nl) :: ar, ai
130
- real (kind (0d0 )), dimension (nl) :: scale
134
+ integer , intent (in ) :: nm, nl
135
+ real (kind (0d0 )), dimension (nm, nl), intent (inout ) :: ar, ai
136
+ integer , intent (out ) :: low, igh
137
+ real (kind (0d0 )), dimension (nl), intent (out ) :: scale
138
+
139
+ integer :: i, j, k, l, ml, jj, iexc
131
140
real (kind (0d0 )) :: c, f, g, r, s, b2, radix
132
- logical noconv
141
+ logical :: noconv
133
142
134
143
radix = 16.0d0
135
144
@@ -295,12 +304,13 @@ subroutine corth(nm, nl, low, igh, ar, ai, ortr, orti)
295
304
! this version dated august 1983.
296
305
!
297
306
! ------------------------------------------------------------------
298
- integer i, j, ml, nl, ii, jj, la, mp, nm, igh, kp1, low
299
- real (kind (0d0 )), dimension (nm, nl) :: ar, ai
300
- real (kind (0d0 )), dimension (igh) :: ortr, orti
307
+ integer , intent (in ) :: nm, nl, low, igh
308
+ real (kind (0d0 )), dimension (nm, nl), intent (inout ) :: ar, ai
309
+ real (kind (0d0 )), dimension (igh), intent (out ) :: ortr, orti
310
+
311
+ integer :: i, j, ml, ii, jj, la, mp, kp1, mll
301
312
real (kind (0d0 )) :: f, g, h, fi, fr, scale, c
302
313
303
- integer mll
304
314
mll = 6
305
315
306
316
la = igh - 1
@@ -460,11 +470,14 @@ subroutine comqr2(nm, nl, low, igh, ortr, orti, hr, hi, wr, wi, zr, zi, ierr)
460
470
! this version dated october 1989.
461
471
!
462
472
! ------------------------------------------------------------------
463
- integer i, j, k, l, ml, nl, en, ii, jj, ll, nm, nn, igh, ip1, &
464
- itn, its, low, lp1, enm1, iend, ierr
465
- real (kind (0d0 )), dimension (nm, nl) :: hr, hi, zr, zi
466
- real (kind (0d0 )), dimension (nl) :: wr, wi
467
- real (kind (0d0 )), dimension (igh) :: ortr, orti
473
+ integer , intent (in ) :: nm, nl, low, igh
474
+ real (kind (0d0 )), dimension (nm, nl), intent (inout ) :: hr, hi
475
+ real (kind (0d0 )), dimension (nl), intent (out ) :: wr, wi
476
+ real (kind (0d0 )), dimension (nm, nl), intent (out ) :: zr, zi
477
+ real (kind (0d0 )), dimension (igh), intent (inout ) :: ortr, orti
478
+ integer , intent (out ) :: ierr
479
+
480
+ integer :: i, j, k, l, ml, en, ii, jj, ll, nn, ip1, itn, its, lp1, enm1, iend
468
481
real (kind (0d0 )) :: si, sr, ti, tr, xi, xr, yi, yr, zzi, zzr, &
469
482
norm, tst1, tst2, c, d
470
483
!
@@ -843,9 +856,13 @@ subroutine cbabk2(nm, nl, low, igh, scale, ml, zr, zi)
843
856
!
844
857
! ------------------------------------------------------------------
845
858
!
846
- integer i, j, k, ml, nl, ii, nm, igh, low
847
- double precision scale (nl), zr(nm, ml), zi(nm, ml)
848
- double precision s
859
+ integer , intent (in ) :: nm, nl, low, igh
860
+ double precision , intent (in ) :: scale (nl)
861
+ integer , intent (in ) :: ml
862
+ double precision , intent (inout ) :: zr(nm, ml), zi(nm, ml)
863
+
864
+ integer :: i, j, k, ii
865
+ double precision :: s
849
866
850
867
if (ml == 0 ) go to 200
851
868
if (igh == low) go to 120
@@ -885,7 +902,8 @@ subroutine cbabk2(nm, nl, low, igh, scale, ml, zr, zi)
885
902
end subroutine cbabk2
886
903
887
904
subroutine csroot (xr , xi , yr , yi )
888
- real (kind (0d0 )) :: xr, xi, yr, yi
905
+ real (kind (0d0 )), intent (in ) :: xr, xi
906
+ real (kind (0d0 )), intent (out ) :: yr, yi
889
907
!
890
908
! (yr,yi) = complex dsqrt(xr,xi)
891
909
! branch chosen so that yr .ge. 0.0 and sign(yi) .eq. sign(xi)
@@ -904,7 +922,8 @@ subroutine csroot(xr, xi, yr, yi)
904
922
end subroutine csroot
905
923
906
924
subroutine cdiv (ar , ai , br , bi , cr , ci )
907
- real (kind (0d0 )) :: ar, ai, br, bi, cr, ci
925
+ real (kind (0d0 )), intent (in ) :: ar, ai, br, bi
926
+ real (kind (0d0 )), intent (out ) :: cr, ci
908
927
!
909
928
! complex division, (cr,ci) = (ar,ai)/(br,bi)
910
929
!
@@ -921,7 +940,8 @@ subroutine cdiv(ar, ai, br, bi, cr, ci)
921
940
end subroutine cdiv
922
941
923
942
subroutine pythag (a , b , c )
924
- real (kind (0d0 )) :: a, b, c
943
+ real (kind (0d0 )), intent (in ) :: a, b
944
+ real (kind (0d0 )), intent (out ) :: c
925
945
!
926
946
! finds dsqrt(a**2+b**2) without overflow or destructive underflow
927
947
!
0 commit comments