./66_nonlocal/cont3.F90

ABINIT/cont3

ABINIT/cont3 [ Functions ]

NAME

 cont3

FUNCTION

 Compute several specialized contractions needed for the
 l=3 part of the stress tensor.

 Copyright (C) 1998-2018 ABINIT group (DCA, XG, GMR)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .
 For the initials of contributors, see ~abinit/doc/developers/contributors.txt.

INPUTS

  gxa(2,10)=complex symmetric rank 3 tensor
  gmet(3,3)=usual metric tensor, a symmetric matrix stored in
            full storage mode (bohr^-2)

OUTPUT

  rank2(6)=2*Re[contraction] given by
   2*Re[(15/2)*r3(a,i,j)*r3(b,j,i)-3*r1(i)*r3(a,b,i)-(3/2)*r1(a)*r1(b)]
   where r3(a,i,j)=gmet(j,k) gxa(a,i,k) and r1(a)=gmet(i,j) gxa(i,j,a).
  rank2 is stored in the compressed form 11 22 33 32 31 21.

NOTES

 Input gxa is a completely symmetric rank 3 tensor (complex)
 in compressed storage: 111 221 331 321 311 211 222 332 322 333.
 The output tensor is completely symmetric rank 2, real, and is given by
  $2 Re[{15 \over 2} r3(a,i,j) r3(b,j,i) - 3 r1(i) r3(a,b,i) - {3 \over 2} r1(a) r1(b)]$
  where $r3(a,i,j)=gmet(j,k) gxa(a,i,k)$ and $r1(a)=gmet(i,j) gxa(i,j,a)$.
  rank2 is stored in the compressed form 11 22 33 32 31 21.

PARENTS

      nonlop_pl

CHILDREN

SOURCE

 43 #if defined HAVE_CONFIG_H
 44 #include "config.h"
 45 #endif
 46 
 47 #include "abi_common.h"
 48 
 49 
 50 subroutine cont3(gxa,gmet,rank2)
 51 
 52  use defs_basis
 53 
 54 !This section has been created automatically by the script Abilint (TD).
 55 !Do not modify the following lines by hand.
 56 #undef ABI_FUNC
 57 #define ABI_FUNC 'cont3'
 58 !End of the abilint section
 59 
 60  implicit none
 61 
 62 !Arguments ------------------------------------
 63 !arrays
 64  real(dp),intent(in) :: gmet(3,3),gxa(2,10)
 65  real(dp),intent(out) :: rank2(6)
 66 
 67 !Local variables-------------------------------
 68 !scalars
 69  integer,parameter :: im=2,re=1
 70  integer :: ii
 71 !arrays
 72  real(dp) :: r1(2,3),r3(2,18),s13(6),s33(6)
 73 
 74 ! *************************************************************************
 75 
 76 !Compute r1(a) = gmet(i,j) gxa(i,j,a)
 77 
 78 !Write out components for 3 distinct terms, Re and Im
 79  do ii=1,2
 80    r1(ii,1)=gmet(1,1)*gxa(ii,1)+gmet(2,2)*gxa(ii,2)+&
 81 &   gmet(3,3)*gxa(ii,3)+2.d0*(&
 82 &   gmet(3,2)*gxa(ii,4)+gmet(3,1)*gxa(ii,5)+&
 83 &   gmet(2,1)*gxa(ii,6))
 84    r1(ii,2)=gmet(1,1)*gxa(ii,6)+gmet(2,2)*gxa(ii,7)+&
 85 &   gmet(3,3)*gxa(ii,8)+2.d0*(&
 86 &   gmet(3,2)*gxa(ii,9)+gmet(3,1)*gxa(ii,4)+&
 87 &   gmet(2,1)*gxa(ii,2))
 88    r1(ii,3)=gmet(1,1)*gxa(ii,5)+gmet(2,2)*gxa(ii,9)+&
 89 &   gmet(3,3)*gxa(ii,10)+2.d0*(&
 90 &   gmet(3,2)*gxa(ii,8)+gmet(3,1)*gxa(ii,3)+&
 91 &   gmet(2,1)*gxa(ii,4))
 92  end do
 93 
 94 !Compute r3(a,b,k)=gmet(k,n) gxa(a,b,n)
 95 
 96 !Write out components for 18 distinct terms, Re and Im
 97 !(symmetric in first two indices, not in all permutations)
 98 !store as 111 221 331 321 311 211
 99 !112 222 332 322 312 212
