Skip to content

HTTPS clone URL

Subversion checkout URL

You can clone with
or
.
Download ZIP
Browse files

BUG: sparse/isolve: fix bicgstab, cgs, and qmr for complex-valued types

These routines contained calls <_t>AXPY(N, ONE, ...), which will fail
for complex-valued data types, because the argument 'ONE' is specified
as double-precision real number and it is not automatically cast to
complex by the Fortran compiler.
  • Loading branch information...
commit e00f70c57bdcd48bf917927cc2f3a3aef9011ca3 1 parent 575262f
@pv pv authored
View
5 scipy/sparse/linalg/isolve/iterative/BiCGSTABREVCOM.f.src
@@ -116,7 +116,7 @@
<rt> TOL, BNRM2,
$ RHOTOL, OMEGATOL, <sdsd=s,d,s,d>GETBREAK,
$ <rc=s,d,sc,dz>NRM2
- <_t> ALPHA, BETA, RHO, RHO1, OMEGA,
+ <_t> ALPHA, BETA, RHO, RHO1, OMEGA, TMPVAL,
$ <xdot=sdot,ddot,cdotc,zdotc>
* indicates where to resume from. Only valid when IJOB = 2!
INTEGER RLBL
@@ -274,7 +274,8 @@
BETA = ( RHO / RHO1 ) * ( ALPHA / OMEGA )
CALL <_c>AXPY( N, -OMEGA, WORK(1,V), 1, WORK(1,P), 1 )
CALL <_c>SCAL( N, BETA, WORK(1,P), 1 )
- CALL <_c>AXPY( N, ONE, WORK(1,R), 1, WORK(1,P), 1 )
+ TMPVAL = ONE
+ CALL <_c>AXPY( N, TMPVAL, WORK(1,R), 1, WORK(1,P), 1 )
ELSE
CALL <_c>COPY( N, WORK(1,R), 1, WORK(1,P), 1 )
ENDIF
View
8 scipy/sparse/linalg/isolve/iterative/CGSREVCOM.f.src
@@ -117,7 +117,7 @@
$ <sdsd=s,d,s,d>GETBREAK,
$ <rc=s,d,sc,dz>NRM2
- <_t> ALPHA, BETA, RHO, RHO1,
+ <_t> ALPHA, BETA, RHO, RHO1, TMPVAL,
$ <xdot=sdot,ddot,cdotc,zdotc>
* ..
* indicates where to resume from. Only valid when IJOB = 2!
@@ -291,7 +291,8 @@
*
CALL <_c>SCAL( N, BETA**2, WORK(1,P), 1 )
CALL <_c>AXPY( N, BETA, WORK(1,Q), 1, WORK(1,P), 1 )
- CALL <_c>AXPY( N, ONE, WORK(1,U), 1, WORK(1,P), 1 )
+ TMPVAL = ONE
+ CALL <_c>AXPY( N, TMPVAL, WORK(1,U), 1, WORK(1,P), 1 )
ELSE
CALL <_c>COPY( N, WORK(1,R), 1, WORK(1,U), 1 )
CALL <_c>COPY( N, WORK(1,U), 1, WORK(1,P), 1 )
@@ -336,7 +337,8 @@
* PHAT is being used as temporary storage here.
*
CALL <_c>COPY( N, WORK(1,Q), 1, WORK(1,PHAT), 1 )
- CALL <_c>AXPY( N, ONE, WORK(1,U), 1, WORK(1,PHAT), 1 )
+ TMPVAL = ONE
+ CALL <_c>AXPY( N, TMPVAL, WORK(1,U), 1, WORK(1,PHAT), 1 )
*********CALL PSOLVE( WORK(1,UHAT), WORK(1,PHAT) )
*
NDX1 = ((UHAT - 1) * LDW) + 1
View
4 scipy/sparse/linalg/isolve/iterative/GMRESREVCOM.f.src
@@ -415,6 +415,7 @@
$ <rc=s,d,sc,dz>NRM2, ONE
PARAMETER ( ONE = 1.0D+0 )
<_t> <xdot=sdot,ddot,cdotc,zdotc>
+ <_t> TMPVAL
EXTERNAL <_c>AXPY, <_c>COPY, <xdot>, <rc>NRM2, <_c>SCAL
*
DO 10 K = 1, I
@@ -423,7 +424,8 @@
10 CONTINUE
H( I+1 ) = <rc>NRM2( N, W, 1 )
CALL <_c>COPY( N, W, 1, V( 1,I+1 ), 1 )
- CALL <_c>SCAL( N, ONE / H( I+1 ), V( 1,I+1 ), 1 )
+ TMPVAL = ONE / H( I+1 )
+ CALL <_c>SCAL( N, TMPVAL, V( 1,I+1 ), 1 )
*
RETURN
*
View
15 scipy/sparse/linalg/isolve/iterative/QMRREVCOM.f.src
@@ -127,7 +127,7 @@
<_t> BETA, GAMMA, GAMMA1, DELTA, EPS, ETA, XI,
- $ RHO, RHO1, THETA, THETA1, C1,
+ $ RHO, RHO1, THETA, THETA1, C1, TMPVAL,
$ <xdot=sdot,ddot,cdotc,zdotc>,
$ toz
*
@@ -347,12 +347,14 @@
$ GO TO 25
*
CALL <_c>COPY( N, WORK(1,VTLD), 1, WORK(1,V), 1 )
- CALL <_c>SCAL( N, ONE / RHO, WORK(1,V), 1 )
- CALL <_c>SCAL( N, ONE / RHO, WORK(1,Y), 1 )
+ TMPVAL = ONE / RHO
+ CALL <_c>SCAL( N, TMPVAL, WORK(1,V), 1 )
+ CALL <_c>SCAL( N, TMPVAL, WORK(1,Y), 1 )
*
+ TMPVAL = ONE / XI
CALL <_c>COPY( N, WORK(1,WTLD), 1, WORK(1,W), 1 )
- CALL <_c>SCAL( N, ONE / XI, WORK(1,W), 1 )
- CALL <_c>SCAL( N, ONE / XI, WORK(1,Z), 1 )
+ CALL <_c>SCAL( N, TMPVAL, WORK(1,W), 1 )
+ CALL <_c>SCAL( N, TMPVAL, WORK(1,Z), 1 )
*
DELTA = <xdot>( N, WORK(1,Z), 1, WORK(1,Y), 1 )
IF ( ABS( DELTA ).LT.DELTATOL ) GO TO 25
@@ -482,7 +484,8 @@
*
* Compute current solution vector x.
*
- CALL <_c>AXPY( N, ONE, WORK(1,D), 1, X, 1 )
+ TMPVAL = ONE
+ CALL <_c>AXPY( N, TMPVAL, WORK(1,D), 1, X, 1 )
*
* Compute residual vector rk, find norm,
* then check for tolerance.
Please sign in to comment.
Something went wrong with that request. Please try again.