1      function t2r(r, dist, invdist, a1, a2, pol1, pol2) result (tr)
2      implicit none
3#include "errquit.fh"
4#include "inp.fh"
5#include "rtdb.fh"
6#include "stdio.fh"
7#include "nwc_const.fh"
8#include "mafdecls.fh"
9#include "global.fh"
10#include "dimqm_constants.fh"
11c
12c     Input variables
13      integer, intent(in) ::  a1, a2
14      double precision, intent(in) :: dist, invdist
15      double precision, intent(in) :: pol1, pol2
16      double precision, intent(in) :: r(3)
17c
18c     Local variables
19      double precision dist3
20      double precision invdist3
21      double precision invdist5
22      double precision temp1, temp2
23      double precision Rp1, Rp2
24      double precision Rpp, invRpp
25      double precision aa, ab, ag
26      double precision bb, bg
27      double precision gg
28c
29c     Returns
30      double precision tr(3,3)
31      double precision util_erf
32      external util_erf
33c
34c     Compute powers and inverses once
35      dist3 = dist * dist * dist
36      invdist3 = invdist  * invdist * invdist
37      invdist5 = invdist3 * invdist * invdist
38c
39c     Calculate R values
40      Rp1   = (SQRTTWOINVPI * (pol1 * THIRD))**THIRD
41      Rp2   = (SQRTTWOINVPI * (pol2 * THIRD))**THIRD
42      Rpp   = SQRT(Rp1**2 + Rp2**2)
43      invRpp = ONE / Rpp
44c
45c     Calculate screening values
46      temp1 = util_erf(dist * invRpp) - TWOINVSQRTPI
47     $         * (dist * invRpp) * EXP(-(dist * invRpp)**2)
48      temp2 = FOUR / (SQRTPI * Rpp**3) * EXP(-(dist * invRpp)**2)
49c     Include these prefactors into the inverse distance terms
50      invdist5 = invdist5 * (temp1 * THREE - temp2 * dist3)
51      invdist3 = invdist3 * temp1
52c     Predetermine all possible values allowed by symmetry
53      aa = r(1)**2 * invdist5 - invdist3
54      bb = r(2)**2 * invdist5 - invdist3
55      gg = r(3)**2 * invdist5 - invdist3
56c
57      ab = r(1) * r(2) * invdist5
58      ag = r(1) * r(3) * invdist5
59      bg = r(2) * r(3) * invdist5
60c     Fill the temp array
61      tr(1,1) = aa
62      tr(2,1) = ab
63      tr(3,1) = ag
64      tr(1,2) = ab
65      tr(2,2) = bb
66      tr(3,2) = bg
67      tr(1,3) = ag
68      tr(2,3) = bg
69      tr(3,3) = gg
70      end function t2r
71
72      function t2c(r, dist, invdist, a1, a2, pol1, pol2) result (tc)
73      implicit none
74#include "errquit.fh"
75#include "inp.fh"
76#include "rtdb.fh"
77#include "stdio.fh"
78#include "nwc_const.fh"
79#include "mafdecls.fh"
80#include "global.fh"
81#include "dimqm_constants.fh"
82c
83c     Input variables
84      integer, intent(in) ::  a1, a2
85      double precision, intent(in) :: dist, invdist
86      double complex, intent(in) :: pol1, pol2
87      double precision, intent(in) :: r(3)
88c
89c     Local variables
90      double precision dist3
91      double precision invdist3
92      double precision invdist5
93      double precision temp1, temp2
94      double precision Rp1, Rp2
95      double precision Rpp, invRpp
96      double precision aa, ab, ag
97      double precision bb, bg
98      double precision gg
99      double precision util_erf
100      external util_erf
101c
102c     Returns
103      double complex tc(3,3)
104c
105c     Compute powers and inverses once
106      dist3 = dist * dist * dist
107      invdist3 = invdist  * invdist * invdist
108      invdist5 = invdist3 * invdist * invdist
109c
110c     Calculate R values
111      Rp1   = (SQRTTWOINVPI * (pol1 * THIRD))**THIRD
112      Rp2   = (SQRTTWOINVPI * (pol2 * THIRD))**THIRD
113      Rpp   = SQRT(Rp1**2 + Rp2**2)
114      invRpp = ONE / Rpp
115c
116c     Calculate screening values
117      temp1 = util_erf(dist * invRpp) - TWOINVSQRTPI
118     $         * (dist * invRpp) * EXP(-(dist * invRpp)**2)
119      temp2 = FOUR / (SQRTPI * Rpp**3) * EXP(-(dist * invRpp)**2)
120c     Include these prefactors into the inverse distance terms
121      invdist5 = invdist5 * (temp1 * THREE - temp2 * dist3)
122      invdist3 = invdist3 * temp1
123c     Predetermine all possible values allowed by symmetry
124      aa = r(1)**2 * invdist5 - invdist3
125      bb = r(2)**2 * invdist5 - invdist3
126      gg = r(3)**2 * invdist5 - invdist3
127c
128      ab = r(1) * r(2) * invdist5
129      ag = r(1) * r(3) * invdist5
130      bg = r(2) * r(3) * invdist5
131c     Fill the temp array
132      tc(1,1) = aa
133      tc(2,1) = ab
134      tc(3,1) = ag
135      tc(1,2) = ab
136      tc(2,2) = bb
137      tc(3,2) = bg
138      tc(1,3) = ag
139      tc(2,3) = bg
140      tc(3,3) = gg
141      end function t2c
142
143