1      subroutine dimqm_printAtomicDipoles(rtdb, lpol, liscmplx, label)
2      implicit none
3#include "stdio.fh"
4#include "util.fh"
5#include "rtdb.fh"
6#include "errquit.fh"
7#include "mafdecls.fh"
8#include "dimqm_constants.fh"
9#include "dimqm.fh"
10#include "global.fh"
11c
12c   Input variables
13      integer rtdb
14      logical lpol
15      logical liscmplx
16      character(len=*) label
17c
18c   Local variables
19      double precision diptot(3)
20      character elements(2,nDIMTypes)
21      integer tx(nDIM)
22      integer l_muind, k_muind
23      character*256 f
24      integer i, i3, j, itype
25      integer icmplx
26      logical stat
27      character*50 dd
28      character*50 d
29      integer l_name, k_name
30      character*2 ele
31      character aele(2)
32      double precision mu(3)
33      double precision mu_c(3)
34      character*16 mu_str
35c      character*16 q_str
36c
37c   Return without execution if this table is not requested
38      if(.not.latmdip) return
39      if(ga_nodeid().eq..0) then
40        stat = rtdb_parallel(.false.)
41c
42c   Check if real or complex
43      icmplx = 1
44      if(liscmplx) icmplx = 2
45      i3 = nDIM * 3
46c
47c   Determine which dipoles to grab
48      mu_str = 'dimqm:muind'
49c
50c   ===============
51c   Allocate Arrays
52c   ===============
53c
54c   Allocate memory for induced dipoles
55      if(.not.ma_push_get(mt_dbl,i3*icmplx,'dimqm muind',
56     $                                  l_muind,k_muind))
57     $  call errquit('printAtomicDipoles malloc muind failed',1,MA_ERR)
58c   Allocate memory for atom names
59c   Read in atom names
60      if(.not.rtdb_get(rtdb,'dimpar:typeindex',mt_int,nDIM,tx))
61     $  call errquit('printAtomicDipoles get tx failed',1,RTDB_ERR)
62      if(.not.rtdb_get(rtdb, 'dimpar:name', mt_byte, nDIMTypes*2,
63     $                 elements))
64     $  call errquit('printAtomicDipoles get names failed',1,RTDB_ERR)
65
66      dd =
67     $ '=================================================='
68      d =
69     $ '--------------------------------------------------'
70c
71c -------------------------------------------------
72c Print out the dipoles for each atom
73c -------------------------------------------------
74c
75      if(liscmplx) then
76
77c   Read in induced dipoles
78      if(.not.rtdb_get(rtdb,mu_str,mt_dbl,i3*icmplx,dbl_mb(k_muind)))
79     $  call errquit('dimqm_print get muind(c) failed',1, RTDB_ERR)
80
81        f = '(1X,A)'
82        write (luout,f) dd
83        write (luout,f) 'Induced dipoles for each DIM atom :'
84        if (label .ne. '') write (luout,f) label
85        write (luout,f) d
86        f = '(5X,A,17X,A)'
87        write (luout,f) 'ATOM', 'Dipole'
88        f = '(17X,A1,11X,A1,11X,A1)'
89        write (luout,f) 'X', 'Y', 'Z'
90        f = '(1X,A)'
91        write (luout,f) d
92        f = '(1X,I5,1X,2A1,3F12.5,1X,A3)'
93        do i = 1, nDIM
94          itype = tx(i)
95          do j = 1, 3
96            mu(j)   = dbl_mb(k_muind+3*(i-1)+j-1)
97            mu_c(j) = dbl_mb(k_muind+3*(i-1)+j-1+i3)
98          end do
99          write (luout,f) i, elements(1:2,itype),  mu(:), '(R)'
100          write (luout,'(9X,3F12.5,1X,A3)')  mu_c(:), '(I)'
101        end do
102        f = '(1X,A)'
103        write (luout,f) dd
104        write (luout,*)
105
106      else
107
108c   Read in induced dipoles
109      if(.not.rtdb_get(rtdb,mu_str,mt_dbl,i3,dbl_mb(k_muind)))
110     $  call errquit('dimqm_print get muind(r) failed',1, RTDB_ERR)
111
112        f = '(1X,A)'
113        write (luout,f) dd
114        write (luout,f) 'Induced dipoles for each DIM atom :'
115        if (label .ne. '') write (luout,f) label
116        write (luout,f) d
117        f = '(5X,A,17X,A)'
118        write (luout,f) 'ATOM', 'Dipole'
119        f = '(17X,A1,11X,A1,11X,A1)'
120        write (luout,f) 'X', 'Y', 'Z'
121        f = '(1X,A)'
122        write (luout,f) d
123        f = '(1X,I5,1X,2A1,3F12.5)'
124        do i =1, nDIM
125          itype = tx(i)
126          do j = 1, 3
127            mu(j) = dbl_mb(k_muind+3*(i-1)+j-1)
128          end do
129          write (luout,f) i, elements(1:2,itype), mu(:)
130        end do
131        f = '(1X,A)'
132        write (luout,f) dd
133        write (luout,*)
134
135      end if
136      call util_flush(LuOut)
137      stat = ma_chop_stack(l_muind)
138      stat = rtdb_parallel(.true.)
139      end if
140
141      end subroutine dimqm_printAtomicDipoles
142