!************************************************************
 !*  This program calculates the determinant of a complex    *
 !* square matrix using Function CFindDet (see ref. below).  *
 !* -------------------------------------------------------- *
 !* SAMPLE RUN:                                              *
 !*                                                          *
 !* Calculate the deteminant of square matrix:               *
 !*                                                          *
 !*           ! 0  -1  -1  -1 !                              *
 !*           ! 1   0  -1  -1 !                              *
 !*           ! 1   1   0  -1 !                              *
 !*           ! 1   1   1   0 !                              *
 !*                                                          *
 !*  Determinant =  (1.000000,0.000000)                      *
 !*                                                          *
 !************************************************************ 
 Program TCFindDet

 integer, parameter :: N=4

 Complex  A(N,N)
 Complex  D
 
 data  A  / (0,0),(1,0),(1,0),(1,0), &  
           (-1,0),(0,0),(1,0),(1,0), &
	   (-1,0),(-1,0),(0,0),(1,0), &
	   (-1,0),(-1,0),(-1,0),(0,0) /

 D = CFindDet(A, N)

 print *,' '
 print *,' Determinant =', D
 print *,' '

 stop

 End
 
 !----------------------------------------------------------------------------------------------------
 !Function to find the determinant of a square matrix
 !Author : Louisda16th a.k.a Ashwith J. Rego
 !Description: The subroutine is based on two key points:
 !1] A determinant is unaltered when row operations are performed: Hence, using this principle,
 !row operations (column operations would work as well) are used
 !to convert the matrix into upper traingular form
 !2]The determinant of a triangular matrix is obtained by finding the product of the diagonal elements
 !----------------------------------------------------------------------------------------------------
 REAL FUNCTION CFindDet(matrix, n)
     IMPLICIT NONE
     INTEGER, INTENT(IN) :: n
     COMPLEX, DIMENSION(n,n) :: matrix
     COMPLEX :: m, temp, CONE, CZERO
     INTEGER :: i, j, k, l
     LOGICAL :: DetExists = .TRUE.
     l = 1
	 CZERO=CMPLX(0,0); CONE=CMPLX(1,0)
     !Convert to upper triangular form
     DO k = 1, n-1
         IF (matrix(k,k) == CZERO) THEN
             DetExists = .FALSE.
             DO i = k+1, n
                 IF (matrix(i,k) /= CZERO) THEN
                     DO j = 1, n
                         temp = matrix(i,j)
                         matrix(i,j)= matrix(k,j)
                         matrix(k,j) = temp
                     END DO
                     DetExists = .TRUE.
                     l=-l
                     EXIT
                 ENDIF
             END DO
             IF (DetExists .EQV. .FALSE.) THEN
                 CFindDet = CZERO
                 return
             END IF
         ENDIF
         DO j = k+1, n
             m = matrix(j,k)/matrix(k,k)
             DO i = k+1, n
                 matrix(j,i) = matrix(j,i) - m*matrix(k,i)
             END DO
         END DO
     END DO
     
     !Calculate determinant by finding product of diagonal elements
     CFindDet = CONE
     DO i = 1, n
         CFindDet = CFindDet * matrix(i,i)
     END DO
     
 END FUNCTION CFindDet

!end of file cfinddet.f90