subroutine cgefa(a,lda,n,ipvt,info)
integer lda,n,ipvt(1),info
complex a(lda,1)
c
c cgefa factors a complex matrix by gaussian elimination.
c
c cgefa is usually called by cgeco, but it can be called
c directly with a saving in time if rcond is not needed.
c (time for cgeco) = (1 + 9/n)*(time for cgefa) .
c
c on entry
c
c a complex(lda, n)
c the matrix to be factored.
c
c lda integer
c the leading dimension of the array a .
c
c n integer
c the order of the matrix a .
c
c on return
c
c a an upper triangular matrix and the multipliers
c 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 cgesl or cgedi will divide by zero
c if called. use rcond in cgeco for a reliable
c indication of singularity.
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 caxpy,cscal,icamax
c fortran abs,aimag,real
c
c internal variables
c
complex t
integer icamax,j,k,kp1,l,nm1
c
complex zdum
real cabs1
cabs1(zdum) = abs(real(zdum)) + abs(aimag(zdum))
c
c gaussian elimination with partial pivoting
c
info = 0
nm1 = n - 1
if (nm1 .lt. 1) go to 70
do 60 k = 1, nm1
kp1 = k + 1
c
c find l = pivot index
c
l = icamax(n-k+1,a(k,k),1) + k - 1
ipvt(k) = l
c
c zero pivot implies this column already triangularized
c
if (cabs1(a(l,k)) .eq. 0.0e0) go to 40
c
c interchange if necessary
c
if (l .eq. k) go to 10
t = a(l,k)
a(l,k) = a(k,k)
a(k,k) = t
10 continue
c
c compute multipliers
c
t = -(1.0e0,0.0e0)/a(k,k)
call cscal(n-k,t,a(k+1,k),1)
c
c row elimination with column indexing
c
do 30 j = kp1, n
t = a(l,j)
if (l .eq. k) go to 20
a(l,j) = a(k,j)
a(k,j) = t
20 continue
call caxpy(n-k,t,a(k+1,k),1,a(k+1,j),1)
30 continue
go to 50
40 continue
info = k
50 continue
60 continue
70 continue
ipvt(n) = n
if (cabs1(a(n,n)) .eq. 0.0e0) info = n
return
end