TABLE OF CONTENTS
ABINIT/m_opernlc_ylm_allwf_cpu [ Modules ]
NAME
m_opernlc_ylm_allwf_cpu
FUNCTION
COPYRIGHT
Copyright (C) 2008-2024 ABINIT group (MT) This file is distributed under the terms of the GNU General Public License, see ~abinit/COPYING or http://www.gnu.org/copyleft/gpl.txt .
PARENTS
CHILDREN
SOURCE
19 #if defined HAVE_CONFIG_H 20 #include "config.h" 21 #endif 22 23 #include "abi_common.h" 24 25 module m_opernlc_ylm_allwf_cpu 26 27 use defs_basis 28 use m_errors 29 use m_abicore 30 use m_xmpi 31 32 use defs_abitypes, only : MPI_type 33 34 implicit none 35 36 private
ABINIT/opernlc_ylm_allwf_cpu [ Functions ]
NAME
opernlc_ylm_allwf_cpu
FUNCTION
* Operate with the non-local part of the hamiltonian, in order to reduce projected scalars * Operate with the non-local projectors and the overlap matrix, in order to reduce projected scalars
INPUTS
atindx1(natom)=index table for atoms (gives the absolute index of
an atom from its rank in a block of atoms)
cplex=1 if <p_lmn|c> scalars are real (equivalent to istwfk>1)
2 if <p_lmn|c> scalars are complex
cplex_dgxdt(ndgxdt) = used only when cplex = 1
cplex_dgxdt(i) = 1 if dgxdt(1,i,:,:) is real, 2 if it is pure imaginary
cplex_enl=1 if enl factors are real, 2 if they are complex
cplex_fac=1 if gxfac scalars are real, 2 if gxfac scalars are complex
dgxdt(cplex,ndgxdt,nlmn,nincat)=grads of projected scalars (only if optder>0)
dimenl1,dimenl2=dimensions of enl (see enl)
dimekbq=1 if enl factors do not contain a exp(-iqR) phase, 2 is they do
enl(cplex_enl*dimenl1,dimenl2,nspinortot**2,dimekbq)=
->Norm conserving : ==== when paw_opt=0 ====
(Real) Kleinman-Bylander energies (hartree)
dimenl1=lmnmax - dimenl2=ntypat
dimekbq is 2 if Enl contains a exp(-iqR) phase, 1 otherwise
->PAW : ==== when paw_opt=1, 2 or 4 ====
(Real or complex, hermitian) Dij coefs to connect projectors
dimenl1=cplex_enl*lmnmax*(lmnmax+1)/2 - dimenl2=natom
These are complex numbers if cplex_enl=2
enl(:,:,1) contains Dij^up-up
enl(:,:,2) contains Dij^dn-dn
enl(:,:,3) contains Dij^up-dn (only if nspinor=2)
enl(:,:,4) contains Dij^dn-up (only if nspinor=2)
dimekbq is 2 if Dij contains a exp(-iqR) phase, 1 otherwise
gx(cplex,nlmn,nincat*abs(enl_opt))= projected scalars
iatm=absolute rank of first atom of the current block of atoms
indlmn(6,nlmn)= array giving l,m,n,lm,ln,s for i=lmn
itypat=type of atoms
lambda=factor to be used when computing (Vln-lambda.S) - only for paw_opt=2
mpi_enreg=information about MPI parallelization
natom=number of atoms in cell
ndgxdt=second dimension of dgxdt
ndgxdtfac=second dimension of dgxdtfac
nincat=number of atoms in the subset here treated
nlmn=number of (l,m,n) numbers for current type of atom
nspinor= number of spinorial components of the wavefunctions (on current proc)
nspinortot=total number of spinorial components of the wavefunctions
optder=0=only gxfac is computed, 1=both gxfac and dgxdtfac are computed
2=gxfac, dgxdtfac and d2gxdtfac are computed
paw_opt= define the nonlocal operator concerned with:
paw_opt=0 : Norm-conserving Vnl (use of Kleinman-Bylander ener.)
paw_opt=1 : PAW nonlocal part of H (use of Dij coeffs)
paw_opt=2 : PAW: (Vnl-lambda.Sij) (Sij=overlap matrix)
paw_opt=3 : PAW overlap matrix (Sij)
paw_opt=4 : both PAW nonlocal part of H (Dij) and overlap matrix (Sij)
sij(nlm*(nlmn+1)/2)=overlap matrix components (only if paw_opt=2, 3 or 4)
OUTPUT
if (paw_opt=0, 1, 2 or 4)
gxfac(cplex_fac,nlmn,nincat,nspinor)= reduced projected scalars related to Vnl (NL operator)
if (paw_opt=3 or 4)
gxfac_sij(cplex,nlmn,nincat,nspinor)= reduced projected scalars related to Sij (overlap)
if (optder==1.and.paw_opt=0, 1, 2 or 4)
dgxdtfac(cplex_fac,ndgxdtfac,nlmn,nincat,nspinor)= gradients of gxfac related to Vnl (NL operator)
if (optder==1.and.paw_opt=3 or 4)
dgxdtfac_sij(cplex,ndgxdtfac,nlmn,nincat,nspinor)= gradients of gxfac related to Sij (overlap)
NOTES
This routine operates for one type of atom, and within this given type of atom, for a subset of at most nincat atoms. About the non-local factors symmetry: - The lower triangular part of the Dij matrix can be deduced from the upper one with the following relation: D^s2s1_ji = (D^s1s2_ij)^* where s1,s2 are spinor components - The Dij factors can contain a exp(-iqR) phase This phase does not have to be included in the symmetry rule For that reason, we first apply the real part (cos(qR).D^s1s2_ij) then, we apply the imaginary part (-sin(qR).D^s1s2_ij)
PARENTS
m_gemm_nonlop,m_nonlop_ylm
CHILDREN
xmpi_sum
SOURCE
135 subroutine opernlc_ylm_allwf_cpu(atindx1,cplex,cplex_enl,cplex_fac, & 136 & dimenl1,dimenl2,dimekbq,enl, & 137 & gx,gxfac,gxfac_sij, & 138 & iatm,indlmn,itypat,ndat,lambda,mpi_enreg,natom, & 139 & nincat,nlmn,nspinor,nspinortot,paw_opt,sij) 140 141 !Arguments ------------------------------------ 142 !scalars 143 integer,intent(in) :: cplex,cplex_enl,cplex_fac,dimenl1,dimenl2,dimekbq,iatm,itypat 144 integer,intent(in) :: natom,nincat,nspinor,nspinortot,paw_opt 145 integer,intent(inout) :: nlmn 146 type(MPI_type) , intent(in) :: mpi_enreg 147 integer,intent(in) :: ndat 148 !arrays 149 integer, intent(in) :: atindx1(natom),indlmn(6,nlmn) 150 real(dp),intent(in),target :: enl(dimenl1,dimenl2,nspinortot**2,dimekbq) 151 real(dp),intent(inout) :: gx(cplex,nlmn,nincat,nspinor*ndat) 152 real(dp),intent(in) :: sij(((paw_opt+1)/3)*nlmn*(nlmn+1)/2) 153 real(dp),intent(out),target :: gxfac(cplex_fac,nlmn,nincat,nspinor*ndat) 154 real(dp),intent(out) :: gxfac_sij(cplex,nlmn,nincat,nspinor*ndat*(paw_opt/3)) 155 real(dp),intent(in) :: lambda(ndat) 156 157 !Local variables------------------------------- 158 !Arrays 159 !scalars 160 integer :: ia, idat, idat_ispinor, idat_jspinor 161 integer :: ijlmn,ilm,ilmn,i0lmn,iln,index_enl,iphase,ispinor,ispinor_index 162 integer :: j0lmn,jjlmn,jlm,jlmn,jspinor,shift 163 !integer :: ierr,ijspin,jilmn,jispin,jspinor_index FIXME Unneeded ? 164 real(dp) :: sijr 165 !arrays 166 real(dp) :: enl_(2), gxi(cplex), gxj(cplex) 167 !real(dp),allocatable :: gxfac_offdiag(:,:,:,:) !FIXME Unneeded 168 real(dp),pointer :: gxfac_(:,:,:,:) 169 real(dp),pointer :: enl_ptr(:,:,:) 170 171 ! ************************************************************************* 172 173 DBG_ENTER("COLL") 174 175 !Parallelization over spinors treatment 176 shift=0 177 if (mpi_enreg%paral_spinor==1) then 178 shift = mpi_enreg%me_spinor 179 end if 180 181 !When Enl factors contain a exp(-iqR) phase: 182 ! - We loop over the real and imaginary parts 183 ! - We need an additional memory space 184 do iphase=1,dimekbq 185 186 if (paw_opt==3) cycle 187 188 if (iphase==1) then 189 gxfac_ => gxfac 190 else 191 ABI_CHECK(cplex_fac==2,"BUG: invalid cplex_fac==1 when dimekbq=2!") 192 ABI_MALLOC(gxfac_,(cplex_fac,nlmn,nincat,ndat*nspinor)) 193 end if 194 195 gxfac_ = zero 196 enl_ptr => enl(:,:,:,iphase) 197 198 199 !Accumulate gxfac related to non-local operator (Norm-conserving) 200 !------------------------------------------------------------------- 201 if (paw_opt==0) then 202 !Enl is E(Kleinman-Bylander) 203 ABI_CHECK(cplex_enl/=2,"BUG: invalid cplex_enl=2!") 204 ABI_CHECK(cplex_fac==cplex,"BUG: invalid cplex_fac/=cplex!") 205 do idat = 1,ndat 206 do ispinor = 1,nspinor 207 do ia = 1,nincat 208 do ilmn = 1,nlmn 209 ispinor_index = ispinor+shift 210 idat_ispinor = nspinor*(idat-1) + ispinor 211 iln = indlmn(5,ilmn) 212 enl_(1) = enl_ptr(iln,itypat,ispinor_index) 213 gxfac_(1:cplex,ilmn,ia,idat_ispinor) = enl_(1)*gx(1:cplex,ilmn,ia,idat_ispinor) 214 end do ! ilmn 215 end do ! ia 216 end do ! ispinor 217 end do ! idat 218 end if ! paw_opt == 0 219 220 !Accumulate gxfac related to nonlocal operator (PAW) 221 !------------------------------------------------------------------- 222 if (paw_opt==1.or.paw_opt==2.or.paw_opt==4) then 223 224 !Enl is psp strength Dij or (Dij-lambda.Sij) 225 226 ! === Diagonal term(s) (up-up, down-down) 227 228 if (cplex_enl==1) then 229 ! 1-Enl is real 230 231 do idat = 1,ndat 232 do ispinor=1,nspinor 233 do ia=1,nincat 234 do jlmn=1,nlmn 235 236 ispinor_index = ispinor+shift 237 index_enl = atindx1(iatm+ia) 238 j0lmn = jlmn*(jlmn-1)/2 239 jjlmn = j0lmn+jlmn 240 enl_(1) = enl_ptr(jjlmn,index_enl,ispinor_index) 241 idat_ispinor = nspinor*(idat-1) + ispinor 242 243 if (paw_opt==2) enl_(1) = enl_(1)-lambda(idat)*sij(jjlmn) 244 gxj (1:cplex) = gx (1:cplex,jlmn,ia,idat_ispinor) 245 gxfac_(1:cplex,jlmn,ia,idat_ispinor) = gxfac_(1:cplex,jlmn,ia,idat_ispinor) + enl_(1)*gxj(1:cplex) 246 247 do ilmn=1,nlmn 248 if (ilmn < jlmn) then 249 ijlmn = j0lmn+ilmn 250 enl_(1) = enl_ptr(ijlmn,index_enl,ispinor_index) 251 if (paw_opt==2) enl_(1)=enl_(1)-lambda(idat)*sij(ijlmn) 252 gxi (1:cplex) = gx (1:cplex,ilmn,ia,idat_ispinor) 253 gxfac_(1:cplex,jlmn,ia,idat_ispinor) = gxfac_(1:cplex,jlmn,ia,idat_ispinor) + enl_(1)*gxi(1:cplex) 254 else if (ilmn > jlmn) then 255 if(jlmn<nlmn) then 256 i0lmn = (ilmn*(ilmn-1)/2) 257 ijlmn = i0lmn+jlmn 258 enl_(1) = enl_ptr(ijlmn,index_enl,ispinor_index) 259 if (paw_opt==2) enl_(1)=enl_(1)-lambda(idat)*sij(ijlmn) 260 gxi (1:cplex) = gx (1:cplex,ilmn,ia,idat_ispinor) 261 gxfac_(1:cplex,jlmn,ia,idat_ispinor) = gxfac_(1:cplex,jlmn,ia,idat_ispinor) + enl_(1)*gxi(1:cplex) 262 end if 263 end if 264 end do ! ilmn 265 266 end do ! jlmn 267 end do ! ia 268 end do ! ispinor 269 end do ! idat 270 271 else 272 ! 2-Enl is complex ===== D^ss'_ij=D^s's_ji^* 273 274 ABI_CHECK(cplex_fac==cplex_enl,"BUG: invalid cplex_fac/=cplex_enl!") 275 276 if (nspinortot==1) then ! -------------> NO SPINORS 277 278 do idat = 1,ndat 279 280 ! ispinor is 1 281 ! nspinor should be 1 too (maybe an ABI_CHECK could ensure that) 282 idat_ispinor = nspinor*(idat-1) + 1 283 284 do ia=1,nincat 285 286 index_enl = atindx1(iatm+ia) 287 288 do jlmn=1,nlmn 289 290 j0lmn = jlmn*(jlmn-1)/2 291 jjlmn = j0lmn+jlmn 292 293 enl_(1) = enl_ptr(2*jjlmn-1,index_enl,1) 294 if (paw_opt==2) then 295 enl_(1) = enl_(1)-lambda(idat)*sij(jjlmn) 296 end if 297 298 gxj(1:cplex) = gx(1:cplex,jlmn,ia,idat_ispinor) 299 300 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) + enl_(1)*gxj(1) 301 if (cplex==2) then 302 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(1)*gxj(2) 303 end if 304 305 do ilmn=1,jlmn-1 306 307 if (ilmn < jlmn) then 308 309 ijlmn = j0lmn + ilmn 310 311 enl_(1:2) = enl_ptr(2*ijlmn-1:2*ijlmn,index_enl,1) 312 if (paw_opt==2) enl_(1)=enl_(1)-lambda(idat)*sij(ijlmn) 313 314 gxi(1:cplex) = gx(1:cplex,ilmn,ia,idat_ispinor) 315 316 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) + enl_(1)*gxi(1) 317 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) - enl_(2)*gxi(1) 318 319 if (cplex==2) then 320 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) + enl_(2)*gxi(2) 321 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(1)*gxi(2) 322 end if 323 324 else if (ilmn > jlmn) then 325 326 if(jlmn<nlmn) then 327 i0lmn=ilmn*(ilmn-1)/2 328 ijlmn=i0lmn+jlmn 329 330 enl_(1:2)=enl_ptr(2*ijlmn-1:2*ijlmn,index_enl,1) 331 if (paw_opt==2) enl_(1)=enl_(1)-lambda(idat)*sij(ijlmn) 332 333 gxi(1:cplex) = gx(1:cplex,ilmn,ia,idat_ispinor) 334 335 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) + enl_(1)*gxi(1) 336 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(2)*gxi(1) 337 338 if (cplex==2) then 339 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) - enl_(2)*gxi(2) 340 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(1)*gxi(2) 341 end if ! cplex == 2 342 343 end if ! jlmn 344 345 end if ! ilmn 346 347 end do ! ilmn 348 349 end do ! jlmn 350 351 end do ! ia 352 353 end do ! idat 354 355 else 356 ! -------------> SPINORIAL CASE 357 358 ! === Diagonal term(s) (up-up, down-down) 359 360 do idat = 1,ndat 361 362 do ispinor = 1,nspinor 363 364 ispinor_index = ispinor+shift 365 idat_ispinor = nspinor*(idat-1) + ispinor 366 367 do ia=1,nincat 368 369 index_enl = atindx1(iatm+ia) 370 371 do jlmn = 1,nlmn 372 373 j0lmn = jlmn*(jlmn-1)/2 374 jjlmn = j0lmn+jlmn 375 376 enl_(1) = enl_ptr(2*jjlmn-1,index_enl,ispinor_index) 377 if (paw_opt==2) then 378 enl_(1) = enl_(1)-lambda(idat)*sij(jjlmn) 379 end if 380 381 gxj(1:cplex) = gx(1:cplex,jlmn,ia,idat_ispinor) 382 383 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) + enl_(1)*gxj(1) 384 if (cplex==2) then 385 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(1)*gxj(2) 386 end if 387 388 do ilmn = 1,nlmn 389 390 if (ilmn < jlmn) then 391 ijlmn = j0lmn+ilmn 392 393 enl_(1:2) = enl_ptr(2*ijlmn-1:2*ijlmn,index_enl,ispinor_index) 394 if (paw_opt==2) then 395 enl_(1) = enl_(1)-lambda(idat)*sij(ijlmn) 396 end if 397 398 gxi(1:cplex) = gx(1:cplex,ilmn,ia,idat_ispinor) 399 400 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) + enl_(1)*gxi(1) 401 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) - enl_(2)*gxi(1) 402 403 if (cplex==2) then 404 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) + enl_(2)*gxi(2) 405 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(1)*gxi(2) 406 end if 407 408 else if (ilmn > jlmn) then 409 410 if(jlmn<nlmn) then 411 i0lmn=ilmn*(ilmn-1)/2 412 ijlmn=i0lmn+jlmn 413 enl_(1:2) = enl_ptr(2*ijlmn-1:2*ijlmn,index_enl,ispinor_index) 414 415 if (paw_opt==2) enl_(1)=enl_(1)-lambda(idat)*sij(ijlmn) 416 gxi(1:cplex) = gx(1:cplex,ilmn,ia,idat_ispinor) 417 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) + enl_(1)*gxi(1) 418 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(2)*gxi(1) 419 if (cplex==2) then 420 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) - enl_(2)*gxi(2) 421 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(1)*gxi(2) 422 end if ! cplex == 2 423 end if 424 425 end if ! ilmn 426 427 end do ! ilmn 428 429 end do ! jlmn 430 431 end do ! ia 432 433 end do ! ispinor 434 435 end do ! idat 436 437 end if ! nspinortot 438 439 end if ! complex_enl 440 441 ! === Off-diagonal term(s) (up-down, down-up) 442 443 ! --- No parallelization over spinors --- 444 if (nspinortot==2 .and. nspinor==nspinortot) then 445 446 ABI_CHECK(cplex_enl==2,"BUG: invalid cplex_enl/=2!") 447 ABI_CHECK(cplex_fac==cplex,"BUG: invalid cplex_fac/=cplex)!") 448 449 do idat = 1,ndat 450 451 do ispinor=1,nspinortot 452 jspinor = 3-ispinor 453 454 idat_ispinor = nspinor*(idat-1) + ispinor 455 456 ! watch out : ispinor replaced by jspinor 457 idat_jspinor = nspinor*(idat-1) + jspinor 458 459 do ia=1,nincat 460 index_enl=atindx1(iatm+ia) 461 462 do jlmn=1,nlmn 463 j0lmn=jlmn*(jlmn-1)/2 464 jjlmn=j0lmn+jlmn 465 enl_(1:2) = enl_ptr(2*jjlmn-1:2*jjlmn,index_enl,2+ispinor ) 466 gxi(1:cplex) = gx(1:cplex,jlmn,ia,idat_ispinor) 467 gxfac_(1,jlmn,ia,idat_jspinor) = gxfac_(1,jlmn,ia,idat_jspinor) + enl_(1)*gxi(1) 468 gxfac_(2,jlmn,ia,idat_jspinor) = gxfac_(2,jlmn,ia,idat_jspinor) - enl_(2)*gxi(1) 469 if (cplex==2) then 470 gxfac_(1,jlmn,ia,idat_jspinor) = gxfac_(1,jlmn,ia,idat_jspinor) + enl_(2)*gxi(2) 471 gxfac_(2,jlmn,ia,idat_jspinor) = gxfac_(2,jlmn,ia,idat_jspinor) + enl_(1)*gxi(2) 472 end if 473 474 do ilmn=1,nlmn 475 476 if (ilmn < jlmn) then 477 478 ijlmn=j0lmn+ilmn 479 enl_(1:2)=enl_ptr(2*ijlmn-1:2*ijlmn,index_enl,2+ispinor) 480 gxi(1:cplex) = gx(1:cplex,ilmn,ia,idat_ispinor) 481 gxfac_(1,jlmn,ia,idat_jspinor) = gxfac_(1,jlmn,ia,idat_jspinor) + enl_(1)*gxi(1) 482 gxfac_(2,jlmn,ia,idat_jspinor) = gxfac_(2,jlmn,ia,idat_jspinor) - enl_(2)*gxi(1) 483 if (cplex==2) then 484 gxfac_(1,jlmn,ia,idat_jspinor) = gxfac_(1,jlmn,ia,idat_jspinor) + enl_(2)*gxi(2) 485 gxfac_(2,jlmn,ia,idat_jspinor) = gxfac_(2,jlmn,ia,idat_jspinor) + enl_(1)*gxi(2) 486 end if 487 488 else if (ilmn > jlmn) then 489 490 if(jlmn<nlmn) then 491 i0lmn=ilmn*(ilmn-1)/2 492 ijlmn=i0lmn+jlmn 493 enl_(1:2)=enl_ptr(2*ijlmn-1:2*ijlmn,index_enl,2+ispinor) 494 gxi(1:cplex)=gx(1:cplex,ilmn,ia,idat_jspinor) 495 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) + enl_(1)*gxi(1) 496 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(2)*gxi(1) 497 if (cplex==2) then 498 gxfac_(1,jlmn,ia,idat_ispinor) = gxfac_(1,jlmn,ia,idat_ispinor) - enl_(2)*gxi(2) 499 gxfac_(2,jlmn,ia,idat_ispinor) = gxfac_(2,jlmn,ia,idat_ispinor) + enl_(1)*gxi(2) 500 end if ! cplex==2 501 502 end if ! if jlmn 503 504 end if ! if ilmn 505 506 end do ! ilmn 507 508 end do ! jlmn 509 510 end do ! ia 511 512 end do ! ispinor 513 514 end do ! idat 515 516 ! --- Parallelization over spinors --- 517 518 ! 519 ! this case is removed for GPU parallelization, npspinor>1 is forbidden when gpu_option is ON 520 ! 521 522 ! else if (nspinortot==2 .and. nspinor/=nspinortot) then 523 524 ! ABI_CHECK(cplex_enl==2,"BUG: invalid cplex_enl/=2!") 525 ! ABI_CHECK(cplex_fac==2,"BUG: invalid cplex_fac/=2!") 526 ! ABI_MALLOC(gxfac_offdiag,(cplex_fac,nlmn,nincat,nspinortot)) 527 528 ! gxfac_offdiag(:,:,:,:) = zero 529 530 ! ispinor_index = mpi_enreg%me_spinor+1 531 ! jspinor_index = 3-ispinor_index 532 ! if (ispinor_index==1) then 533 ! ijspin=3 534 ! jispin=4 535 ! else 536 ! ijspin=4 537 ! jispin=3 538 ! end if 539 540 ! do ia=1,nincat 541 ! index_enl = atindx1(iatm+ia) 542 543 ! do jlmn=1,nlmn 544 ! j0lmn=jlmn*(jlmn-1)/2 545 546 ! do ilmn=1,nlmn 547 ! i0lmn=ilmn*(ilmn-1)/2 548 549 ! if (ilmn<=jlmn) then 550 ! ijlmn=j0lmn+ilmn 551 ! enl_(1) = enl_ptr(2*ijlmn-1,index_enl,ijspin) 552 ! enl_(2) = -enl_ptr(2*ijlmn ,index_enl,ijspin) 553 ! else 554 ! jilmn=i0lmn+jlmn 555 ! enl_(1) = enl_ptr(2*jilmn-1,index_enl,jispin) 556 ! enl_(2) = enl_ptr(2*jilmn ,index_enl,jispin) 557 ! end if 558 559 ! gxi(1:cplex) = gx(1:cplex,ilmn,ia,1) 560 561 ! gxfac_offdiag(1,jlmn,ia,jspinor_index) = & 562 ! & gxfac_offdiag(1,jlmn,ia,jspinor_index) + enl_(1)*gxi(1) 563 564 ! gxfac_offdiag(2,jlmn,ia,jspinor_index) = & 565 ! & gxfac_offdiag(2,jlmn,ia,jspinor_index) + enl_(2)*gxi(1) 566 567 ! if (cplex==2) then 568 ! gxfac_offdiag(1,jlmn,ia,jspinor_index) = & 569 ! & gxfac_offdiag(1,jlmn,ia,jspinor_index) - enl_(2)*gxi(2) 570 ! gxfac_offdiag(2,jlmn,ia,jspinor_index) = & 571 ! & gxfac_offdiag(2,jlmn,ia,jspinor_index) + enl_(1)*gxi(2) 572 ! end if 573 574 ! end do !ilmn 575 576 ! end do !jlmn 577 578 ! end do !iat 579 580 ! !call xmpi_sum(gxfac_offdiag,mpi_enreg%comm_spinor,ierr) 581 582 ! gxfac_(:,:,:,1) = gxfac_(:,:,:,1) + gxfac_offdiag(:,:,:,ispinor_index) 583 ! ABI_FREE(gxfac_offdiag) 584 585 end if 586 587 end if !paw_opt 588 589 !End of loop when a exp(-iqR) phase is present 590 !------------------------------------------- ------------------------ 591 592 !When iphase=1, gxfac and gxfac_ point to the same memory space 593 !When iphase=2, we add i.gxfac_ to gxfac 594 if (iphase==2) then 595 do idat = 1,ndat 596 do ispinor=1,nspinor 597 idat_ispinor = nspinor*(idat-1) + ispinor 598 do ia=1,nincat 599 do ilmn=1,nlmn 600 gxfac(1,ilmn,ia,idat_ispinor) = gxfac(1,ilmn,ia,idat_ispinor) - gxfac_(2,ilmn,ia,idat_ispinor) 601 gxfac(2,ilmn,ia,idat_ispinor) = gxfac(2,ilmn,ia,idat_ispinor) + gxfac_(1,ilmn,ia,idat_ispinor) 602 end do ! ilmn 603 end do ! ia 604 end do ! ispinor 605 end do ! idat 606 ABI_FREE(gxfac_) 607 end if ! iphase==2 608 609 !End loop over real/imaginary part of the exp(-iqR) phase 610 end do 611 612 613 !Accumulate gxfac related to overlap (Sij) (PAW) 614 !------------------------------------------- ------------------------ 615 if (paw_opt==3.or.paw_opt==4) then ! Use Sij, overlap contribution 616 617 gxfac_sij(1:cplex,1:nlmn,1:nincat,1:nspinor) = zero 618 619 do idat = 1,ndat 620 621 do ispinor = 1,nspinor 622 623 idat_ispinor = nspinor*(idat-1) + ispinor 624 625 do ia = 1,nincat 626 627 do jlmn = 1,nlmn 628 629 j0lmn=jlmn*(jlmn-1)/2 630 jjlmn=j0lmn+jlmn 631 jlm=indlmn(4,jlmn) 632 sijr=sij(jjlmn) 633 gxj(1:cplex) = gx(1:cplex,jlmn,ia,idat_ispinor) 634 gxfac_sij(1:cplex,jlmn,ia,idat_ispinor) = gxfac_sij(1:cplex,jlmn,ia,idat_ispinor) + sijr*gxj(1:cplex) 635 636 do ilmn = 1,jlmn-1 637 ilm=indlmn(4,ilmn) 638 !if (ilm==jlm) then 639 ijlmn=j0lmn+ilmn 640 sijr=sij(ijlmn) 641 gxi(1:cplex) = gx(1:cplex,ilmn,ia,idat_ispinor) 642 gxfac_sij(1:cplex,jlmn,ia,idat_ispinor) = gxfac_sij(1:cplex,jlmn,ia,idat_ispinor) + sijr*gxi(1:cplex) 643 !end if 644 end do 645 646 if(jlmn<nlmn) then 647 648 do ilmn = jlmn+1,nlmn 649 ilm=indlmn(4,ilmn) 650 !if (ilm==jlm) then 651 i0lmn=ilmn*(ilmn-1)/2 652 ijlmn=i0lmn+jlmn 653 sijr=sij(ijlmn) 654 gxi(1:cplex) = gx(1:cplex,ilmn,ia,idat_ispinor) 655 gxfac_sij(1:cplex,jlmn,ia,idat_ispinor) = gxfac_sij(1:cplex,jlmn,ia,idat_ispinor) + sijr*gxi(1:cplex) 656 !end if 657 end do ! ilmn 658 659 end if ! jlmn 660 661 end do ! jlmn 662 663 end do ! ia 664 665 end do ! ispinor 666 667 end do ! idat 668 669 end if ! paw_opt 670 671 end subroutine opernlc_ylm_allwf_cpu