      SUBROUTINE ComplexGaussJordan (coefary, constary, numcol,
     +         solcoef, invary, det)
      INCLUDE 'STDHDR.FOR'
      COMPLEX coefary(0:maxc, 0:maxc), constary(0:maxc)
      COMPLEX solcoef(0:maxc), invary(0:maxc, 0:maxc)
      REAL det
      REAL inv[ALLOCATABLE](:,:),coef[ALLOCATABLE](:,:)
      REAL sol[ALLOCATABLE](:),cons[ALLOCATABLE](:)
      INTEGER i, j, nn, ii, jj,iErr

      ALLOCATE(inv(0:maxc,0:maxc),coef(0:maxc, 0:maxc), STAT=iErr)
      CALL CheckMem(iErr)
      ALLOCATE(sol(0:maxc),cons(0:maxc),STAT=iErr)
      CALL CheckMem(iErr)

      DO i = 0, numcol - 1
        ii = 2 * i
        DO j = 0, numcol - 1
          jj = 2 * j
          coef(ii, jj) = REAL(coefary(i, j))
          coef(ii, jj  + 1) = -IMAG(coefary(i, j))
          coef(ii + 1, jj) = IMAG(coefary(i, j))
          coef(ii + 1, jj + 1) = REAL(coefary(i, j))
        END DO
        cons(ii) = REAL(constary(i))
        cons(ii + 1) = IMAG(constary(i))
      END DO
      nn = 2 * numcol
      CALL GaussJordan(coef, cons, nn, sol, inv, det)

      DO i = 0, numcol - 1
        ii = 2 * i
        DO j = 0, numcol - 1
          jj = 2 * j
          invary(i, j) = CMPLX(inv(ii, jj),inv(ii + 1, jj))
        END DO
      solcoef(i) = CMPLX(sol(ii),sol(ii + 1))
      END DO
      DEALLOCATE(inv,coef, STAT=iErr)
      CALL CheckDealloc(iErr)
      DEALLOCATE(sol,cons,STAT=iErr)
      CALL CheckDealloc(iErr)
      END !SUBROUTINE ComplexGaussJordan