100 !113 223 333 323 313 213
101  do ii=1,2
102    r3(ii, 1)=gmet(1,1)*gxa(ii,1)+gmet(2,1)*gxa(ii,6)+&
103 &   gmet(3,1)*gxa(ii,5)
104    r3(ii, 2)=gmet(1,1)*gxa(ii,2)+gmet(2,1)*gxa(ii,7)+&
105 &   gmet(3,1)*gxa(ii,9)
106    r3(ii, 3)=gmet(1,1)*gxa(ii,3)+gmet(2,1)*gxa(ii,8)+&
107 &   gmet(3,1)*gxa(ii,10)
108    r3(ii, 4)=gmet(1,1)*gxa(ii,4)+gmet(2,1)*gxa(ii,9)+&
109 &   gmet(3,1)*gxa(ii,8)
110    r3(ii, 5)=gmet(1,1)*gxa(ii,5)+gmet(2,1)*gxa(ii,4)+&
111 &   gmet(3,1)*gxa(ii,3)
112    r3(ii, 6)=gmet(1,1)*gxa(ii,6)+gmet(2,1)*gxa(ii,2)+&
113 &   gmet(3,1)*gxa(ii,4)
114    r3(ii, 7)=gmet(2,1)*gxa(ii,1)+gmet(2,2)*gxa(ii,6)+&
115 &   gmet(3,2)*gxa(ii,5)
116    r3(ii, 8)=gmet(2,1)*gxa(ii,2)+gmet(2,2)*gxa(ii,7)+&
117 &   gmet(3,2)*gxa(ii,9)
118    r3(ii, 9)=gmet(2,1)*gxa(ii,3)+gmet(2,2)*gxa(ii,8)+&
119 &   gmet(3,2)*gxa(ii,10)
120    r3(ii,10)=gmet(2,1)*gxa(ii,4)+gmet(2,2)*gxa(ii,9)+&
121 &   gmet(3,2)*gxa(ii,8)
122    r3(ii,11)=gmet(2,1)*gxa(ii,5)+gmet(2,2)*gxa(ii,4)+&
123 &   gmet(3,2)*gxa(ii,3)
124    r3(ii,12)=gmet(2,1)*gxa(ii,6)+gmet(2,2)*gxa(ii,2)+&
125 &   gmet(3,2)*gxa(ii,4)
126    r3(ii,13)=gmet(3,1)*gxa(ii,1)+gmet(3,2)*gxa(ii,6)+&
127 &   gmet(3,3)*gxa(ii,5)
128    r3(ii,14)=gmet(3,1)*gxa(ii,2)+gmet(3,2)*gxa(ii,7)+&
129 &   gmet(3,3)*gxa(ii,9)
130    r3(ii,15)=gmet(3,1)*gxa(ii,3)+gmet(3,2)*gxa(ii,8)+&
131 &   gmet(3,3)*gxa(ii,10)
132    r3(ii,16)=gmet(3,1)*gxa(ii,4)+gmet(3,2)*gxa(ii,9)+&
133 &   gmet(3,3)*gxa(ii,8)
134    r3(ii,17)=gmet(3,1)*gxa(ii,5)+gmet(3,2)*gxa(ii,4)+&
135 &   gmet(3,3)*gxa(ii,3)
136    r3(ii,18)=gmet(3,1)*gxa(ii,6)+gmet(3,2)*gxa(ii,2)+&
137 &   gmet(3,3)*gxa(ii,4)
138 
139  end do
140 
141 !Now need
142 !2*Re[(15/2)*r3(a,i,j)*r3(b,j,i)-3*r1(i)*r3(a,b,i)-(3/2)*r1(a)*r1(b)].
143 
144 !Write out s33(a,b)=2*Re[r3(a,i,j)*r3(b,j,i)]
145 
146  s33(1)=2.d0*(r3(re, 1)*r3(re, 1)+r3(im, 1)*r3(im, 1)+&
147 & r3(re,12)*r3(re,12)+r3(im,12)*r3(im,12)+&
148 & r3(re,17)*r3(re,17)+r3(im,17)*r3(im,17)+&
149 & r3(re,11)*r3(re,18)+r3(im,11)*r3(im,18)+&
150 & r3(re,18)*r3(re,11)+r3(im,18)*r3(im,11)+&
151 & r3(re, 5)*r3(re,13)+r3(im, 5)*r3(im,13)+&
152 & r3(re,13)*r3(re, 5)+r3(im,13)*r3(im, 5)+&
153 & r3(re, 6)*r3(re, 7)+r3(im, 6)*r3(im, 7)+&
154 & r3(re, 7)*r3(re, 6)+r3(im, 7)*r3(im, 6))
155 
156  s33(2)=2.d0*(r3(re, 6)*r3(re, 6)+r3(im, 6)*r3(im, 6)+&
157 & r3(re, 8)*r3(re, 8)+r3(im, 8)*r3(im, 8)+&
158 & r3(re,16)*r3(re,16)+r3(im,16)*r3(im,16)+&
159 & r3(re,10)*r3(re,14)+r3(im,10)*r3(im,14)+&
160 & r3(re,14)*r3(re,10)+r3(im,14)*r3(im,10)+&
161 & r3(re, 4)*r3(re,18)+r3(im, 4)*r3(im,18)+&
162 & r3(re,18)*r3(re, 4)+r3(im,18)*r3(im, 4)+&
163 & r3(re, 2)*r3(re,12)+r3(im, 2)*r3(im,12)+&
164 & r3(re,12)*r3(re, 2)+r3(im,12)*r3(im, 2))
165 
166  s33(3)=2.d0*(r3(re, 5)*r3(re, 5)+r3(im, 5)*r3(im, 5)+&
167 & r3(re,10)*r3(re,10)+r3(im,10)*r3(im,10)+&
168 & r3(re,15)*r3(re,15)+r3(im,15)*r3(im,15)+&
169 & r3(re, 9)*r3(re,16)+r3(im, 9)*r3(im,16)+&
170 & r3(re,16)*r3(re, 9)+r3(im,16)*r3(im, 9)+&
171 & r3(re, 3)*r3(re,17)+r3(im, 3)*r3(im,17)+&
172 & r3(re,17)*r3(re, 3)+r3(im,17)*r3(im, 3)+&
173 & r3(re, 4)*r3(re,11)+r3(im, 4)*r3(im,11)+&
174 & r3(re,11)*r3(re, 4)+r3(im,11)*r3(im, 4))
175 
176  s33(4)=2.d0*(r3(re, 5)*r3(re, 6)+r3(im, 5)*r3(im, 6)+&
177 & r3(re,10)*r3(re, 8)+r3(im,10)*r3(im, 8)+&
178 & r3(re,15)*r3(re,16)+r3(im,15)*r3(im,16)+&
179 & r3(re, 9)*r3(re,14)+r3(im, 9)*r3(im,14)+&
180 & r3(re,16)*r3(re,10)+r3(im,16)*r3(im,10)+&
181 & r3(re, 3)*r3(re,18)+r3(im, 3)*r3(im,18)+&
182 & r3(re,17)*r3(re, 4)+r3(im,17)*r3(im, 4)+&
183 & r3(re, 4)*r3(re,12)+r3(im, 4)*r3(im,12)+&
184 & r3(re,11)*r3(re, 2)+r3(im,11)*r3(im, 2))
185 
186  s33(5)=2.d0*(r3(re, 5)*r3(re, 1)+r3(im, 5)*r3(im, 1)+&
187 & r3(re,10)*r3(re,12)+r3(im,10)*r3(im,12)+&
188 & r3(re,15)*r3(re,17)+r3(im,15)*r3(im,17)+&
189 & r3(re, 9)*r3(re,18)+r3(im, 9)*r3(im,18)+&
190 & r3(re,16)*r3(re,11)+r3(im,16)*r3(im,11)+&
191 & r3(re, 3)*r3(re,13)+r3(im, 3)*r3(im,13)+&
192 & r3(re,17)*r3(re, 5)+r3(im,17)*r3(im, 5)+&
193 & r3(re, 4)*r3(re, 7)+r3(im, 4)*r3(im, 7)+&
194 & r3(re,11)*r3(re, 6)+r3(im,11)*r3(im, 6))
195 
196  s33(6)=2.d0*(r3(re, 6)*r3(re, 1)+r3(im, 6)*r3(im, 1)+&
197 & r3(re, 8)*r3(re,12)+r3(im, 8)*r3(im,12)+&
198 & r3(re,16)*r3(re,17)+r3(im,16)*r3(im,17)+&
199 & r3(re,10)*r3(re,18)+r3(im,10)*r3(im,18)+&
200 & r3(re,14)*r3(re,11)+r3(im,14)*r3(im,11)+&
201 & r3(re, 4)*r3(re,13)+r3(im, 4)*r3(im,13)+&
202 & r3(re,18)*r3(re, 5)+r3(im,18)*r3(im, 5)+&
203 & r3(re, 2)*r3(re, 7)+r3(im, 2)*r3(im, 7)+&
204 & r3(re,12)*r3(re, 6)+r3(im,12)*r3(im, 6))
205 
206 
207 !Write out s13(a,b)=2*Re[r1(i)*r3(a,b,i)]
208 
209  s13(1)=2.d0*(r1(re,1)*r3(re, 1)+r1(im,1)*r3(im, 1)+&
210 & r1(re,2)*r3(re, 7)+r1(im,2)*r3(im, 7)+&
211 & r1(re,3)*r3(re,13)+r1(im,3)*r3(im,13))
212  s13(2)=2.d0*(r1(re,1)*r3(re, 2)+r1(im,1)*r3(im, 2)+&
213 & r1(re,2)*r3(re, 8)+r1(im,2)*r3(im, 8)+&
214 & r1(re,3)*r3(re,14)+r1(im,3)*r3(im,14))
215  s13(3)=2.d0*(r1(re,1)*r3(re, 3)+r1(im,1)*r3(im, 3)+&
216 & r1(re,2)*r3(re, 9)+r1(im,2)*r3(im, 9)+&
217 & r1(re,3)*r3(re,15)+r1(im,3)*r3(im,15))
218  s13(4)=2.d0*(r1(re,1)*r3(re, 4)+r1(im,1)*r3(im, 4)+&
219 & r1(re,2)*r3(re,10)+r1(im,2)*r3(im,10)+&
220 & r1(re,3)*r3(re,16)+r1(im,3)*r3(im,16))
221  s13(5)=2.d0*(r1(re,1)*r3(re, 5)+r1(im,1)*r3(im, 5)+&
222 & r1(re,2)*r3(re,11)+r1(im,2)*r3(im,11)+&
223 & r1(re,3)*r3(re,17)+r1(im,3)*r3(im,17))
224  s13(6)=2.d0*(r1(re,1)*r3(re, 6)+r1(im,1)*r3(im, 6)+&
225 & r1(re,2)*r3(re,12)+r1(im,2)*r3(im,12)+&
226 & r1(re,3)*r3(re,18)+r1(im,3)*r3(im,18))
227 
228 !Finally, write out the six terms as final answer
229 !rank2(a,b)=(15/2)*s33(a,b)-3*s13(a,b)-(3/2)*2*Re[r1(a)*r1(b)]
230 
231  rank2(1)=7.5d0*s33(1)-3.d0*s13(1)&
232 & -3.d0*(r1(re,1)*r1(re,1)+r1(im,1)*r1(im,1))
233  rank2(2)=7.5d0*s33(2)-3.d0*s13(2)&
234 & -3.d0*(r1(re,2)*r1(re,2)+r1(im,2)*r1(im,2))
235  rank2(3)=7.5d0*s33(3)-3.d0*s13(3)&
236 & -3.d0*(r1(re,3)*r1(re,3)+r1(im,3)*r1(im,3))
237  rank2(4)=7.5d0*s33(4)-3.d0*s13(4)&
238 & -3.d0*(r1(re,3)*r1(re,2)+r1(im,3)*r1(im,2))
239  rank2(5)=7.5d0*s33(5)-3.d0*s13(5)&
240 & -3.d0*(r1(re,3)*r1(re,1)+r1(im,3)*r1(im,1))
241  rank2(6)=7.5d0*s33(6)-3.d0*s13(6)&
242 & -3.d0*(r1(re,2)*r1(re,1)+r1(im,2)*r1(im,1))
243 
244 end subroutine cont3

TABLE OF CONTENTS

ABINIT/cont3 [ Functions ]