1c $Id$
2      subroutine et_det(op,handle,n,detval)
3#include "rtdb.fh"
4#include "mafdecls.fh"
5#include "inp.fh"
6#include "util.fh"
7#include "global.fh"
8#include "apiP.fh"
9#include "bas.fh"
10#include "geom.fh"
11#include "cscf.fh"
12c
13c returns the determinant of the square matrix with file handle 'handle'.
14c n is the dimension of the array, detval is the returned value of the
15c determinant. If op='N', a=a. If op='T', a=a'.
16c
17      character*1 op
18      integer handle
19      integer n,i,j,info
20      integer l_ipvt,k_ipvt,l_a,k_a,l_work,k_work
21      double precision detval,det(2),ten
22c
23      ten = 10.0d0
24      if (.not.ma_push_get(mt_int,n,'ipvt',
25     1  l_ipvt,k_ipvt)) call errquit
26     2  ('et_det: failed to allocate ipvt',0,0)
27      if (.not.ma_push_get(mt_dbl,n*n,'a',
28     1  l_a,k_a)) call errquit
29     2  ('et_det: failed to allocate a',0,0)
30      if (.not.ma_push_get(mt_dbl,n,'work',
31     1  l_work,k_work)) call errquit
32     2  ('et_det: failed to allocate work',0,0)
33c
34      do i = 1,n
35      do j = 1,n
36       if(op.eq.'N') call ga_get(handle,i,i,j,j,
37     1 dbl_mb(k_a+i-1+(j-1)*n),1)
38       if(op.eq.'T') call ga_get(handle,j,j,i,i,
39     1 dbl_mb(k_a+i-1+(j-1)*n),1)
40      enddo
41      enddo
42c
43      call dgefa(dbl_mb(k_a),n,n,int_mb(k_ipvt),info)
44      call dgedi(dbl_mb(k_a),n,n,int_mb(k_ipvt),det,dbl_mb(k_work),10)
45      detval =  det(1) * ten**det(2)
46c
47      if (.not.ma_pop_stack(l_work)) call errquit
48     1  ('et_det: failed to deallocate work',0,0)
49      if (.not.ma_pop_stack(l_a)) call errquit
50     1  ('et_det: failed to deallocate a',0,0)
51      if (.not.ma_pop_stack(l_ipvt)) call errquit
52     1  ('et_det: failed to deallocate ipvt',0,0)
53      call ga_sync()
54      return
55      end
56
57