diff --git a/ext/SConscript b/ext/SConscript index eaf59bd58..16c045727 100644 --- a/ext/SConscript +++ b/ext/SConscript @@ -6,6 +6,11 @@ localenv = env.Clone() def prep_default(env): return env.Clone() +def prep_fortran(env): + localenv = env.Clone() + localenv.Append(CPPPATH='#include/cantera/base') # for config.h + return localenv + def prep_f2c(env): localenv = env.Clone() @@ -73,7 +78,7 @@ if env['build_with_f2c']: libs.append(('f2c_lapack', ['c'], prep_f2c)) else: - libs.append(('math', ['cpp','c','f'], prep_default)) + libs.append(('math', ['cpp','c','f'], prep_fortran)) if env['BUILD_BLAS_LAPACK']: libs.append(('blas', ['f'], prep_default)) diff --git a/ext/blas/lsame.f b/ext/blas/lsame.f deleted file mode 100755 index f89517460..000000000 --- a/ext/blas/lsame.f +++ /dev/null @@ -1,87 +0,0 @@ - LOGICAL FUNCTION LSAME( CA, CB ) -* -* -- LAPACK auxiliary routine (version 2.0) -- -* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., -* Courant Institute, Argonne National Lab, and Rice University -* January 31, 1994 -* -* .. Scalar Arguments .. - CHARACTER CA, CB -* .. -* -* Purpose -* ======= -* -* LSAME returns .TRUE. if CA is the same letter as CB regardless of -* case. -* -* Arguments -* ========= -* -* CA (input) CHARACTER*1 -* CB (input) CHARACTER*1 -* CA and CB specify the single characters to be compared. -* -* ===================================================================== -* -* .. Intrinsic Functions .. - INTRINSIC ICHAR -* .. -* .. Local Scalars .. - INTEGER INTA, INTB, ZCODE -* .. -* .. Executable Statements .. -* -* Test if the characters are equal -* - LSAME = CA.EQ.CB - IF( LSAME ) - $ RETURN -* -* Now test for equivalence if both characters are alphabetic. -* - ZCODE = ICHAR( 'Z' ) -* -* Use 'Z' rather than 'A' so that ASCII can be detected on Prime -* machines, on which ICHAR returns a value with bit 8 set. -* ICHAR('A') on Prime machines returns 193 which is the same as -* ICHAR('A') on an EBCDIC machine. -* - INTA = ICHAR( CA ) - INTB = ICHAR( CB ) -* - IF( ZCODE.EQ.90 .OR. ZCODE.EQ.122 ) THEN -* -* ASCII is assumed - ZCODE is the ASCII code of either lower or -* upper case 'Z'. -* - IF( INTA.GE.97 .AND. INTA.LE.122 ) INTA = INTA - 32 - IF( INTB.GE.97 .AND. INTB.LE.122 ) INTB = INTB - 32 -* - ELSE IF( ZCODE.EQ.233 .OR. ZCODE.EQ.169 ) THEN -* -* EBCDIC is assumed - ZCODE is the EBCDIC code of either lower or -* upper case 'Z'. -* - IF( INTA.GE.129 .AND. INTA.LE.137 .OR. - $ INTA.GE.145 .AND. INTA.LE.153 .OR. - $ INTA.GE.162 .AND. INTA.LE.169 ) INTA = INTA + 64 - IF( INTB.GE.129 .AND. INTB.LE.137 .OR. - $ INTB.GE.145 .AND. INTB.LE.153 .OR. - $ INTB.GE.162 .AND. INTB.LE.169 ) INTB = INTB + 64 -* - ELSE IF( ZCODE.EQ.218 .OR. ZCODE.EQ.250 ) THEN -* -* ASCII is assumed, on Prime machines - ZCODE is the ASCII code -* plus 128 of either lower or upper case 'Z'. -* - IF( INTA.GE.225 .AND. INTA.LE.250 ) INTA = INTA - 32 - IF( INTB.GE.225 .AND. INTB.LE.250 ) INTB = INTB - 32 - END IF - LSAME = INTA.EQ.INTB -* -* RETURN -* -* End of LSAME -* - END diff --git a/ext/lapack/xerbla.f b/ext/lapack/xerbla.f deleted file mode 100755 index 618dfcf97..000000000 --- a/ext/lapack/xerbla.f +++ /dev/null @@ -1,46 +0,0 @@ - SUBROUTINE XERBLA( SRNAME, INFO ) -* -* -- LAPACK auxiliary routine (version 2.0) -- -* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., -* Courant Institute, Argonne National Lab, and Rice University -* September 30, 1994 -* -* .. Scalar Arguments .. - CHARACTER*6 SRNAME - INTEGER INFO -* .. -* -* Purpose -* ======= -* -* XERBLA is an error handler for the LAPACK routines. -* It is called by an LAPACK routine if an input parameter has an -* invalid value. A message is printed and execution stops. -* -* Installers may consider modifying the STOP statement in order to -* call system-specific exception-handling facilities. -* -* Arguments -* ========= -* -* SRNAME (input) CHARACTER*6 -* The name of the routine which called XERBLA. -* -* INFO (input) INTEGER -* The position of the invalid parameter in the parameter list -* of the calling routine. -* -* ===================================================================== -* -* .. Executable Statements .. -* - WRITE( *, FMT = 9999 )SRNAME, INFO -* - STOP -* - 9999 FORMAT( ' ** On entry to ', A6, ' parameter number ', I2, ' had ', - $ 'an illegal value' ) -* -* End of XERBLA -* - END diff --git a/ext/math/dgbefa.f b/ext/math/dgbefa.f deleted file mode 100644 index c26e6f579..000000000 --- a/ext/math/dgbefa.f +++ /dev/null @@ -1,174 +0,0 @@ - subroutine dgbfa(abd,lda,n,ml,mu,ipvt,info) - integer lda,n,ml,mu,ipvt(1),info - double precision abd(lda,1) -c -c dgbfa factors a double precision band matrix by elimination. -c -c dgbfa is usually called by dgbco, but it can be called -c directly with a saving in time if rcond is not needed. -c -c on entry -c -c abd double precision(lda, n) -c contains the matrix in band storage. the columns -c of the matrix are stored in the columns of abd and -c the diagonals of the matrix are stored in rows -c ml+1 through 2*ml+mu+1 of abd . -c see the comments below for details. -c -c lda integer -c the leading dimension of the array abd . -c lda must be .ge. 2*ml + mu + 1 . -c -c n integer -c the order of the original matrix. -c -c ml integer -c number of diagonals below the main diagonal. -c 0 .le. ml .lt. n . -c -c mu integer -c number of diagonals above the main diagonal. -c 0 .le. mu .lt. n . -c more efficient if ml .le. mu . -c on return -c -c abd an upper triangular matrix in band storage and -c the multipliers which were used to obtain it. -c the factorization can be written a = l*u where -c l is a product of permutation and unit lower -c triangular matrices and u is upper triangular. -c -c ipvt integer(n) -c an integer vector of pivot indices. -c -c info integer -c = 0 normal value. -c = k if u(k,k) .eq. 0.0 . this is not an error -c condition for this subroutine, but it does -c indicate that dgbsl will divide by zero if -c called. use rcond in dgbco for a reliable -c indication of singularity. -c -c band storage -c -c if a is a band matrix, the following program segment -c will set up the input. -c -c ml = (band width below the diagonal) -c mu = (band width above the diagonal) -c m = ml + mu + 1 -c do 20 j = 1, n -c i1 = max0(1, j-mu) -c i2 = min0(n, j+ml) -c do 10 i = i1, i2 -c k = i - j + m -c abd(k,j) = a(i,j) -c 10 continue -c 20 continue -c -c this uses rows ml+1 through 2*ml+mu+1 of abd . -c in addition, the first ml rows in abd are used for -c elements generated during the triangularization. -c the total number of rows needed in abd is 2*ml+mu+1 . -c the ml+mu by ml+mu upper left triangle and the -c ml by ml lower right triangle are not referenced. -c -c linpack. this version dated 08/14/78 . -c cleve moler, university of new mexico, argonne national lab. -c -c subroutines and functions -c -c blas daxpy,dscal,idamax -c fortran max0,min0 -c -c internal variables -c - double precision t - integer i,idamax,i0,j,ju,jz,j0,j1,k,kp1,l,lm,m,mm,nm1 -c -c - m = ml + mu + 1 - info = 0 -c -c zero initial fill-in columns -c - j0 = mu + 2 - j1 = min0(n,m) - 1 - if (j1 .lt. j0) go to 30 - do 20 jz = j0, j1 - i0 = m + 1 - jz - do 10 i = i0, ml - abd(i,jz) = 0.0d0 - 10 continue - 20 continue - 30 continue - jz = j1 - ju = 0 -c -c gaussian elimination with partial pivoting -c - nm1 = n - 1 - if (nm1 .lt. 1) go to 130 - do 120 k = 1, nm1 - kp1 = k + 1 -c -c zero next fill-in column -c - jz = jz + 1 - if (jz .gt. n) go to 50 - if (ml .lt. 1) go to 50 - do 40 i = 1, ml - abd(i,jz) = 0.0d0 - 40 continue - 50 continue -c -c find l = pivot index -c - lm = min0(ml,n-k) - l = idamax(lm+1,abd(m,k),1) + m - 1 - ipvt(k) = l + k - m -c -c zero pivot implies this column already triangularized -c - if (abd(l,k) .eq. 0.0d0) go to 100 -c -c interchange if necessary -c - if (l .eq. m) go to 60 - t = abd(l,k) - abd(l,k) = abd(m,k) - abd(m,k) = t - 60 continue -c -c compute multipliers -c - t = -1.0d0/abd(m,k) - call dscal(lm,t,abd(m+1,k),1) -c -c row elimination with column indexing -c - ju = min0(max0(ju,mu+ipvt(k)),n) - mm = m - if (ju .lt. kp1) go to 90 - do 80 j = kp1, ju - l = l - 1 - mm = mm - 1 - t = abd(l,j) - if (l .eq. mm) go to 70 - abd(l,j) = abd(mm,j) - abd(mm,j) = t - 70 continue - call daxpy(lm,t,abd(m+1,k),1,abd(mm+1,j),1) - 80 continue - 90 continue - go to 110 - 100 continue - info = k - 110 continue - 120 continue - 130 continue - ipvt(n) = n - if (abd(m,n) .eq. 0.0d0) info = n - return - end diff --git a/ext/math/idamax.f b/ext/math/idamax.f deleted file mode 100755 index eb0de47a7..000000000 --- a/ext/math/idamax.f +++ /dev/null @@ -1,39 +0,0 @@ - - integer function idamax(n,dx,incx) -c -c finds the index of element having max. absolute value. -c jack dongarra, linpack, 3/11/78. -c modified 3/93 to return if incx .le. 0. -c - double precision dx(1),dmax - integer i,incx,ix,n -c - idamax = 0 - if( n.lt.1 .or. incx.le.0 ) return - idamax = 1 - if(n.eq.1)return - if(incx.eq.1)go to 20 -c -c code for increment not equal to 1 -c - ix = 1 - dmax = dabs(dx(1)) - ix = ix + incx - do 10 i = 2,n - if(dabs(dx(ix)).le.dmax) go to 5 - idamax = i - dmax = dabs(dx(ix)) - 5 ix = ix + incx - 10 continue - return -c -c code for increment equal to 1 -c - 20 dmax = dabs(dx(1)) - do 30 i = 2,n - if(dabs(dx(i)).le.dmax) go to 30 - idamax = i - dmax = dabs(dx(i)) - 30 continue - return - end