phynexis-deps/lapack
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

fix a few bugs

Browse Source
This commit is contained in:
thijs
2021-02-14 11:25:07 +01:00
parent 4a26b58426
commit 2d276d9e02
3 changed files with 40 additions and 24 deletions
Download Patch File Download Diff File Expand all files Collapse all files
SRC/dlaqz0.f
+17 -15
View File
@@ -409,8 +409,8 @@
* rcost = ilaenv( 17,'DLAQZ0',jbcmpz,n,ilo,ihi,lwork )
rcost = 10
itemp1 = int(nshifts/sqrt( 1+2*nshifts/(dble(rcost)/100*n) ))
itemp1 = ((k-1)/4)*4+4
itemp1 = int(nsr/sqrt( 1+2*nsr/(dble(rcost)/100*n) ))
itemp1 = ((itemp1-1)/4)*4+4
nbr = nsr+itemp1
if( n .lt. nmin ) then
@@ -442,7 +442,7 @@
info =-19
end if
if( info.NE.0 ) then
call xerbla( 'DLAQZ0',-info )
call xerbla( 'DLAQZ0',info )
return
end if
*
@@ -462,7 +462,7 @@
istop = ihi
maxit = 3*(ihi-ilo+1)
ld = 0
do iiter = 1,maxit
if(iiter .ge. maxit) then
info = istop
@@ -473,28 +473,30 @@
end if
* Check deflations at the end
if (abs(A(istop-1,istop-2)) .le. ulp*(abs(A(istop-1,
$ istop-1))+abs(A(istop-2,istop-2)))) then
if (abs(A(istop-1,istop-2)) .le. max(smlnum,ulp*(abs(A(istop-1,
$ istop-1))+abs(A(istop-2,istop-2))))) then
A(istop-1,istop-2) = zero
istop = istop-2
ld = 0
eshift = zero
else if (abs(A(istop,istop-1)) .le. ulp*(abs(A(istop,
$ istop))+abs(A(istop-1,istop-1)))) then
else if (abs(A(istop,istop-1)) .le. max(smlnum,
$ ulp*(abs(A(istop,istop))+abs(A(istop-1,istop-1))))) then
A(istop,istop-1) = zero
istop = istop-1
ld = 0
eshift = zero
end if
* Check deflations at the start
if (abs(A(istart+2,istart+1)) .le. ulp*(abs(A(istart+1,
$ istart+1))+abs(A(istart+2,istart+2)))) then
if (abs(A(istart+2,istart+1)) .le. max(smlnum,
$ ulp*(abs(A(istart+1,istart+1))+abs(A(istart+2,
$ istart+2))))) then
A(istart+2,istart+1) = zero
istart = istart+2
ld = 0
eshift = zero
else if (abs(A(istart+1,istart)) .le. ulp*(abs(A(istart,
$ istart))+abs(A(istart+1,istart+1)))) then
else if (abs(A(istart+1,istart)) .le. max(smlnum,
$ ulp*(abs(A(istart,istart))+abs(A(istart+1,
$ istart+1))))) then
A(istart+1,istart) = zero
istart = istart+1
ld = 0
@@ -508,8 +510,8 @@
* Check interior deflations
istart2 = istart
do k = istop,istart+1,-1
if (abs(A(k,k-1)) .le. ulp*(abs(A(k,k))+abs(A(k-1,
$ k-1)))) then
if (abs(A(k,k-1)) .le. max(smlnum,ulp*(abs(A(k,k))+abs(A(k-
$ 1,k-1))))) then
A(k,k-1) = zero
istart2 = k
exit
@@ -671,7 +673,7 @@
end do
80 call dlaqz3(ilschur,ilq,ilz,n,ilo,ihi,A,ldA,B,ldB,Q,ldQ,Z,ldZ,
80 call dlaqz3(ilschur,ilq,ilz,n,1,n,A,ldA,B,ldB,Q,ldQ,Z,ldZ,
$alphar,alphai,beta,norm_info)
end subroutine
SRC/dlaqz1.f
+15 -1
View File
@@ -119,7 +119,13 @@
parameter(zero=0.0d0,one=1.0d0,half=0.5d0)
* Local scalars
double precision :: w(2)
double precision :: w(2),safmin,safmax
* External Functions
double precision,external :: dlamch
safmin = dlamch('SAFE MINIMUM')
safmax = one/safmin
* Calculate first shifted vector
w(1) = beta1*A(1,1)-sr1*B(1,1)
@@ -140,4 +146,12 @@
* Account for imaginary part
v(1) = v(1)+si*si*B(1,1)
if( abs(v(1)).gt.safmax .or. abs(v(2)) .gt. safmax .or. abs(v(3
$ )).gt.safmax .or. isnan(v(1)) .or. isnan(v(2)) .or. isnan(v(
$ 3)) ) then
v(1) = zero
v(2) = zero
v(3) = zero
end if
end subroutine
SRC/dlaqz6.f
+8 -8
View File
@@ -324,7 +324,7 @@
info =-17
end if
if( info.NE.0 ) then
call xerbla( 'DLAQZ0',-info )
call xerbla( 'DLAQZ6',-info )
return
end if
@@ -343,22 +343,22 @@
jbcmpz(2:2)=wantQ
jbcmpz(3:3)=wantZ
nmin = ilaenv( 12,'DLAQZ0',jbcmpz,n,ilo,ihi,lwork )
nmin = ilaenv( 12,'DLAQZ6',jbcmpz,n,ilo,ihi,lwork )
nwr = ilaenv( 13,'DLAQZ0',jbcmpz,n,ilo,ihi,lwork )
nwr = ilaenv( 13,'DLAQZ6',jbcmpz,n,ilo,ihi,lwork )
nwr = max( 2,nwr )
nwr = min( ihi-ilo+1,( n-1 ) / 3,nwr )
nibble = ilaenv( 14,'DLAQZ0',jbcmpz,n,ilo,ihi,lwork )
nibble = ilaenv( 14,'DLAQZ6',jbcmpz,n,ilo,ihi,lwork )
nsr = ilaenv( 15,'DLAQZ0',jbcmpz,n,ilo,ihi,lwork )
nsr = ilaenv( 15,'DLAQZ6',jbcmpz,n,ilo,ihi,lwork )
nsr = min( nsr,( n+6 ) / 9,ihi-ilo )
nsr = max( 2,nsr-mod( nsr,2 ) )
* rcost = ilaenv( 17,'DLAQZ0',jbcmpz,n,ilo,ihi,lwork )
rcost = 10
itemp1 = int(nshifts/sqrt( 1+2*nshifts/(dble(rcost)/100*n) ))
itemp1 = ((k-1)/4)*4+4
itemp1 = int(nsr/sqrt( 1+2*nsr/(dble(rcost)/100*n) ))
itemp1 = ((itemp1-1)/4)*4+4
nbr = nsr+itemp1
if( n .lt. nmin ) then
@@ -390,7 +390,7 @@
info =-19
end if
if( info.NE.0 ) then
call xerbla( 'DLAQZ0',-info )
call xerbla( 'DLAQZ6',info )
return
end if
*