Line data Source code
1 : !!****m* abinit/m_pawcprj
2 : !! NAME
3 : !! m_pawcprj
4 : !!
5 : !! FUNCTION
6 : !! This module contains functions used to manipulate variables of
7 : !! structured datatype pawcprj_type.
8 : !! pawcprj_type variables are <p_lmn|Cnk> projected quantities,
9 : !! where |p_lmn> are non-local projectors
10 : !! |Cnk> are wave functions
11 : !!
12 : !! COPYRIGHT
13 : !! Copyright (C) 2012-2026 ABINIT group (MT,JWZ)
14 : !! This file is distributed under the terms of the
15 : !! GNU General Public License, see ~abinit/COPYING
16 : !! or http://www.gnu.org/copyleft/gpl.txt .
17 : !!
18 : !! NOTES
19 : !! FOR DEVELOPPERS: in order to preserve the portability of libPAW library,
20 : !! please consult ~abinit/src/??_libpaw/libpaw-coding-rules.txt
21 : !!
22 : !! SOURCE
23 :
24 : #include "libpaw.h"
25 :
26 : module m_pawcprj
27 :
28 : USE_DEFS
29 : USE_MSG_HANDLING
30 : USE_MPI_WRAPPERS
31 : USE_MEMORY_PROFILING
32 :
33 : use m_pawtab, only : pawtab_type
34 :
35 : implicit none
36 :
37 : private
38 : !!***
39 :
40 : !!****t* m_pawcprj/pawcprj_type
41 : !! NAME
42 : !! pawcprj_type
43 : !!
44 : !! FUNCTION
45 : !! This structured datatype contains <p_lmn|Cnk> projected scalars and derivatives
46 : !! where |p_lmn> are non-local projectors for a given atom
47 : !! |Cnk> is a wave function
48 : !! Used only for PAW calculations.
49 : !!
50 : !! SOURCE
51 :
52 : type,public :: pawcprj_type
53 :
54 : !Integer scalars
55 :
56 : integer :: ncpgr=0
57 : ! Number of gradients of cp=<p_lmn|Cnk>
58 :
59 : integer :: nlmn=0
60 : ! Number of (l,m,n) non-local projectors
61 :
62 : !Real (real(dp)) arrays
63 :
64 : real(dp), allocatable :: cp (:,:)
65 : ! cp(2,nlmn)
66 : ! <p_lmn|Cnk> projected scalars for a given atom and wave function
67 :
68 : real(dp), allocatable :: dcp (:,:,:)
69 : ! dcp(2,ncpgr,nlmn)
70 : ! derivatives of <p_lmn|Cnk> projected scalars for a given atom and wave function
71 :
72 : end type pawcprj_type
73 :
74 : !public procedures.
75 : public :: pawcprj_alloc ! Allocation
76 : public :: pawcprj_free ! Deallocation
77 : public :: pawcprj_set_zero ! Set to zero all arrays in a cprj datastructure
78 : public :: pawcprj_copy ! Copy a cprj datastructure into another
79 : public :: pawcprj_axpby ! cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:)
80 : public :: pawcprj_zaxpby ! cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:), alpha and beta are COMPLEX scalars
81 : public :: pawcprj_projbd ! cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:), alpha and beta are COMPLEX scalars
82 : public :: pawcprj_conjg ! cprj(:,:) <- conjugate(cprj(:,:))
83 : public :: pawcprj_symkn ! construct cprj from that at a symmetry related k point
84 : public :: pawcprj_lincom ! Compute a LINear COMbination of cprj datastructure:
85 : public :: pawcprj_output ! Output a cprj. Useful for debugging.
86 : public :: pawcprj_get ! Read the cprj for a given k-point from memory or from a temporary file
87 : public :: pawcprj_put ! Write the cprj for a given set of (n,k) into memory or into a temporary file
88 : public :: pawcprj_reorder ! Change the order of a cprj datastructure
89 : public :: pawcprj_mpi_allgather ! Perform MPI_ALLGATHER on a pawcprj_type inside a MPI communicator.
90 : public :: pawcprj_bcast ! Broadcast a pawcprj_type from master to all nodes inside a MPI communicator.
91 : public :: pawcprj_transpose ! Transpose a cprj datastructure FOR A GIVEN (K,SPIN)
92 : public :: pawcprj_gather_spin ! Collect spin distributed cprjs.
93 : public :: pawcprj_mpi_exch ! Exchange a pawcprj_type between two processors inside a MPI communicator.
94 : public :: pawcprj_mpi_send ! Send a pawcprj_type inside a MPI communicator.
95 : public :: pawcprj_mpi_recv ! Receive a pawcprj_type inside a MPI communicator.
96 : public :: pawcprj_mpi_sum ! Perform MPI_SUM on a pawcprj_type inside a MPI communicator.
97 : public :: pawcprj_getdim ! Returns the number of lmn components in the <p_{lmn}^i|\psi> for the i-th atom.
98 : public :: paw_overlap ! Compute the onsite contribution to the overlap between two states.
99 : public :: pawcprj_pack ! Copy data from a cprj to a simple real buffer
100 : public :: pawcprj_unpack ! Copy data from a simple real buffer to a cprj
101 : !!***
102 :
103 : CONTAINS
104 :
105 : !===========================================================
106 : !!***
107 :
108 : !!****f* m_pawcprj/pawcprj_alloc
109 : !! NAME
110 : !! pawcprj_alloc
111 : !!
112 : !! FUNCTION
113 : !! Allocation of a cprj datastructure
114 : !!
115 : !! INPUTS
116 : !! ncpgr=number of gradients to be allocated
117 : !! nlmn(:)=sizes of cprj%cp
118 : !!
119 : !! SIDE EFFECTS
120 : !! cprj(:,:) <type(pawcprj_type)>= cprj datastructure
121 : !!
122 : !! SOURCE
123 :
124 4135848 : subroutine pawcprj_alloc(cprj,ncpgr,nlmn)
125 :
126 : !Arguments ------------------------------------
127 : !scalars
128 : integer,intent(in) :: ncpgr
129 : !arrays
130 : integer,intent(in) :: nlmn(:)
131 : type(pawcprj_type),intent(inout) :: cprj(:,:)
132 :
133 : !Local variables-------------------------------
134 : !scalars
135 : integer :: ii,jj,n1dim,n2dim,nn
136 : character(len=500) :: msg
137 :
138 : ! *************************************************************************
139 :
140 4135848 : n1dim=size(cprj,dim=1);n2dim=size(cprj,dim=2);nn=size(nlmn,dim=1)
141 4135848 : if (nn/=n1dim) then
142 0 : write(msg,*) 'wrong sizes (pawcprj_alloc)! :',nn,n1dim
143 0 : LIBPAW_ERROR(msg)
144 : end if
145 :
146 12693911 : do jj=1,n2dim
147 30658472 : do ii=1,n1dim
148 17964561 : if (allocated(cprj(ii,jj)%cp)) then
149 0 : LIBPAW_DEALLOCATE(cprj(ii,jj)%cp)
150 : end if
151 17964561 : if (allocated(cprj(ii,jj)%dcp)) then
152 0 : LIBPAW_DEALLOCATE(cprj(ii,jj)%dcp)
153 : end if
154 17964561 : nn=nlmn(ii)
155 17964561 : cprj(ii,jj)%nlmn=nn
156 53893683 : LIBPAW_ALLOCATE(cprj(ii,jj)%cp,(2,nn))
157 520304496 : cprj(ii,jj)%cp=zero
158 17964561 : cprj(ii,jj)%ncpgr=ncpgr
159 26522624 : if (ncpgr>0) then
160 8669316 : LIBPAW_ALLOCATE(cprj(ii,jj)%dcp,(2,ncpgr,nn))
161 110329299 : cprj(ii,jj)%dcp=zero
162 : end if
163 : end do
164 : end do
165 :
166 4135848 : end subroutine pawcprj_alloc
167 : !!***
168 :
169 : !----------------------------------------------------------------------
170 :
171 : !!****f* m_pawcprj/pawcprj_free
172 : !! NAME
173 : !! pawcprj_free
174 : !!
175 : !! FUNCTION
176 : !! Deallocation of a cprj datastructure
177 : !!
178 : !! SIDE EFFECTS
179 : !! cprj(:,:) <type(pawcprj_type)>= cprj datastructure
180 : !!
181 : !! SOURCE
182 :
183 4144011 : subroutine pawcprj_free(cprj)
184 :
185 : !Arguments ------------------------------------
186 : !scalars
187 : !arrays
188 : type(pawcprj_type),intent(inout) :: cprj(:,:)
189 :
190 : !Local variables-------------------------------
191 : !scalars
192 : integer :: ii,jj,n1dim,n2dim
193 :
194 : ! *************************************************************************
195 :
196 4144011 : n1dim=size(cprj,dim=1);n2dim=size(cprj,dim=2)
197 :
198 12705242 : do jj=1,n2dim
199 30675851 : do ii=1,n1dim
200 17970609 : if (allocated(cprj(ii,jj)%cp)) then
201 17964561 : LIBPAW_DEALLOCATE(cprj(ii,jj)%cp)
202 : end if
203 26531840 : if (allocated(cprj(ii,jj)%dcp)) then
204 2167329 : LIBPAW_DEALLOCATE(cprj(ii,jj)%dcp)
205 : end if
206 : end do
207 : end do
208 :
209 4144011 : end subroutine pawcprj_free
210 : !!***
211 :
212 : !----------------------------------------------------------------------
213 :
214 : !!****f* m_pawcprj/pawcprj_set_zero
215 : !! NAME
216 : !! pawcprj_set_zero
217 : !!
218 : !! FUNCTION
219 : !! Set to zero all arrays in a cprj datastructure
220 : !!
221 : !! SIDE EFFECTS
222 : !! cprj(:,:) <type(pawcprj_type)>= cprj datastructure
223 : !!
224 : !! SOURCE
225 :
226 621131 : subroutine pawcprj_set_zero(cprj)
227 :
228 : !Arguments ------------------------------------
229 : !scalars
230 : !arrays
231 : type(pawcprj_type),intent(inout) :: cprj(:,:)
232 :
233 : !Local variables-------------------------------
234 : !scalars
235 : integer :: ii,jj,n1dim,n2dim
236 :
237 : ! *************************************************************************
238 :
239 621131 : n1dim=size(cprj,dim=1);n2dim=size(cprj,dim=2)
240 :
241 1529999 : do jj=1,n2dim
242 3435917 : do ii=1,n1dim
243 51753192 : if (cprj(ii,jj)%nlmn>0) cprj(ii,jj)%cp(:,:)=zero
244 24000494 : if (cprj(ii,jj)%ncpgr>0) cprj(ii,jj)%dcp(:,:,:)=zero
245 : end do
246 : end do
247 :
248 621131 : end subroutine pawcprj_set_zero
249 : !!***
250 :
251 : !----------------------------------------------------------------------
252 :
253 : !!****f* m_pawcprj/pawcprj_copy
254 : !! NAME
255 : !! pawcprj_copy
256 : !!
257 : !! FUNCTION
258 : !! Copy a cprj datastructure into another
259 : !!
260 : !! INPUTS
261 : !! icpgr= (optional argument) if present, only component icpgr of
262 : !! input cprj gradient is copied into output cprj
263 : !! Not used if cprj(:,:)%ncpgr<icpgr
264 : !! -1 only copy cp
265 : !! cprj_in(:,:) <type(pawcprj_type)>= input cprj datastructure
266 : !!
267 : !! OUTPUT
268 : !! cprj_out(:,:) <type(pawcprj_type)>= output cprj datastructure
269 : !!
270 : !! NOTES
271 : !! MG: What about an option to report a pointer to cprj_in?
272 : !!
273 : !! SOURCE
274 :
275 5550184 : subroutine pawcprj_copy(cprj_in,cprj_out,&
276 : & icpgr) ! optional argument
277 :
278 : !Arguments ------------------------------------
279 : !scalars
280 : integer,intent(in),optional :: icpgr
281 : !arrays
282 : type(pawcprj_type),intent(in) :: cprj_in(:,:)
283 : type(pawcprj_type),intent(inout) :: cprj_out(:,:)
284 :
285 : !Local variables-------------------------------
286 : !scalars
287 : integer :: ii,jj,kk,n1dim_in,n1dim_out,n2dim_in,n2dim_out,ncpgr_in,ncpgr_out,nlmn
288 : logical :: has_icpgr,copy_dcp
289 : character(len=500) :: msg
290 :
291 : ! *************************************************************************
292 :
293 5550184 : n1dim_in=size(cprj_in,dim=1); n1dim_out=size(cprj_out,dim=1)
294 5550184 : n2dim_in=size(cprj_in,dim=2); n2dim_out=size(cprj_out,dim=2)
295 5550184 : ncpgr_in=cprj_in(1,1)%ncpgr; ncpgr_out=cprj_out(1,1)%ncpgr
296 :
297 5550184 : if (n1dim_in/=n1dim_out) then
298 0 : write(msg,'(a,2(1x,i0))')" Error in pawcprj_copy: n1 wrong sizes ",n1dim_in,n1dim_out
299 0 : LIBPAW_ERROR(msg)
300 : end if
301 5550184 : if (n2dim_in/=n2dim_out) then
302 0 : write(msg,'(a,2(1x,i0))')" Error in pawcprj_copy: n2 wrong sizes ",n2dim_in,n2dim_out
303 0 : LIBPAW_ERROR(msg)
304 : end if
305 5550184 : if (ncpgr_in<ncpgr_out) then
306 0 : write(msg,'(a,2(1x,i0))')" Error in pawcprj_copy: ncpgr wrong sizes ",ncpgr_in,ncpgr_out
307 0 : LIBPAW_ERROR(msg)
308 : end if
309 :
310 : !Check if icgr is present and if dcp have to be copy
311 5550184 : has_icpgr=present(icpgr)
312 5550184 : copy_dcp = .TRUE.
313 5550184 : if(has_icpgr)then
314 309168 : copy_dcp = icpgr>=0
315 : end if
316 :
317 20515720 : do jj=1,n2dim_in
318 50023004 : do ii=1,n1dim_in
319 29507284 : nlmn=cprj_in(ii,jj)%nlmn
320 29507284 : cprj_out(ii,jj)%nlmn =nlmn
321 282354026 : do kk=1,nlmn
322 743150902 : cprj_out(ii,jj)%cp(1:2,kk)=cprj_in(ii,jj)%cp(1:2,kk)
323 : end do
324 : end do
325 : end do
326 :
327 5550184 : if (ncpgr_in>0.and.copy_dcp) then
328 1491805 : if (has_icpgr) has_icpgr=(ncpgr_out>0.and.icpgr>0.or.icpgr<=ncpgr_in)
329 :
330 : if (has_icpgr) then
331 400896 : do jj=1,n2dim_in
332 838464 : do ii=1,n1dim_in
333 437568 : nlmn=cprj_in(ii,jj)%nlmn
334 4325520 : do kk=1,nlmn
335 11500080 : cprj_out(ii,jj)%dcp(1:2,1,kk)=cprj_in(ii,jj)%dcp(1:2,icpgr,kk)
336 : end do
337 : end do
338 : end do
339 : else
340 1291357 : if (ncpgr_out>=ncpgr_in) then
341 4383459 : do jj=1,n2dim_in
342 9879677 : do ii=1,n1dim_in
343 5496218 : nlmn=cprj_in(ii,jj)%nlmn
344 52835827 : do kk=1,nlmn
345 1113446776 : cprj_out(ii,jj)%dcp(1:2,1:ncpgr_in,kk)=cprj_in(ii,jj)%dcp(1:2,1:ncpgr_in,kk)
346 : end do
347 : end do
348 : end do
349 : end if
350 : end if
351 : end if
352 :
353 5550184 : end subroutine pawcprj_copy
354 : !!***
355 :
356 : !----------------------------------------------------------------------
357 :
358 : !!****f* m_pawcprj/pawcprj_axpby
359 : !! NAME
360 : !! pawcprj_axpby
361 : !!
362 : !! FUNCTION
363 : !! Apply AXPBY (blas-like) operation with 2 cprj datastructures:
364 : !! cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:)
365 : !! alpha and beta are REAL scalars
366 : !!
367 : !! INPUTS
368 : !! alpha,beta= alpha,beta REAL factors
369 : !! cprjx(:,:) <type(pawcprj_type)>= input cprjx datastructure
370 : !!
371 : !! SIDE EFFECTS
372 : !! cprjy(:,:) <type(pawcprj_type)>= input/output cprjy datastructure
373 : !!
374 : !! SOURCE
375 :
376 2104994 : subroutine pawcprj_axpby(alpha,beta,cprjx,cprjy)
377 :
378 : !Arguments ------------------------------------
379 : !scalars
380 : real(dp),intent(in) :: alpha,beta
381 : !arrays
382 : type(pawcprj_type),intent(in) :: cprjx(:,:)
383 : type(pawcprj_type),intent(inout) :: cprjy(:,:)
384 :
385 : !Local variables-------------------------------
386 : !scalars
387 : integer :: ii,jj,kk,n1dimx,n1dimy,n2dimx,n2dimy,ncpgrx,ncpgry,nlmn
388 : character(len=500) :: msg
389 :
390 : ! *************************************************************************
391 :
392 2104994 : n1dimy=size(cprjy,dim=1);n2dimy=size(cprjy,dim=2);ncpgry=cprjy(1,1)%ncpgr
393 2104994 : if (abs(alpha)>tol16) then
394 2104942 : n1dimx=size(cprjx,dim=1);n2dimx=size(cprjx,dim=2);ncpgrx=cprjx(1,1)%ncpgr
395 2104942 : msg = ""
396 2104942 : if (n1dimx/=n1dimy) msg = TRIM(msg)//"Error in pawcprj_axpby: n1 wrong sizes !"//ch10
397 2104942 : if (n2dimx/=n2dimy) msg = TRIM(msg)//"Error in pawcprj_axpby: n2 wrong sizes !"//ch10
398 2104942 : if (ncpgrx/=ncpgry) msg = TRIM(msg)//"Error in pawcprj_axpby: ncpgr wrong sizes !"//ch10
399 2104942 : if (LEN_TRIM(msg) > 0) then
400 0 : LIBPAW_ERROR(msg)
401 : end if
402 : else
403 : n1dimx=0;n2dimx=0;ncpgrx=0
404 : end if
405 :
406 2104994 : if (abs(alpha)<=tol16) then
407 104 : do jj=1,n2dimy
408 158 : do ii=1,n1dimy
409 54 : nlmn=cprjy(ii,jj)%nlmn
410 538 : do kk=1,nlmn
411 1350 : cprjy(ii,jj)%cp(1:2,kk)=beta*cprjy(ii,jj)%cp(1:2,kk)
412 : end do
413 : end do
414 : end do
415 52 : if (ncpgry>0) then
416 0 : do jj=1,n2dimy
417 0 : do ii=1,n1dimy
418 0 : nlmn=cprjy(ii,jj)%nlmn
419 0 : do kk=1,nlmn
420 0 : cprjy(ii,jj)%dcp(1:2,1:ncpgry,kk)=beta*cprjy(ii,jj)%dcp(1:2,1:ncpgry,kk)
421 : end do
422 : end do
423 : end do
424 : end if
425 2104942 : else if (abs(beta)<=tol16) then
426 562728 : do jj=1,n2dimx
427 1256172 : do ii=1,n1dimx
428 693444 : nlmn=cprjx(ii,jj)%nlmn
429 693444 : cprjy(ii,jj)%nlmn=nlmn
430 6746760 : do kk=1,nlmn
431 18009300 : cprjy(ii,jj)%cp(1:2,kk)=alpha*cprjx(ii,jj)%cp(1:2,kk)
432 : end do
433 : end do
434 : end do
435 281364 : if (ncpgrx>0) then
436 0 : do jj=1,n2dimx
437 0 : do ii=1,n1dimx
438 0 : nlmn=cprjx(ii,jj)%nlmn
439 0 : do kk=1,nlmn
440 0 : cprjy(ii,jj)%dcp(1:2,1:ncpgrx,kk)=alpha*cprjx(ii,jj)%dcp(1:2,1:ncpgrx,kk)
441 : end do
442 : end do
443 : end do
444 : end if
445 : else ! alpha/=0 and beta/=0
446 4391878 : do jj=1,n2dimx
447 9907365 : do ii=1,n1dimx
448 5515487 : nlmn=cprjx(ii,jj)%nlmn
449 5515487 : cprjy(ii,jj)%nlmn=nlmn
450 55724323 : do kk=1,nlmn
451 : cprjy(ii,jj)%cp(1:2,kk)=alpha*cprjx(ii,jj)%cp(1:2,kk) &
452 148437095 : & +beta *cprjy(ii,jj)%cp(1:2,kk)
453 : end do
454 : end do
455 : end do
456 1823578 : if (ncpgrx>0) then
457 0 : do jj=1,n2dimx
458 0 : do ii=1,n1dimx
459 0 : nlmn=cprjx(ii,jj)%nlmn
460 0 : do kk=1,nlmn
461 : cprjy(ii,jj)%dcp(1:2,1:ncpgrx,kk)=alpha*cprjx(ii,jj)%dcp(1:2,1:ncpgrx,kk) &
462 0 : & +beta *cprjy(ii,jj)%dcp(1:2,1:ncpgrx,kk)
463 : end do
464 : end do
465 : end do
466 : end if
467 : end if
468 :
469 2104994 : end subroutine pawcprj_axpby
470 : !!***
471 :
472 : !----------------------------------------------------------------------
473 :
474 : !!****f* m_pawcprj/pawcprj_zaxpby
475 : !! NAME
476 : !! pawcprj_zaxpby
477 : !!
478 : !! FUNCTION
479 : !! Apply ZAXPBY (blas-like) operation with 2 cprj datastructures:
480 : !! cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:)
481 : !! alpha and beta are COMPLEX scalars
482 : !!
483 : !! INPUTS
484 : !! alpha(2),beta(2)= alpha,beta COMPLEX factors
485 : !! cprjx(:,:) <type(pawcprj_type)>= input cprjx datastructure
486 : !!
487 : !! SIDE EFFECTS
488 : !! cprjy(:,:) <type(pawcprj_type)>= input/output cprjy datastructure
489 : !!
490 : !! SOURCE
491 :
492 1107685 : subroutine pawcprj_zaxpby(alpha,beta,cprjx,cprjy)
493 :
494 : !Arguments ------------------------------------
495 : !scalars
496 : real(dp),intent(in) :: alpha(2),beta(2)
497 : !arrays
498 : type(pawcprj_type),intent(in) :: cprjx(:,:)
499 : type(pawcprj_type),intent(inout) :: cprjy(:,:)
500 :
501 : !Local variables-------------------------------
502 : !scalars
503 : integer :: ii,jj,kk,ll,n1dimx,n1dimy,n2dimx,n2dimy,ncpgrx,ncpgry,nlmn
504 : real(dp) :: cp1,cp2,norma,normb
505 : character(len=500) :: msg
506 :
507 : ! *************************************************************************
508 :
509 1107685 : norma=alpha(1)**2+alpha(2)**2
510 1107685 : normb=beta(1) **2+beta(2) **2
511 1107685 : n1dimy=size(cprjy,dim=1);n2dimy=size(cprjy,dim=2);ncpgry=cprjy(1,1)%ncpgr
512 1107685 : if (norma>tol16*tol16) then
513 1107606 : n1dimx=size(cprjx,dim=1);n2dimx=size(cprjx,dim=2);ncpgrx=cprjx(1,1)%ncpgr
514 1107606 : msg = ""
515 1107606 : if (n1dimx/=n1dimy) msg = TRIM(msg)//"Error in pawcprj_zaxpby: n1 wrong sizes !"//ch10
516 1107606 : if (n2dimx/=n2dimy) msg = TRIM(msg)//"Error in pawcprj_zaxpby: n2 wrong sizes !"//ch10
517 1107606 : if (ncpgrx/=ncpgry) msg = TRIM(msg)//"Error in pawcprj_zaxpby: ncpgr wrong sizes !"//ch10
518 1107606 : if (LEN_TRIM(msg) > 0) then
519 0 : LIBPAW_ERROR(msg)
520 : end if
521 : end if
522 :
523 1107685 : if (norma<=tol16*tol16) then
524 158 : do jj=1,n2dimy
525 438 : do ii=1,n1dimy
526 280 : nlmn=cprjy(ii,jj)%nlmn
527 2719 : do kk=1,nlmn
528 2360 : cp1=beta(1)*cprjy(ii,jj)%cp(1,kk)-beta(2)*cprjy(ii,jj)%cp(2,kk)
529 2360 : cp2=beta(1)*cprjy(ii,jj)%cp(2,kk)+beta(2)*cprjy(ii,jj)%cp(1,kk)
530 2360 : cprjy(ii,jj)%cp(1,kk)=cp1
531 2640 : cprjy(ii,jj)%cp(2,kk)=cp2
532 : end do
533 : end do
534 : end do
535 79 : if (ncpgry>0) then
536 0 : do jj=1,n2dimy
537 0 : do ii=1,n1dimy
538 0 : nlmn=cprjy(ii,jj)%nlmn
539 0 : do kk=1,nlmn
540 0 : do ll=1,ncpgry
541 0 : cp1=beta(1)*cprjy(ii,jj)%dcp(1,ll,kk)-beta(2)*cprjy(ii,jj)%dcp(2,ll,kk)
542 0 : cp2=beta(1)*cprjy(ii,jj)%dcp(2,ll,kk)+beta(2)*cprjy(ii,jj)%dcp(1,ll,kk)
543 0 : cprjy(ii,jj)%dcp(1,ll,kk)=cp1
544 0 : cprjy(ii,jj)%dcp(2,ll,kk)=cp2
545 : end do
546 : end do
547 : end do
548 : end do
549 : end if
550 1107606 : else if (normb<=tol16*tol16) then
551 0 : do jj=1,n2dimx
552 0 : do ii=1,n1dimx
553 0 : nlmn=cprjx(ii,jj)%nlmn
554 0 : cprjy(ii,jj)%nlmn=nlmn
555 0 : do kk=1,nlmn
556 0 : cprjy(ii,jj)%cp(1,kk)=alpha(1)*cprjx(ii,jj)%cp(1,kk)-alpha(2)*cprjx(ii,jj)%cp(2,kk)
557 0 : cprjy(ii,jj)%cp(2,kk)=alpha(1)*cprjx(ii,jj)%cp(2,kk)+alpha(2)*cprjx(ii,jj)%cp(1,kk)
558 : end do
559 : end do
560 : end do
561 0 : if (ncpgrx>0) then
562 0 : do jj=1,n2dimx
563 0 : do ii=1,n1dimx
564 0 : nlmn=cprjx(ii,jj)%nlmn
565 0 : do kk=1,nlmn
566 : cprjy(ii,jj)%dcp(1,1:ncpgrx,kk)=alpha(1)*cprjx(ii,jj)%dcp(1,1:ncpgrx,kk) &
567 0 : & -alpha(2)*cprjx(ii,jj)%dcp(2,1:ncpgrx,kk)
568 : cprjy(ii,jj)%dcp(2,1:ncpgrx,kk)=alpha(1)*cprjx(ii,jj)%dcp(2,1:ncpgrx,kk) &
569 0 : & +alpha(2)*cprjx(ii,jj)%dcp(1,1:ncpgrx,kk)
570 : end do
571 : end do
572 : end do
573 : end if
574 : ! else if (abs(beta(1)-one)<tol16.and.abs(beta(2))<tol16) then
575 : ! do jj=1,n2dimx
576 : ! do ii=1,n1dimx
577 : ! nlmn=cprjx(ii,jj)%nlmn
578 : ! cprjy(ii,jj)%nlmn =nlmn
579 : ! do kk=1,nlmn
580 : ! cp1=cprjy(ii,jj)%cp(1,kk)
581 : ! cp2=cprjy(ii,jj)%cp(2,kk)
582 : ! cp1=cp1+alpha(1)*cprjx(ii,jj)%cp(1,kk)-alpha(2)*cprjx(ii,jj)%cp(2,kk)
583 : ! cp2=cp2+alpha(1)*cprjx(ii,jj)%cp(2,kk)+alpha(2)*cprjx(ii,jj)%cp(1,kk)
584 : ! cprjy(ii,jj)%cp(1,kk)=cp1
585 : ! cprjy(ii,jj)%cp(2,kk)=cp2
586 : ! end do
587 : ! end do
588 : ! end do
589 : ! if (ncpgrx>0) then
590 : ! do jj=1,n2dimx
591 : ! do ii=1,n1dimx
592 : ! nlmn=cprjx(ii,jj)%nlmn
593 : ! do kk=1,nlmn
594 : ! do ll=1,ncpgrx
595 : ! cp1=cprjy(ii,jj)%dcp(1,ll,kk)
596 : ! cp2=cprjy(ii,jj)%dcp(2,ll,kk)
597 : ! cp1=cp1+alpha(1)*cprjx(ii,jj)%dcp(1,ll,kk)-alpha(2)*cprjx(ii,jj)%dcp(2,ll,kk)
598 : ! cp2=cp2+alpha(1)*cprjx(ii,jj)%dcp(2,ll,kk)+alpha(2)*cprjx(ii,jj)%dcp(1,ll,kk)
599 : ! cprjy(ii,jj)%dcp(1,ll,kk)=cp1
600 : ! cprjy(ii,jj)%dcp(2,ll,kk)=cp2
601 : ! end do
602 : ! end do
603 : ! end do
604 : ! end do
605 : ! end if
606 : else
607 2257512 : do jj=1,n2dimx
608 4605817 : do ii=1,n1dimx
609 2348305 : nlmn=cprjx(ii,jj)%nlmn
610 2348305 : cprjy(ii,jj)%nlmn =nlmn
611 24723531 : do kk=1,nlmn
612 : cp1=alpha(1)*cprjx(ii,jj)%cp(1,kk)-alpha(2)*cprjx(ii,jj)%cp(2,kk) &
613 21225320 : & +beta(1) *cprjy(ii,jj)%cp(1,kk)-beta(2) *cprjy(ii,jj)%cp(2,kk)
614 : cp2=alpha(1)*cprjx(ii,jj)%cp(2,kk)+alpha(2)*cprjx(ii,jj)%cp(1,kk) &
615 21225320 : & +beta(1) *cprjy(ii,jj)%cp(2,kk)+beta(2) *cprjy(ii,jj)%cp(1,kk)
616 : ! cp1=beta(1) *cprjy(ii,jj)%cp(1,kk)-beta(2) *cprjy(ii,jj)%cp(2,kk)
617 : ! cp1=cp1+alpha(1)*cprjx(ii,jj)%cp(1,kk)-alpha(2)*cprjx(ii,jj)%cp(2,kk)
618 : ! cp2=beta(1) *cprjy(ii,jj)%cp(2,kk)+beta(2) *cprjy(ii,jj)%cp(1,kk)
619 : ! cp2=cp2+alpha(1)*cprjx(ii,jj)%cp(2,kk)+alpha(2)*cprjx(ii,jj)%cp(1,kk)
620 21225320 : cprjy(ii,jj)%cp(1,kk)=cp1
621 23573625 : cprjy(ii,jj)%cp(2,kk)=cp2
622 : end do
623 : end do
624 : end do
625 1107606 : if (ncpgrx>0) then
626 0 : do jj=1,n2dimx
627 0 : do ii=1,n1dimx
628 0 : nlmn=cprjx(ii,jj)%nlmn
629 0 : do kk=1,nlmn
630 0 : do ll=1,ncpgrx
631 : cp1=alpha(1)*cprjx(ii,jj)%dcp(1,ll,kk)-alpha(2)*cprjx(ii,jj)%dcp(2,ll,kk) &
632 0 : & +beta(1) *cprjy(ii,jj)%dcp(1,ll,kk)-beta(2) *cprjy(ii,jj)%dcp(2,ll,kk)
633 : cp2=alpha(1)*cprjx(ii,jj)%dcp(2,ll,kk)+alpha(2)*cprjx(ii,jj)%dcp(1,ll,kk) &
634 0 : & +beta(1) *cprjy(ii,jj)%dcp(2,ll,kk)+beta(2) *cprjy(ii,jj)%dcp(1,ll,kk)
635 : ! cp1=beta(1) *cprjy(ii,jj)%dcp(1,ll,kk)-beta(2) *cprjy(ii,jj)%dcp(2,ll,kk)
636 : ! cp1=cp1+alpha(1)*cprjx(ii,jj)%dcp(1,ll,kk)-alpha(2)*cprjx(ii,jj)%dcp(2,ll,kk)
637 : ! cp2=beta(1) *cprjy(ii,jj)%dcp(2,ll,kk)+beta(2) *cprjy(ii,jj)%dcp(1,ll,kk)
638 : ! cp2=cp2+alpha(1)*cprjx(ii,jj)%dcp(2,ll,kk)+alpha(2)*cprjx(ii,jj)%dcp(1,ll,kk)
639 0 : cprjy(ii,jj)%dcp(1,ll,kk)=cp1
640 0 : cprjy(ii,jj)%dcp(2,ll,kk)=cp2
641 : end do
642 : end do
643 : end do
644 : end do
645 : end if
646 : end if
647 :
648 1107685 : end subroutine pawcprj_zaxpby
649 : !!***
650 :
651 : !----------------------------------------------------------------------
652 :
653 : !!****f* m_pawcprj/pawcprj_projbd
654 : !! NAME
655 : !! pawcprj_projbd
656 : !!
657 : !! FUNCTION
658 : !! Apply ZAXPBY (blas-like) operation with 2 cprj datastructures:
659 : !! cprjy(:,:) <- alpha.cprjx(:,:)+beta.cprjy(:,:)
660 : !! alpha and beta are COMPLEX scalars
661 : !!
662 : !! INPUTS
663 : !! alpha(2),beta(2)= alpha,beta COMPLEX factors
664 : !! cprjx(:,:) <type(pawcprj_type)>= input cprjx datastructure
665 : !!
666 : !! SIDE EFFECTS
667 : !! cprjy(:,:) <type(pawcprj_type)>= input/output cprjy datastructure
668 : !!
669 : !! SOURCE
670 :
671 0 : subroutine pawcprj_projbd(alpha,cprjx,cprjy)
672 :
673 : !Arguments ------------------------------------
674 : !scalars
675 : real(dp),intent(in) :: alpha(:,:)
676 : !arrays
677 : type(pawcprj_type),intent(in) :: cprjx(:,:)
678 : type(pawcprj_type),intent(inout) :: cprjy(:,:)
679 :
680 : !Local variables-------------------------------
681 : !scalars
682 : integer :: ia,ii,jj,kk,ll,n1dima,n1dimx,n1dimy,n2dimx,n2dimy,n2dima,ncpgrx,ncpgry,nlmn
683 : real(dp) :: cp1,cp2,norma
684 : character(len=500) :: msg
685 :
686 : ! *************************************************************************
687 :
688 0 : n1dimy=size(cprjy,dim=1);n2dimy=size(cprjy,dim=2);ncpgry=cprjy(1,1)%ncpgr
689 0 : n1dimx=size(cprjx,dim=1);n2dimx=size(cprjx,dim=2);ncpgrx=cprjx(1,1)%ncpgr
690 0 : n1dima=size(alpha,dim=1);n2dima=size(alpha,dim=2)
691 0 : msg = ""
692 0 : if (n1dima/=2) msg = TRIM(msg)//"Error in pawcprj_projbd: alpha n1 wrong sizes !"//ch10
693 0 : if (n1dimx/=n1dimy) msg = TRIM(msg)//"Error in pawcprj_projbd: n1 wrong sizes !"//ch10
694 0 : if (n2dimx/=n2dimy*n2dima) msg = TRIM(msg)//"Error in pawcprj_projbd: n2 wrong sizes !"//ch10
695 0 : if (ncpgrx/=ncpgry) msg = TRIM(msg)//"Error in pawcprj_projbd: ncpgr wrong sizes !"//ch10
696 0 : if (LEN_TRIM(msg) > 0) then
697 0 : LIBPAW_ERROR(msg)
698 : end if
699 :
700 0 : do ia=1,n2dima
701 0 : norma=alpha(1,ia)**2+alpha(2,ia)**2
702 0 : if (norma>tol16*tol16) then
703 0 : do jj=1,n2dimy
704 0 : do ii=1,n1dimx
705 0 : nlmn=cprjy(ii,jj)%nlmn
706 : cprjy(ii,jj)%nlmn =nlmn
707 0 : do kk=1,nlmn
708 0 : cp1=alpha(1,ia)*cprjx(ii,jj+(ia-1)*n2dimy)%cp(1,kk)-alpha(2,ia)*cprjx(ii,jj+(ia-1)*n2dimy)%cp(2,kk)
709 0 : cp2=alpha(1,ia)*cprjx(ii,jj+(ia-1)*n2dimy)%cp(2,kk)+alpha(2,ia)*cprjx(ii,jj+(ia-1)*n2dimy)%cp(1,kk)
710 0 : cprjy(ii,jj)%cp(1,kk)=cprjy(ii,jj)%cp(1,kk)+cp1
711 0 : cprjy(ii,jj)%cp(2,kk)=cprjy(ii,jj)%cp(2,kk)+cp2
712 : end do
713 : end do
714 : end do
715 0 : if (ncpgrx>0) then
716 0 : do jj=1,n2dimy
717 0 : do ii=1,n1dimx
718 0 : nlmn=cprjy(ii,jj)%nlmn
719 0 : do kk=1,nlmn
720 0 : do ll=1,ncpgrx
721 0 : cp1=alpha(1,ia)*cprjx(ii,jj+(ia-1)*n2dimy)%dcp(1,ll,kk)-alpha(2,ia)*cprjx(ii,jj+(ia-1)*n2dimy)%dcp(2,ll,kk)
722 0 : cp2=alpha(1,ia)*cprjx(ii,jj+(ia-1)*n2dimy)%dcp(2,ll,kk)+alpha(2,ia)*cprjx(ii,jj+(ia-1)*n2dimy)%dcp(1,ll,kk)
723 0 : cprjy(ii,jj)%dcp(1,ll,kk)=cprjy(ii,jj)%dcp(1,ll,kk)+cp1
724 0 : cprjy(ii,jj)%dcp(2,ll,kk)=cprjy(ii,jj)%dcp(2,ll,kk)+cp2
725 : end do
726 : end do
727 : end do
728 : end do
729 : end if
730 : end if
731 : end do
732 :
733 0 : end subroutine pawcprj_projbd
734 : !!***
735 :
736 : !----------------------------------------------------------------------
737 :
738 : !!****f* m_pawcprj/pawcprj_symkn
739 : !! NAME
740 : !! pawcprj_symkn
741 : !!
742 : !! FUNCTION
743 : !! compute cprj for a given band and k point based on cprj at a symmetry-related
744 : !! k point.
745 : !!
746 : !! INPUTS
747 : !! cprj_ikn (pawcprj_type) :: cprj for a single band and k point, typically a k point in the IBZ
748 : !! cprj_sym(4,nsym,natom) :: 1:3 shift, and 4 final atom, of symmetry isym operating on iatom
749 : !! (S^{-1}(R - t) = r0 + L, see symatm.F90
750 : !! dimlmn(natom) :: ln dimension of each atom
751 : !! iband :: number of bands to treat, use -1 to treat all nband bands
752 : !! indlmn(6,lmnmax,ntypat) :: n,l,m dimensions for each atom type (see psps type)
753 : !! isym :: symmetry element used in current application
754 : !! itim :: 1 if time reversal also used, 0 else
755 : !! kpt(3) :: kpt vector used
756 : !! lmax :: max l value
757 : !! lmnmax :: max lmn value
758 : !! mband :: maximum number of bands
759 : !! natom :: number of atoms in cell
760 : !! nband :: number of bands in cprj_ikn
761 : !! nspinor :: number of spinors
762 : !! nsym :: total number of symmetry elements
763 : !! ntypat :: number of types of atoms
764 : !! typat(natom) :: type of each atom
765 : !! zarot(2*lmax+1,2*lmax+1,lmax+1,nsym) :: elements of rotation matrix for angular momentum states
766 : !! and symmetry operations. See m_paw_sphharm/setsym_ylm.
767 : !!
768 : !! OUTPUT
769 : !! cprj_fkn (pawcprj_type) :: cprj for a single band and k point where the k point is related to
770 : !! the input k point by a symmetry operation
771 : !!
772 : !! SIDE EFFECTS
773 : !!
774 : !! NOTES
775 : !! This routine is based on M. Giantomassi's doctoral dissertation, formula 7.77. It is not clear
776 : !! whether it is implemented correctly for nonsymmorphic symmetries.
777 : !!
778 : !! SOURCE
779 :
780 1305428 : subroutine pawcprj_symkn(cprj_fkn,cprj_ikn,cprj_sym,dimlmn,iband,indlmn,&
781 : & isym,itim,kpt,lmax,lmnmax,mband,natom,nband,nspinor,nsym,ntypat,&
782 1305428 : & typat,zarot)
783 :
784 : !Arguments---------------------------
785 : !scalars
786 : integer,intent(in) :: iband,isym,itim,lmax,lmnmax,mband
787 : integer,intent(in) :: natom,nband,nspinor,nsym,ntypat
788 :
789 : !arrays
790 : integer,intent(in) :: cprj_sym(4,nsym,natom),dimlmn(natom)
791 : integer,intent(in) :: indlmn(6,lmnmax,ntypat),typat(natom)
792 : real(dp),intent(in) :: kpt(3)
793 : real(dp),intent(in) :: zarot(2*lmax+1,2*lmax+1,lmax+1,nsym)
794 : type(pawcprj_type),intent(in) :: cprj_ikn(natom,mband*nspinor)
795 : type(pawcprj_type),intent(inout) :: cprj_fkn(natom,mband*nspinor) !vz_i
796 :
797 : !Local variables---------------------------
798 : !scalars
799 : integer :: iatm,iatom, ibct, ibnd, ibsp, ibst, icpgr, iin, il, il0, im
800 : integer :: ilmn, iln, iln0, ilpm, indexi, ispinor, itypat, jatm,jatom, mm, nlmn
801 : real(dp) :: kdotL, phr, phi
802 : !arrays
803 : real(dp) :: rl(3), t1(2), t2(2)
804 :
805 : ! *************************************************************************
806 :
807 1305428 : if (iband == -1) then
808 1302958 : ibst = 1
809 1302958 : ibnd = nband
810 : else
811 : ibst = iband
812 : ibnd = iband
813 : end if
814 :
815 3916284 : do iatom = 1, natom
816 2610856 : iatm=iatom
817 2610856 : itypat = typat(iatom)
818 2610856 : nlmn = dimlmn(iatm)
819 2610856 : jatom = cprj_sym(4,isym,iatom)
820 2610856 : jatm=jatom
821 10443424 : rl(:) = cprj_sym(1:3,isym,iatom)
822 10443424 : kdotL = dot_product(rl,kpt)
823 2610856 : phr = cos(two_pi*kdotL)
824 2610856 : phi = sin(two_pi*kdotL)
825 :
826 2610856 : il0 = -1; iln0 = -1; indexi = 1
827 24799952 : do ilmn = 1, nlmn
828 :
829 20883668 : il = indlmn(1,ilmn,itypat)
830 20883668 : im = indlmn(2,ilmn,itypat)
831 20883668 : iin = indlmn(3,ilmn,itypat)
832 20883668 : iln = indlmn(5,ilmn,itypat)
833 20883668 : ilpm = 1 + il + im
834 20883668 : if (iln /= iln0) indexi = indexi + 2*il0 + 1
835 :
836 104309320 : do ibct = ibst, ibnd
837 :
838 187734972 : do ispinor = 1, nspinor
839 :
840 83425652 : ibsp = nspinor*(ibct-1) + ispinor
841 :
842 83425652 : t1(:) = zero
843 291988192 : do mm = 1, 2*il+1
844 208562540 : t1(1) = t1(1) + zarot(mm,ilpm,il+1,isym)*cprj_ikn(jatm,ibsp)%cp(1,indexi+mm)
845 291988192 : t1(2) = t1(2) + zarot(mm,ilpm,il+1,isym)*cprj_ikn(jatm,ibsp)%cp(2,indexi+mm)
846 : end do
847 83425652 : t2(1) = t1(1)*phr - t1(2)*phi
848 83425652 : t2(2) = t1(2)*phr + t1(1)*phi
849 :
850 83425652 : if (itim == 1) t2(2) = -t2(2)
851 :
852 83425652 : cprj_fkn(iatm,ibsp)%cp(1,ilmn) = t2(1)
853 83425652 : cprj_fkn(iatm,ibsp)%cp(2,ilmn) = t2(2)
854 :
855 : ! do same transformations for gradients of cprj_ikn
856 : ! note that ncpgr = 0 if no gradients present so this loop will not be executed
857 : ! in this case
858 :
859 339022084 : do icpgr = 1, cprj_ikn(jatom,ibsp)%ncpgr
860 : t1(:) = zero
861 :
862 602592960 : do mm = 1, 2*il+1
863 430422180 : t1(1) = t1(1) + zarot(mm,ilpm,il+1,isym)*cprj_ikn(jatm,ibsp)%dcp(1,icpgr,indexi+mm)
864 602592960 : t1(2) = t1(2) + zarot(mm,ilpm,il+1,isym)*cprj_ikn(jatm,ibsp)%dcp(2,icpgr,indexi+mm)
865 : end do
866 :
867 172170780 : t2(1) = t1(1)*phr - t1(2)*phi
868 172170780 : t2(2) = t1(2)*phr + t1(1)*phi
869 :
870 172170780 : if (itim == 1) t2(2) = -t2(2)
871 :
872 172170780 : cprj_fkn(iatm,ibsp)%dcp(1,icpgr,ilmn) = t2(1)
873 255596432 : cprj_fkn(iatm,ibsp)%dcp(2,icpgr,ilmn) = t2(2)
874 :
875 : end do ! end loop over ncpgr
876 :
877 : end do ! end loop over nspinor
878 :
879 : end do ! end loop over bands
880 :
881 23494524 : il0 = il; iln0 = iln
882 : end do ! end loop over ilmn
883 : end do ! end loop over atoms
884 :
885 1305428 : end subroutine pawcprj_symkn
886 : !!***
887 :
888 : !----------------------------------------------------------------------
889 :
890 : !!****f* m_pawcprj/pawcprj_conjg
891 : !! NAME
892 : !! pawcprj_conjg
893 : !!
894 : !! FUNCTION
895 : !! conjugate a cprj datastructures:
896 : !! cprj(:,:) <- conjugate(cprj(:,:))
897 : !!
898 : !! INPUTS
899 : !!
900 : !! SIDE EFFECTS
901 : !! cprj(:,:) <type(pawcprj_type)>= input/output cprj datastructure
902 : !!
903 : !! SOURCE
904 :
905 0 : subroutine pawcprj_conjg(cprj)
906 :
907 : !Arguments ------------------------------------
908 : !scalars
909 : !arrays
910 : type(pawcprj_type),intent(inout) :: cprj(:,:)
911 :
912 : !Local variables-------------------------------
913 : !scalars
914 : integer :: ii,jj,kk,n1dim,n2dim,ncpgr,nlmn
915 :
916 : ! *************************************************************************
917 :
918 :
919 0 : n1dim=size(cprj,dim=1);n2dim=size(cprj,dim=2);ncpgr=cprj(1,1)%ncpgr
920 :
921 0 : do jj=1,n2dim
922 0 : do ii=1,n1dim
923 0 : nlmn=cprj(ii,jj)%nlmn
924 0 : do kk=1,nlmn
925 0 : cprj(ii,jj)%cp(2,kk)=-cprj(ii,jj)%cp(2,kk)
926 : end do
927 : end do
928 : end do
929 0 : if (ncpgr>0) then
930 0 : do jj=1,n2dim
931 0 : do ii=1,n1dim
932 0 : nlmn=cprj(ii,jj)%nlmn
933 0 : do kk=1,nlmn
934 0 : cprj(ii,jj)%dcp(2,1:ncpgr,kk)=-cprj(ii,jj)%dcp(2,1:ncpgr,kk)
935 : end do
936 : end do
937 : end do
938 : end if
939 :
940 0 : end subroutine pawcprj_conjg
941 : !!***
942 :
943 : !----------------------------------------------------------------------
944 :
945 : !!****f* m_pawcprj/pawcprj_lincom
946 : !! NAME
947 : !! pawcprj_lincom
948 : !!
949 : !! FUNCTION
950 : !! Compute a LINear COMbination of cprj datastructure:
951 : !! cprj_out(:,:) <--- Sum_i [ alpha_i . cprj_i(:,:) ]
952 : !! alpha_i are COMPLEX scalars
953 : !!
954 : !! INPUTS
955 : !! alpha(2,nn)= alpha COMPLEX factors
956 : !! cprj_in(:,:) <type(pawcprj_type)>= input cprj_in datastructure
957 : !! nn= number of cprj involved in the linear combination
958 : !!
959 : !! OUTPUT
960 : !! cprj_out(:,:) <type(pawcprj_type)>= output cprj_out datastructure
961 : !!
962 : !! NOTES
963 : !! cprj_in and cprj_out must be dimensionned as cprj_in(n1,n2*nn) and cprj_in(n1,n2)
964 : !!
965 : !! SOURCE
966 :
967 288348 : subroutine pawcprj_lincom(alpha,cprj_in,cprj_out,nn)
968 :
969 : !Arguments ------------------------------------
970 : !scalars
971 : integer,intent(in) :: nn
972 : real(dp),intent(in) :: alpha(2,nn)
973 : !arrays
974 : type(pawcprj_type),intent(in) :: cprj_in(:,:)
975 : type(pawcprj_type),intent(inout) :: cprj_out(:,:)
976 :
977 : !Local variables-------------------------------
978 : !scalars
979 : integer :: ii,in,jj,jn,kk,ll,n1in,n1out,n2in,n2out,ncpgrin,ncpgrout,nlmn
980 : real(dp) :: cp1,cp2
981 : character(len=500) :: msg
982 :
983 : ! *************************************************************************
984 :
985 288348 : n1in=size(cprj_in,dim=1);n1out=size(cprj_out,dim=1)
986 288348 : n2in=size(cprj_in,dim=2);n2out=size(cprj_out,dim=2)
987 288348 : ncpgrin=cprj_in(1,1)%ncpgr;ncpgrout=cprj_out(1,1)%ncpgr
988 :
989 288348 : msg = ""
990 288348 : if (n1in/=n1out) msg = TRIM(msg)//"Bug in pawcprj_lincom: n1 wrong sizes!"//ch10
991 288348 : if (n2in/=n2out*nn) msg = TRIM(msg)//"Bug in pawcprj_lincom: n2 wrong sizes!"//ch10
992 288348 : if (ncpgrin/=ncpgrout) msg = TRIM(msg)//"Bug in pawcprj_lincom: ncpgr wrong sizes!"//ch10
993 288348 : if (LEN_TRIM(msg) > 0) then
994 0 : LIBPAW_ERROR(msg)
995 : end if
996 :
997 576696 : do jj=1,n2out
998 1324932 : do ii=1,n1out
999 748236 : nlmn=cprj_in(ii,jj)%nlmn
1000 748236 : cprj_out(ii,jj)%nlmn=nlmn
1001 19379100 : cprj_out(ii,jj)%cp(1:2,1:nlmn)=zero
1002 : jn=jj
1003 6315696 : do in=1,nn
1004 50231448 : do kk=1,nlmn
1005 : cp1=cprj_out(ii,jj)%cp(1,kk) &
1006 44952336 : & +alpha(1,in)*cprj_in(ii,jn)%cp(1,kk)-alpha(2,in)*cprj_in(ii,jn)%cp(2,kk)
1007 : cp2=cprj_out(ii,jj)%cp(2,kk) &
1008 44952336 : & +alpha(1,in)*cprj_in(ii,jn)%cp(2,kk)+alpha(2,in)*cprj_in(ii,jn)%cp(1,kk)
1009 44952336 : cprj_out(ii,jj)%cp(1,kk)=cp1
1010 50231448 : cprj_out(ii,jj)%cp(2,kk)=cp2
1011 : end do
1012 6027348 : jn=jn+n2out
1013 : end do
1014 : end do
1015 : end do
1016 :
1017 288348 : if (ncpgrin>0) then
1018 0 : do jj=1,n2out
1019 0 : do ii=1,n1out
1020 0 : nlmn=cprj_in(ii,jj)%nlmn
1021 0 : cprj_out(ii,jj)%dcp(1:2,1:ncpgrin,1:nlmn)=zero
1022 : jn=jj
1023 0 : do in=1,nn
1024 0 : do kk=1,nlmn
1025 0 : do ll=1,ncpgrin
1026 : cp1=cprj_out(ii,jj)%dcp(1,ll,kk) &
1027 : & +alpha(1,in)*cprj_in(ii,jn)%dcp(1,ll,kk) &
1028 0 : & -alpha(2,in)*cprj_in(ii,jn)%dcp(2,ll,kk)
1029 : cp2=cprj_out(ii,jj)%dcp(2,ll,kk) &
1030 : & +alpha(1,in)*cprj_in(ii,jn)%dcp(2,ll,kk) &
1031 0 : +alpha(2,in)*cprj_in(ii,jn)%dcp(1,ll,kk)
1032 0 : cprj_out(ii,jj)%dcp(1,ll,kk)=cp1
1033 0 : cprj_out(ii,jj)%dcp(2,ll,kk)=cp2
1034 : end do
1035 : end do
1036 0 : jn=jn+n2out
1037 : end do
1038 : end do
1039 : end do
1040 : end if
1041 :
1042 288348 : end subroutine pawcprj_lincom
1043 : !!***
1044 :
1045 : !----------------------------------------------------------------------
1046 :
1047 : !!****f* m_pawcprj/pawcprj_output
1048 : !! NAME
1049 : !! pawcprj_output
1050 : !!
1051 : !! FUNCTION
1052 : !! Output a cprj. Useful for debugging.
1053 : !!
1054 : !! INPUTS
1055 : !! cprj(:,:) <type(pawcprj_type)>= cprj datastructure
1056 : !! prtgrads :: optional, 1 to print gradients also
1057 : !!
1058 : !! OUTPUT
1059 : !!
1060 : !! SOURCE
1061 :
1062 0 : subroutine pawcprj_output(cprj,prtgrads)
1063 :
1064 : !Arguments ------------------------------------
1065 : !scalar
1066 : integer,optional :: prtgrads
1067 : !arrays
1068 : type(pawcprj_type),intent(in) :: cprj(:,:)
1069 :
1070 : !Local variables-------------------------------
1071 : !scalar
1072 : integer :: ii,jj,kk,nlmn,n1dim,n2dim
1073 : logical :: gradoutput
1074 :
1075 : ! *************************************************************************
1076 :
1077 0 : n1dim=size(cprj,dim=1)
1078 0 : n2dim=size(cprj,dim=2)
1079 0 : gradoutput = .FALSE.
1080 0 : if(present(prtgrads)) then
1081 0 : gradoutput = (prtgrads .EQ. 1)
1082 : end if
1083 :
1084 0 : write(std_out,'(a)')' pawcprj_output '
1085 :
1086 0 : do jj=1,n2dim
1087 0 : do ii=1,n1dim
1088 0 : write(std_out,'(a,i4,a,i4)')'atom ',ii,' band*k ',jj
1089 0 : nlmn=cprj(ii,jj)%nlmn
1090 0 : do kk=1,nlmn
1091 0 : write(std_out,'(2f12.8)')cprj(ii,jj)%cp(1,kk),cprj(ii,jj)%cp(2,kk)
1092 0 : if(gradoutput) then
1093 0 : write(std_out,'(6f12.8)')cprj(ii,jj)%dcp(1,1,kk),cprj(ii,jj)%dcp(2,1,kk),&
1094 0 : &cprj(ii,jj)%dcp(1,2,kk),cprj(ii,jj)%dcp(2,2,kk),&
1095 0 : &cprj(ii,jj)%dcp(1,3,kk),cprj(ii,jj)%dcp(2,3,kk)
1096 : end if
1097 : end do
1098 : end do
1099 : end do
1100 :
1101 0 : end subroutine pawcprj_output
1102 : !!***
1103 :
1104 : !----------------------------------------------------------------------
1105 :
1106 : !!****f* m_pawcprj/pawcprj_get
1107 : !! NAME
1108 : !! pawcprj_get
1109 : !!
1110 : !! FUNCTION
1111 : !! Read the cprj_k for a given k-point from memory in cprj or from a temporary file
1112 : !!
1113 : !! INPUTS
1114 : !! atind(natom)=index table for atoms (see iorder below)
1115 : !! cprj(dimcp,nspinor*mband*mkmem*nsppol)=input cprj (used if mkmem/=0)
1116 : !! dimcp=first dimension of cprj_k,cprj arrays (1 or natom)
1117 : !! iband1=index of first band in cprj
1118 : !! ibg=shift in cprj array to locate current k-point
1119 : !! [icpgr]= (optional argument) if present, only component icpgr of
1120 : !! input cprj gradient is copied into output cprj
1121 : !! Not used if cprj(:,:)%ncpgr<icpgr (mkmem>0)
1122 : !! or ncpgr(optional)<icpgr (mkmem=0)
1123 : !! ikpt=index of current k-point (only needed for the parallel distribution)
1124 : !! iorder=0 if cprj ordering does not change during reading
1125 : !! 1 if cprj ordering changes during reading, depending on content of atind array:
1126 : !! - if atind=atindx (type-sorted=>unsorted)
1127 : !! - if atind=atindx1 (unsorted=>type-sorted)
1128 : !! isppol=index of current spin component
1129 : !! mband=maximum number of bands
1130 : !! mkmem=number of k points which can fit in memory; set to 0 if use disk
1131 : !! [mpi_comm]=(optional argument) MPI communicator over (k-pts,bands,spins)
1132 : !! Must be used in association with proc_distrb argument
1133 : !! natom=number of atoms in cell
1134 : !! nband=number of bands to import (usually 1 or nband_k)
1135 : !! nband_k=total number of bands for this k-point
1136 : !! [ncpgr]=(optional argument) second dimension of cprj%dcp(2,ncpgr,nlmn
1137 : !! stored in memory (mkmem>0) or present on disk (mkmem=0))
1138 : !! needed only when optional argument icpgr is present
1139 : !! nspinor=number of spinorial components of the wavefunctions (on current proc)
1140 : !! nsppol=1 for unpolarized, 2 for spin-polarized
1141 : !! [proc_distrb(nkpt,nband,nsppol)]=(optional argument) processor distribution
1142 : !! Describe how cprj datastructures are distributed over processors
1143 : !! When present, mpicomm argument must be also present
1144 : !! uncp=unit number for cprj data (used if mkmem=0)
1145 : !!
1146 : !! OUTPUT
1147 : !! cprj_k(dimcp,nspinor*nband) <type(pawcprj_type)>= output cprj datastructure
1148 : !!
1149 : !! SOURCE
1150 :
1151 9937816 : subroutine pawcprj_get(atind,cprj_k,cprj,dimcp,iband1,ibg,ikpt,iorder,isppol,mband,&
1152 : & mkmem,natom,nband,nband_k,nspinor,nsppol,uncp,&
1153 4968908 : & icpgr,ncpgr,mpicomm,proc_distrb) ! optionals arguments
1154 :
1155 : !Arguments ------------------------------------
1156 : !scalars
1157 : integer,intent(in) :: dimcp,iband1,ibg,ikpt,iorder,isppol,mband,mkmem,natom
1158 : integer,intent(in) :: nband,nband_k,nspinor,nsppol,uncp
1159 : integer,intent(in),optional :: icpgr,mpicomm,ncpgr
1160 : !arrays
1161 : integer,intent(in) :: atind(natom)
1162 : integer,intent(in),optional :: proc_distrb(:,:,:)
1163 : type(pawcprj_type),intent(in) :: cprj(dimcp,nspinor*mband*mkmem*nsppol)
1164 : type(pawcprj_type),intent(inout) :: cprj_k(dimcp,nspinor*nband)
1165 :
1166 : !Local variables-------------------------------
1167 : !scalars
1168 : integer :: iatm,iatom,ib,ibsp,icpgr_,isp,ispinor,jband,me,nband0,ncpgr_
1169 : logical :: has_distrb,has_icpgr
1170 : character(len=500) :: msg
1171 : !arrays
1172 4968908 : real(dp),allocatable :: tmp(:,:,:)
1173 :
1174 : ! *************************************************************************
1175 :
1176 4968908 : ncpgr_=cprj_k(1,1)%ncpgr;if (present(ncpgr)) ncpgr_=ncpgr
1177 4968908 : icpgr_=-1;if(present(icpgr)) icpgr_=icpgr
1178 4968908 : has_icpgr=(icpgr_>0.and.icpgr_<=ncpgr_)
1179 4968908 : if (present(icpgr).and.(.not.present(ncpgr))) then
1180 0 : msg='ncpgr must be present when icpgr is present (pawcprj_get)!'
1181 0 : LIBPAW_BUG(msg)
1182 : end if
1183 4968908 : if (has_icpgr.and.cprj_k(1,1)%ncpgr<1) then
1184 0 : msg='cprj_k%ncpgr not consistent with icpgr (pawcprj_get)!'
1185 0 : LIBPAW_BUG(msg)
1186 : end if
1187 :
1188 : !MPI data
1189 4968908 : has_distrb=present(proc_distrb)
1190 4968908 : if (has_distrb) then
1191 4201193 : if (.not.present(mpicomm)) then
1192 0 : msg='mpicomm must be present when proc_distrb is present (pawcprj_get)!'
1193 0 : LIBPAW_BUG(msg)
1194 : end if
1195 4201193 : me=xmpi_comm_rank(mpicomm)
1196 : end if
1197 :
1198 4968908 : if (mkmem==0) then
1199 :
1200 0 : if (iband1==1) then
1201 0 : read(uncp) nband0
1202 0 : if (nband_k/=nband0) then
1203 0 : msg='_PAW file was not created with the right options (pawcprj_get)!'
1204 0 : LIBPAW_BUG(msg)
1205 : end if
1206 : end if
1207 :
1208 0 : isp=0;jband=iband1-1
1209 0 : do ib=1,nband
1210 0 : jband=jband+1
1211 0 : if (has_distrb) then
1212 0 : if (abs(proc_distrb(ikpt,jband,isppol)-me)/=0) then
1213 0 : isp=isp+nspinor
1214 0 : cycle
1215 : end if
1216 : end if
1217 0 : do ispinor=1,nspinor
1218 0 : isp=isp+1
1219 0 : if (iorder==0) then
1220 0 : if (ncpgr_==0) then
1221 0 : do iatom=1,dimcp
1222 0 : read(uncp) cprj_k(iatom,isp)%cp(:,:)
1223 : end do
1224 : else
1225 0 : if (has_icpgr) then
1226 0 : do iatom=1,dimcp
1227 0 : LIBPAW_ALLOCATE(tmp,(2,ncpgr_,cprj_k(iatom,1)%nlmn))
1228 0 : read(uncp) cprj_k(iatom,isp)%cp(:,:),tmp(:,:,:)
1229 0 : cprj_k(iatom,isp)%dcp(:,1,:)=tmp(:,icpgr_,:)
1230 0 : LIBPAW_DEALLOCATE(tmp)
1231 : end do
1232 : else
1233 0 : do iatom=1,dimcp
1234 0 : read(uncp) cprj_k(iatom,isp)%cp(:,:),cprj_k(iatom,isp)%dcp(:,:,:)
1235 : end do
1236 : end if
1237 : end if
1238 : else
1239 0 : if (ncpgr_==0) then
1240 0 : do iatom=1,dimcp
1241 0 : iatm=min(atind(iatom),dimcp)
1242 0 : read(uncp) cprj_k(iatm,isp)%cp(:,:)
1243 : end do
1244 : else
1245 0 : if (has_icpgr) then
1246 0 : do iatom=1,dimcp
1247 0 : iatm=min(atind(iatom),dimcp)
1248 0 : LIBPAW_ALLOCATE(tmp,(2,ncpgr_,cprj_k(iatm,1)%nlmn))
1249 0 : read(uncp) cprj_k(iatm,isp)%cp(:,:),tmp(:,:,:)
1250 0 : cprj_k(iatm,isp)%dcp(:,1,:)=tmp(:,icpgr_,:)
1251 0 : LIBPAW_DEALLOCATE(tmp)
1252 : end do
1253 : else
1254 0 : do iatom=1,dimcp
1255 0 : iatm=min(atind(iatom),dimcp)
1256 0 : read(uncp) cprj_k(iatm,isp)%cp(:,:),cprj_k(iatm,isp)%dcp(:,:,:)
1257 : end do
1258 : end if
1259 : end if
1260 : end if
1261 : end do
1262 : end do
1263 :
1264 : else
1265 :
1266 4968908 : isp=0;ibsp=ibg+nspinor*(iband1-1);jband=iband1-1
1267 14340168 : do ib=1,nband
1268 9371260 : jband=jband+1
1269 9371260 : if (has_distrb) then
1270 8599984 : if (abs(proc_distrb(ikpt,jband,isppol)-me)/=0) then
1271 0 : isp=isp+nspinor;ibsp=ibsp+nspinor
1272 0 : cycle
1273 : end if
1274 : end if
1275 23948545 : do ispinor=1,nspinor
1276 9608377 : isp=isp+1;ibsp=ibsp+1
1277 18979637 : if (iorder==0) then
1278 7666873 : if (ncpgr_==0) then
1279 17242771 : do iatom=1,dimcp
1280 339684562 : cprj_k(iatom,isp)%cp(:,:)=cprj(iatom,ibsp)%cp(:,:)
1281 : end do
1282 : else
1283 2074560 : if (has_icpgr) then
1284 1791010 : do iatom=1,dimcp
1285 31196550 : cprj_k(iatom,isp)%cp(:,:) =cprj(iatom,ibsp)%cp(:,:)
1286 31777522 : cprj_k(iatom,isp)%dcp(:,1,:)=cprj(iatom,ibsp)%dcp(:,icpgr_,:)
1287 : end do
1288 : else
1289 4509012 : do iatom=1,dimcp
1290 76849600 : cprj_k(iatom,isp)%cp(:,:) =cprj(iatom,ibsp)%cp(:,:)
1291 602693188 : cprj_k(iatom,isp)%dcp(:,:,:)=cprj(iatom,ibsp)%dcp(:,:,:)
1292 : end do
1293 : end if
1294 : end if
1295 : else
1296 1941504 : if (ncpgr_==0) then
1297 5824512 : do iatom=1,dimcp
1298 3883008 : iatm=min(atind(iatom),dimcp)
1299 99016704 : cprj_k(iatm,isp)%cp(:,:)=cprj(iatom,ibsp)%cp(:,:)
1300 : end do
1301 : else
1302 0 : if (has_icpgr) then
1303 0 : do iatom=1,dimcp
1304 0 : iatm=min(atind(iatom),dimcp)
1305 0 : cprj_k(iatm,isp)%cp(:,:) =cprj(iatom,ibsp)%cp(:,:)
1306 0 : cprj_k(iatm,isp)%dcp(:,1,:)=cprj(iatom,ibsp)%dcp(:,icpgr_,:)
1307 : end do
1308 : else
1309 0 : do iatom=1,dimcp
1310 0 : iatm=min(atind(iatom),dimcp)
1311 0 : cprj_k(iatm,isp)%cp(:,:) =cprj(iatom,ibsp)%cp(:,:)
1312 0 : cprj_k(iatm,isp)%dcp(:,:,:)=cprj(iatom,ibsp)%dcp(:,:,:)
1313 : end do
1314 : end if
1315 : end if
1316 : end if
1317 : end do
1318 : end do
1319 :
1320 : end if
1321 :
1322 4968908 : end subroutine pawcprj_get
1323 : !!***
1324 :
1325 : !----------------------------------------------------------------------
1326 :
1327 : !!****f* m_pawcprj/pawcprj_put
1328 : !! NAME
1329 : !! pawcprj_put
1330 : !!
1331 : !! FUNCTION
1332 : !! Write cprj_k for a given set of (n,k) into memory in cprj, or into a temporary file
1333 : !!
1334 : !! INPUTS
1335 : !! atind(natom)=index table for atoms (see iorder below)
1336 : !! cprj_k(dimcp,nspinor*nband) <type(pawcprj_type)>= input cprj datastructure
1337 : !! dimcp=first dimension of cprj_k,cprjnk arrays (1 or natom)
1338 : !! iband1=index of first band in cprj
1339 : !! ibg=shift in cprj array to locate current k-point
1340 : !! ikpt=index of current k-point (only needed for the parallel distribution)
1341 : !! iorder=0 if cprj ordering does not change during reading
1342 : !! 1 if cprj ordering changes during writing, depending on content of atind array:
1343 : !! - if atind=atindx (type-sorted->unsorted)
1344 : !! - if atind=atindx1 (unsorted->type-sorted)
1345 : !! isppol=index of current spin component
1346 : !! mband=maximum number of bands
1347 : !! mkmem=number of k points which can fit in memory; set to 0 if use disk
1348 : !! [mpi_comm]=(optional argument) MPI communicator over (k-pts,bands,spins)
1349 : !! Must be used in association with proc_distrb argument
1350 : !! [mpi_comm_band]=(optional argument) MPI communicator over bands
1351 : !! Must be used in association with proc_distrb argument
1352 : !! natom=number of atoms in cell
1353 : !! nband=number of bands to export (usually 1, nband_k or nblockbd)
1354 : !! nband_k=total number of bands for this k-point
1355 : !! nlmn(dimcp)=array of dimensions of cprj_k,cprjnk datastructures
1356 : !! nspinor=number of spinorial components of the wavefunctions (on current proc)
1357 : !! nsppol=1 for unpolarized, 2 for spin-polarized
1358 : !! [proc_distrb(nkpt,nband,nsppol)]=(optional argument) processor distribution
1359 : !! Describe how cprj datastructures are distributed over processors
1360 : !! When present, mpicomm argument must be also present
1361 : !! [to_be_gathered]=(optional argument) TRUE if cprj_k arrays have to be
1362 : !! gathered between procs (band-fft parallelism only)
1363 : !! uncp=unit number for cprj data (used if mkmem=0)
1364 : !!
1365 : !! SIDE EFFECTS
1366 : !! cprj(dimcp,nspinor*mband*mkmem*nsppol)=output cprj (used if mkmem/=0)
1367 : !!
1368 : !! SOURCE
1369 :
1370 4850200 : subroutine pawcprj_put(atind,cprj_k,cprj,dimcp,iband1,ibg,ikpt,iorder,isppol,mband,&
1371 2425100 : & mkmem,natom,nband,nband_k,nlmn,nspinor,nsppol,uncp,&
1372 2425100 : & mpicomm,mpi_comm_band,proc_distrb,to_be_gathered) ! Optional arguments
1373 :
1374 : !Arguments ------------------------------------
1375 : !scalars
1376 : integer,intent(in) :: iband1,ibg,ikpt,iorder,isppol,dimcp,mband,mkmem
1377 : integer,intent(in) :: natom,nband,nband_k,nspinor,nsppol,uncp
1378 : integer,intent(in),optional :: mpicomm,mpi_comm_band
1379 : logical,optional,intent(in) :: to_be_gathered
1380 : !arrays
1381 : integer,intent(in) :: atind(natom),nlmn(dimcp)
1382 : integer,intent(in),optional :: proc_distrb(:,:,:)
1383 : type(pawcprj_type),intent(inout) :: cprj(dimcp,nspinor*mband*mkmem*nsppol)
1384 : type(pawcprj_type),intent(in) :: cprj_k(dimcp,nspinor*nband)
1385 :
1386 : !Local variables-------------------------------
1387 : !scalars
1388 : integer :: iatm,iatom,iband,ibsp,icpgr,ierr,ii,ilmn,isp,ispinor,jband,jj
1389 : integer :: lmndim,me,ncpgr,nproc_band
1390 : logical :: has_distrb,to_be_gathered_
1391 : character(len=500) :: msg
1392 : !arrays
1393 2425100 : real(dp),allocatable :: buffer1(:),buffer2(:)
1394 : ! *************************************************************************
1395 :
1396 2425100 : ncpgr=cprj_k(1,1)%ncpgr
1397 2425100 : to_be_gathered_=.false.;if (present(to_be_gathered)) to_be_gathered_=to_be_gathered
1398 :
1399 : !MPI data
1400 2425100 : nproc_band=1;if (present(mpi_comm_band)) nproc_band=xmpi_comm_size(mpi_comm_band)
1401 2425100 : has_distrb=present(proc_distrb)
1402 2425100 : if (has_distrb) then
1403 1577190 : if (.not.present(mpicomm)) then
1404 0 : msg='mpicomm must be present when proc_distrb is present (pawcprj_put)!'
1405 0 : LIBPAW_BUG(msg)
1406 : end if
1407 1577190 : me=xmpi_comm_rank(mpicomm)
1408 : end if
1409 :
1410 2425100 : if (nproc_band==1.or.(.not.to_be_gathered_)) then
1411 :
1412 2425100 : if (mkmem==0) then
1413 :
1414 0 : if (iband1==1) write(uncp) nband_k
1415 :
1416 0 : isp=0;jband=iband1-1
1417 0 : do iband=1,nband
1418 0 : jband=jband+1
1419 0 : if (has_distrb) then
1420 0 : if (abs(proc_distrb(ikpt,jband,isppol)-me)/=0) then
1421 0 : isp=isp+nspinor
1422 0 : cycle
1423 : end if
1424 : end if
1425 0 : do ispinor=1,nspinor
1426 0 : isp=isp+1
1427 0 : if (iorder==0) then
1428 0 : do iatom=1,dimcp
1429 0 : if (ncpgr==0) then
1430 0 : write(uncp) cprj_k(iatom,isp)%cp(:,:)
1431 : else
1432 0 : write(uncp) cprj_k(iatom,isp)%cp(:,:),cprj_k(iatom,isp)%dcp(:,:,:)
1433 : end if
1434 : end do
1435 : else
1436 0 : do iatom=1,dimcp
1437 0 : iatm=min(atind(iatom),dimcp)
1438 0 : if (ncpgr==0) then
1439 0 : write(uncp) cprj_k(iatm,isp)%cp(:,:)
1440 : else
1441 0 : write(uncp) cprj_k(iatm,isp)%cp(:,:),cprj_k(iatm,isp)%dcp(:,:,:)
1442 : end if
1443 : end do
1444 : end if
1445 : end do
1446 : end do
1447 :
1448 : else
1449 :
1450 2425100 : isp=0;ibsp=ibg+nspinor*(iband1-1);jband=iband1-1
1451 5019434 : do iband=1,nband
1452 2594334 : jband=jband+1
1453 2594334 : if (has_distrb) then
1454 1726456 : if (abs(proc_distrb(ikpt,jband,isppol)-me)/=0) then
1455 0 : isp=isp+nspinor;ibsp=ibsp+nspinor
1456 0 : cycle
1457 : end if
1458 : end if
1459 7831284 : do ispinor=1,nspinor
1460 2811850 : isp=isp+1;ibsp=ibsp+1
1461 5406184 : if (iorder==0) then
1462 7305898 : do iatom=1,dimcp
1463 172422692 : cprj(iatom,ibsp)%cp(:,:)=cprj_k(iatom,isp)%cp(:,:)
1464 40183902 : if (ncpgr>0) cprj(iatom,ibsp)%dcp(:,:,:)=cprj_k(iatom,isp)%dcp(:,:,:)
1465 : end do
1466 : else
1467 1869378 : do iatom=1,dimcp
1468 1243954 : iatm=min(atind(iatom),dimcp)
1469 31256260 : cprj(iatom,ibsp)%cp(:,:)=cprj_k(iatm,isp)%cp(:,:)
1470 10860610 : if (ncpgr>0) cprj(iatom,ibsp)%dcp(:,:,:)=cprj_k(iatm,isp)%dcp(:,:,:)
1471 : end do
1472 : end if
1473 : end do
1474 : end do
1475 :
1476 : end if
1477 :
1478 : else ! np_band>1
1479 :
1480 0 : lmndim=2*sum(nlmn(1:dimcp))*(1+ncpgr)*nspinor
1481 0 : LIBPAW_ALLOCATE(buffer1,(lmndim))
1482 0 : LIBPAW_ALLOCATE(buffer2,(lmndim*nproc_band))
1483 0 : isp=0;ibsp=ibg+nspinor*(iband1-1)
1484 0 : do iband=1,nband ! must be nblockbd for band-fft parallelism
1485 : jj=1
1486 0 : do ispinor=1,nspinor
1487 0 : isp=isp+1
1488 0 : do iatom=1,dimcp
1489 0 : if (iorder==0) then
1490 : iatm=iatom
1491 : else
1492 0 : iatm=min(atind(iatom),dimcp)
1493 : end if
1494 0 : do ilmn=1,nlmn(iatm)
1495 0 : buffer1(jj:jj+1)=cprj_k(iatm,isp)%cp(1:2,ilmn)
1496 0 : jj=jj+2
1497 : end do
1498 0 : if (ncpgr>0) then
1499 0 : do ilmn=1,nlmn(iatm)
1500 0 : do icpgr=1,ncpgr
1501 0 : buffer1(jj:jj+1)=cprj_k(iatm,isp)%dcp(1:2,icpgr,ilmn)
1502 0 : jj=jj+2
1503 : end do
1504 : end do
1505 : end if
1506 : end do !iatom
1507 : end do !ispinor
1508 0 : call xmpi_allgather(buffer1,lmndim,buffer2,mpi_comm_band,ierr)
1509 0 : jj=1
1510 0 : do ii=1,nproc_band
1511 0 : do ispinor=1,nspinor
1512 0 : ibsp=ibsp+1
1513 0 : do iatom=1,dimcp
1514 0 : if (iorder==0) then
1515 : iatm=iatom
1516 : else
1517 0 : iatm=min(atind(iatom),dimcp)
1518 : end if
1519 0 : do ilmn=1,nlmn(iatm)
1520 0 : cprj(iatom,ibsp)%cp(1:2,ilmn)=buffer2(jj:jj+1)
1521 0 : jj=jj+2
1522 : end do
1523 0 : if (ncpgr>0) then
1524 0 : do ilmn=1,nlmn(iatm)
1525 0 : do icpgr=1,ncpgr
1526 0 : cprj(iatom,ibsp)%dcp(1:2,icpgr,ilmn)=buffer2(jj:jj+1)
1527 0 : jj=jj+2
1528 : end do
1529 : end do
1530 : end if
1531 : end do !iatom
1532 : end do !ispinor
1533 : end do !ii=1,nproc_band
1534 : end do !iband
1535 0 : LIBPAW_DEALLOCATE(buffer1)
1536 0 : LIBPAW_DEALLOCATE(buffer2)
1537 :
1538 : end if ! mode_para=b, nband
1539 :
1540 2425100 : end subroutine pawcprj_put
1541 : !!***
1542 :
1543 : !----------------------------------------------------------------------
1544 :
1545 : !!****f* m_pawcprj/pawcprj_reorder
1546 : !! NAME
1547 : !! pawcprj_reorder
1548 : !!
1549 : !! FUNCTION
1550 : !! Change the order of a cprj datastructure
1551 : !! From unsorted cprj to atom-sorted cprj (atm_indx=atindx)
1552 : !! From atom-sorted cprj to unsorted cprj (atm_indx=atindx1)
1553 : !!
1554 : !! INPUTS
1555 : !! atm_indx(natom)=index table for atoms
1556 : !! From unsorted cprj to atom-sorted cprj (atm_indx=atindx)
1557 : !! From atom-sorted cprj to unsorted cprj (atm_indx=atindx1)
1558 : !!
1559 : !! OUTPUT
1560 : !!
1561 : !! SIDE EFFECTS
1562 : !! cprj(:,:) <type(pawcprj_type)>= cprj datastructure
1563 : !!
1564 : !! SOURCE
1565 :
1566 844 : subroutine pawcprj_reorder(cprj,atm_indx)
1567 :
1568 : !Arguments ------------------------------------
1569 : !scalars
1570 : !arrays
1571 : integer,intent(in) :: atm_indx(:)
1572 : type(pawcprj_type),intent(inout) :: cprj(:,:)
1573 :
1574 : !Local variables-------------------------------
1575 : !scalars
1576 : integer :: iexit,ii,jj,kk,n1atindx,n1cprj,n2cprj,ncpgr
1577 : character(len=100) :: msg
1578 : !arrays
1579 844 : integer,allocatable :: nlmn(:)
1580 844 : type(pawcprj_type),allocatable :: cprj_tmp(:,:)
1581 :
1582 : ! *************************************************************************
1583 :
1584 844 : n1cprj=size(cprj,dim=1);n2cprj=size(cprj,dim=2)
1585 844 : n1atindx=size(atm_indx,dim=1)
1586 844 : if (n1cprj==0.or.n2cprj==0.or.n1atindx<=1) return
1587 :
1588 819 : if (n1cprj/=n1atindx) then
1589 0 : msg='wrong sizes (pawcprj_reorder)!'
1590 0 : LIBPAW_BUG(msg)
1591 : end if
1592 :
1593 : !Nothing to do when the atoms are already sorted
1594 819 : iexit=1;ii=0
1595 2501 : do while (iexit==1.and.ii<n1atindx)
1596 1682 : ii=ii+1
1597 2501 : if (atm_indx(ii)/=ii) iexit=0
1598 : end do
1599 819 : if (iexit==1) return
1600 :
1601 0 : LIBPAW_ALLOCATE(nlmn,(n1cprj))
1602 0 : do ii=1,n1cprj
1603 0 : nlmn(ii)=cprj(ii,1)%nlmn
1604 : end do
1605 0 : ncpgr=cprj(1,1)%ncpgr
1606 0 : LIBPAW_DATATYPE_ALLOCATE(cprj_tmp,(n1cprj,n2cprj))
1607 0 : call pawcprj_alloc(cprj_tmp,ncpgr,nlmn)
1608 0 : call pawcprj_copy(cprj,cprj_tmp)
1609 0 : call pawcprj_free(cprj)
1610 :
1611 0 : do jj=1,n2cprj
1612 0 : do ii=1,n1cprj
1613 0 : kk=atm_indx(ii)
1614 0 : cprj(kk,jj)%nlmn=nlmn(ii)
1615 0 : cprj(kk,jj)%ncpgr=ncpgr
1616 0 : LIBPAW_ALLOCATE(cprj(kk,jj)%cp,(2,nlmn(ii)))
1617 0 : cprj(kk,jj)%cp(:,:)=cprj_tmp(ii,jj)%cp(:,:)
1618 0 : if (ncpgr>0) then
1619 0 : LIBPAW_ALLOCATE(cprj(kk,jj)%dcp,(2,ncpgr,nlmn(ii)))
1620 0 : cprj(kk,jj)%dcp(:,:,:)=cprj_tmp(ii,jj)%dcp(:,:,:)
1621 : end if
1622 : end do
1623 : end do
1624 :
1625 0 : call pawcprj_free(cprj_tmp)
1626 0 : LIBPAW_DATATYPE_DEALLOCATE(cprj_tmp)
1627 0 : LIBPAW_DEALLOCATE(nlmn)
1628 :
1629 844 : end subroutine pawcprj_reorder
1630 : !!***
1631 :
1632 : !----------------------------------------------------------------------
1633 :
1634 : !!****f* m_pawcprj/pawcprj_mpi_exch
1635 : !! NAME
1636 : !! pawcprj_mpi_exch
1637 : !!
1638 : !! FUNCTION
1639 : !! Exchange a pawcprj_type between two processors inside a MPI communicator.
1640 : !!
1641 : !! INPUTS
1642 : !! natom=Number of atoms (size of first dimension of Cprj_send and Cprj_recv).
1643 : !! n2dim=Size of the second dimension.
1644 : !! nlmn(natom)=Number of nlm partial waves for each atom.
1645 : !! Cprj_send= The datatype to be transmitted.
1646 : !! receiver=ID of the receiver in spaceComm.
1647 : !! sender=ID of the sender in spaceComm.
1648 : !! spaceComm=MPI Communicator.
1649 : !! mtag= message tag
1650 : !!
1651 : !! OUTPUT
1652 : !! ierr=Error status.
1653 : !! Cprj_recv=The datatype copied on proc. receiver.
1654 : !!
1655 : !! NOTES
1656 : !! If sender==receiver, Cprj_send is copied into Cprj_recv.
1657 : !! It should be easy to avoid this additional copy in the calling routine.
1658 : !!
1659 : !! SOURCE
1660 :
1661 16 : subroutine pawcprj_mpi_exch(natom,n2dim,nlmn,ncpgr,Cprj_send,Cprj_recv,sender,receiver,spaceComm,mtag,ierr)
1662 :
1663 : !Arguments ------------------------------------
1664 : !scalars
1665 : integer,intent(in) :: mtag,natom,n2dim,ncpgr
1666 : integer,intent(in) :: sender,receiver,spaceComm
1667 : integer,intent(out) :: ierr
1668 : !arrays
1669 : integer,intent(in) :: nlmn(natom)
1670 : type(pawcprj_type),intent(in) :: Cprj_send(:,:)
1671 : type(pawcprj_type),intent(inout) :: Cprj_recv(:,:)
1672 :
1673 : !Local variables-------------------------------
1674 : !scalars
1675 : integer :: iat,jj,t2dim,tcpgr,n1dim,nn
1676 : integer :: ntotcp,ipck,rank
1677 : character(len=500) :: msg
1678 : !arrays
1679 16 : real(dp),allocatable :: buffer_cp(:,:),buffer_cpgr(:,:,:)
1680 :
1681 : ! *************************************************************************
1682 :
1683 16 : n1dim=0
1684 16 : t2dim=0
1685 16 : tcpgr=0
1686 16 : ierr=0
1687 16 : if (sender==receiver) then
1688 0 : call pawcprj_copy(Cprj_send,Cprj_recv)
1689 : return
1690 : end if
1691 :
1692 16 : rank = xmpi_comm_rank(spaceComm)
1693 :
1694 16 : nn=size(nlmn,dim=1)
1695 16 : if (rank==sender) then
1696 8 : n1dim=size(Cprj_send,dim=1)
1697 8 : t2dim=size(Cprj_send,dim=2)
1698 8 : tcpgr=Cprj_send(1,1)%ncpgr
1699 : end if
1700 16 : if (rank==receiver) then
1701 8 : n1dim=size(Cprj_recv,dim=1)
1702 8 : t2dim=size(Cprj_recv,dim=2)
1703 8 : tcpgr=Cprj_recv(1,1)%ncpgr
1704 : end if
1705 16 : if (rank/=sender.and.rank/=receiver) then
1706 : write(msg,'(a,3i0)') &
1707 0 : & 'rank is not equal to sender or receiver (pawcprj_mpi_exch): ',rank, sender, receiver
1708 0 : LIBPAW_BUG(msg)
1709 : end if
1710 :
1711 32 : ntotcp=n2dim*SUM(nlmn(:))
1712 :
1713 48 : LIBPAW_ALLOCATE(buffer_cp,(2,ntotcp))
1714 16 : if (ncpgr/=0) then
1715 0 : LIBPAW_ALLOCATE(buffer_cpgr,(2,ncpgr,ntotcp))
1716 : end if
1717 :
1718 : !=== Pack Cprj_send ===
1719 16 : if (rank==sender) then
1720 : ipck=0
1721 168 : do jj=1,n2dim
1722 328 : do iat=1,natom
1723 160 : nn=nlmn(iat)
1724 8800 : buffer_cp(:,ipck+1:ipck+nn)=Cprj_send(iat,jj)%cp(:,1:nn)
1725 160 : if (ncpgr/=0) buffer_cpgr(:,:,ipck+1:ipck+nn)=Cprj_send(iat,jj)%dcp(:,:,1:nn)
1726 320 : ipck=ipck+nn
1727 : end do
1728 : end do
1729 : end if
1730 :
1731 : !=== Transmit data ===
1732 16 : call xmpi_exch(buffer_cp,2*ntotcp,sender,buffer_cp,receiver,spaceComm,2*mtag,ierr)
1733 16 : if (ncpgr/=0) then
1734 0 : call xmpi_exch(buffer_cpgr,2*ncpgr*ntotcp,sender,buffer_cpgr,receiver,spaceComm,2*mtag+1,ierr)
1735 : end if
1736 :
1737 : !=== UnPack buffers into Cprj_recv ===
1738 16 : if (rank==receiver) then
1739 : ipck=0
1740 168 : do jj=1,n2dim
1741 328 : do iat=1,natom
1742 160 : nn=nlmn(iat)
1743 8800 : Cprj_recv(iat,jj)%cp(:,1:nn)=buffer_cp(:,ipck+1:ipck+nn)
1744 160 : if (ncpgr/=0) Cprj_recv(iat,jj)%dcp(:,:,1:nn)=buffer_cpgr(:,:,ipck+1:ipck+nn)
1745 320 : ipck=ipck+nn
1746 : end do
1747 : end do
1748 : end if
1749 :
1750 16 : LIBPAW_DEALLOCATE(buffer_cp)
1751 16 : if (ncpgr/=0) then
1752 0 : LIBPAW_DEALLOCATE(buffer_cpgr)
1753 : end if
1754 :
1755 16 : end subroutine pawcprj_mpi_exch
1756 : !!***
1757 :
1758 : !----------------------------------------------------------------------
1759 :
1760 : !!****f* m_pawcprj/pawcprj_mpi_send
1761 : !! NAME
1762 : !! pawcprj_mpi_send
1763 : !!
1764 : !! FUNCTION
1765 : !! Send a pawcprj_type inside a MPI communicator.
1766 : !!
1767 : !! INPUTS
1768 : !! natom=Number of atoms (size of first dimension of cprj_out).
1769 : !! n2dim=Size of the second dimension.
1770 : !! nlmn(natom)=Number of nlm partial waves for each atom.
1771 : !! ncpgr = number of gradients in cprj_out
1772 : !! cprj_out= The datatype to be transmitted.
1773 : !! receiver=ID of the receiver in spaceComm.
1774 : !! spaceComm=MPI Communicator.
1775 : !!
1776 : !! OUTPUT
1777 : !! ierr=Error status.
1778 : !!
1779 : !! NOTES
1780 : !! perhaps in general it is more efficient to use pawcprj_mpi_exch but it is
1781 : !! convenient for coding to have separate send and recieve routines.
1782 : !!
1783 : !! SOURCE
1784 :
1785 40884 : subroutine pawcprj_mpi_send(natom,n2dim,nlmn,ncpgr,cprj_out,receiver,spaceComm,ierr)
1786 :
1787 : !Arguments ------------------------------------
1788 : !scalars
1789 : integer,intent(in) :: natom,n2dim,ncpgr
1790 : integer,intent(in) :: receiver,spaceComm
1791 : integer,intent(out) :: ierr
1792 : !arrays
1793 : integer,intent(in) :: nlmn(natom)
1794 : type(pawcprj_type),intent(in) :: cprj_out(:,:)
1795 :
1796 : !Local variables-------------------------------
1797 : !scalars
1798 : integer :: iat,jj,t2dim,tcpgr,n1dim,nn
1799 : integer :: ntotcp,ipck,tag
1800 : character(len=100) :: msg
1801 : !arrays
1802 40884 : real(dp),allocatable :: buffer_cp(:,:),buffer_cpgr(:,:,:)
1803 :
1804 : ! *************************************************************************
1805 :
1806 40884 : n1dim=0
1807 40884 : t2dim=0
1808 40884 : tcpgr=0
1809 40884 : ierr=0
1810 :
1811 40884 : nn=size(nlmn,dim=1)
1812 40884 : n1dim=size(cprj_out,dim=1)
1813 40884 : t2dim=size(cprj_out,dim=2)
1814 40884 : tcpgr=cprj_out(1,1)%ncpgr
1815 :
1816 40884 : if (nn/=n1dim) then
1817 0 : msg='size mismatch in natom (pawcprj_mpi_send)!'
1818 0 : LIBPAW_BUG(msg)
1819 : end if
1820 40884 : if (t2dim/=n2dim) then
1821 0 : msg='size mismatch in dim=2 (pawcprj_mpi_send)!'
1822 0 : LIBPAW_BUG(msg)
1823 : end if
1824 40884 : if (tcpgr/=ncpgr) then
1825 0 : msg='size mismatch in ncpgr (pawcprj_mpi_send)!'
1826 0 : LIBPAW_BUG(msg)
1827 : end if
1828 :
1829 122652 : ntotcp=n2dim*SUM(nlmn(:))
1830 :
1831 122652 : LIBPAW_ALLOCATE(buffer_cp,(2,ntotcp))
1832 40884 : if (ncpgr/=0) then
1833 142080 : LIBPAW_ALLOCATE(buffer_cpgr,(2,ncpgr,ntotcp))
1834 : end if
1835 :
1836 : !=== Pack cprj_out ====
1837 : ipck=0
1838 204420 : do jj=1,n2dim
1839 531492 : do iat=1,natom
1840 327072 : nn=nlmn(iat)
1841 8176800 : buffer_cp(:,ipck+1:ipck+nn)=cprj_out(iat,jj)%cp(:,1:nn)
1842 63978912 : if (ncpgr/=0) buffer_cpgr(:,:,ipck+1:ipck+nn)=cprj_out(iat,jj)%dcp(:,:,1:nn)
1843 490608 : ipck=ipck+nn
1844 : end do
1845 : end do
1846 :
1847 : !=== Transmit data ===
1848 40884 : tag = 2*ntotcp
1849 40884 : call xmpi_send(buffer_cp,receiver,tag,spaceComm,ierr)
1850 40884 : if (ncpgr/=0) then
1851 35520 : tag=tag*ncpgr
1852 35520 : call xmpi_send(buffer_cpgr,receiver,tag,spaceComm,ierr)
1853 : end if
1854 :
1855 : !=== Clean up ===
1856 40884 : LIBPAW_DEALLOCATE(buffer_cp)
1857 40884 : if (ncpgr/=0) then
1858 35520 : LIBPAW_DEALLOCATE(buffer_cpgr)
1859 : end if
1860 :
1861 81768 : end subroutine pawcprj_mpi_send
1862 : !!***
1863 :
1864 : !----------------------------------------------------------------------
1865 :
1866 : !!****f* m_pawcprj/pawcprj_mpi_recv
1867 : !! NAME
1868 : !! pawcprj_mpi_recv
1869 : !!
1870 : !! FUNCTION
1871 : !! Receive a pawcprj_type inside a MPI communicator.
1872 : !!
1873 : !! INPUTS
1874 : !! natom=Number of atoms (size of first dimension of Cprj_in).
1875 : !! n2dim=Size of the second dimension.
1876 : !! nlmn(natom)=Number of nlm partial waves for each atom.
1877 : !! ncpgr = number of gradients in cprj_in
1878 : !! sender=ID of the sender in spaceComm.
1879 : !! spaceComm=MPI Communicator.
1880 : !!
1881 : !! OUTPUT
1882 : !! ierr=Error status.
1883 : !! cprj_in=The datatype copied on proc. receiver.
1884 : !!
1885 : !! NOTES
1886 : !! Perhaps in general it is more efficient to use pawcprj_mpi_exch but it is
1887 : !! convenient for coding to have separate send and receive routines.
1888 : !!
1889 : !! SOURCE
1890 :
1891 40884 : subroutine pawcprj_mpi_recv(natom,n2dim,nlmn,ncpgr,cprj_in,sender,spaceComm,ierr)
1892 :
1893 : !Arguments ------------------------------------
1894 : !scalars
1895 : integer,intent(in) :: natom,n2dim,ncpgr
1896 : integer,intent(in) :: sender,spaceComm
1897 : integer,intent(out) :: ierr
1898 : !arrays
1899 : integer,intent(in) :: nlmn(natom)
1900 : type(pawcprj_type),intent(inout) :: cprj_in(:,:)
1901 :
1902 : !Local variables-------------------------------
1903 : !scalars
1904 : integer :: iat,jj,t2dim,tcpgr,n1dim,nn
1905 : integer :: ntotcp,ipck,tag
1906 : character(len=100) :: msg
1907 : !arrays
1908 40884 : real(dp),allocatable :: buffer_cp(:,:),buffer_cpgr(:,:,:)
1909 :
1910 : ! *************************************************************************
1911 :
1912 40884 : n1dim=0
1913 40884 : t2dim=0
1914 40884 : tcpgr=0
1915 40884 : ierr=0
1916 :
1917 40884 : nn=size(nlmn,dim=1)
1918 40884 : n1dim=size(cprj_in,dim=1)
1919 40884 : t2dim=size(cprj_in,dim=2)
1920 40884 : tcpgr=cprj_in(1,1)%ncpgr
1921 :
1922 40884 : if (nn/=n1dim) then
1923 0 : msg='size mismatch in natom (pawcprj_mpi_recv)!'
1924 0 : LIBPAW_BUG(msg)
1925 : end if
1926 40884 : if (t2dim/=n2dim) then
1927 0 : msg='size mismatch in dim=2 (pawcprj_mpi_recv)!'
1928 0 : LIBPAW_BUG(msg)
1929 : end if
1930 40884 : if (tcpgr/=ncpgr) then
1931 0 : msg='size mismatch in ncpgr (pawcprj_mpi_recv)!'
1932 0 : LIBPAW_BUG(msg)
1933 : end if
1934 :
1935 122652 : ntotcp=n2dim*SUM(nlmn(:))
1936 :
1937 122652 : LIBPAW_ALLOCATE(buffer_cp,(2,ntotcp))
1938 40884 : if (ncpgr/=0) then
1939 142080 : LIBPAW_ALLOCATE(buffer_cpgr,(2,ncpgr,ntotcp))
1940 : end if
1941 :
1942 : !=== Receive data ===
1943 40884 : tag = 2*ntotcp
1944 40884 : call xmpi_recv(buffer_cp,sender,tag,spaceComm,ierr)
1945 40884 : if (ncpgr/=0) then
1946 35520 : tag=tag*ncpgr
1947 35520 : call xmpi_recv(buffer_cpgr,sender,tag,spaceComm,ierr)
1948 : end if
1949 :
1950 : !=== UnPack buffers into cprj_in ===
1951 : ipck=0
1952 204420 : do jj=1,n2dim
1953 531492 : do iat=1,natom
1954 327072 : nn=nlmn(iat)
1955 8176800 : cprj_in(iat,jj)%cp(:,1:nn)=buffer_cp(:,ipck+1:ipck+nn)
1956 63978912 : if (ncpgr/=0) cprj_in(iat,jj)%dcp(:,:,1:nn)=buffer_cpgr(:,:,ipck+1:ipck+nn)
1957 490608 : ipck=ipck+nn
1958 : end do
1959 : end do
1960 :
1961 : !=== Clean up ===
1962 40884 : LIBPAW_DEALLOCATE(buffer_cp)
1963 40884 : if (ncpgr/=0) then
1964 35520 : LIBPAW_DEALLOCATE(buffer_cpgr)
1965 : end if
1966 :
1967 81768 : end subroutine pawcprj_mpi_recv
1968 : !!***
1969 :
1970 : !----------------------------------------------------------------------
1971 :
1972 : !!****f* m_pawcprj/pawcprj_mpi_sum
1973 : !! NAME
1974 : !! pawcprj_mpi_sum
1975 : !!
1976 : !! FUNCTION
1977 : !! Perform MPI_SUM on a pawcprj_type inside a MPI communicator.
1978 : !!
1979 : !! INPUTS
1980 : !! spaceComm=MPI Communicator.
1981 : !!
1982 : !! SIDE EFFECTS
1983 : !! cprj=the cprj datastructure
1984 : !! ierr=Error status.
1985 : !!
1986 : !! SOURCE
1987 :
1988 880558 : subroutine pawcprj_mpi_sum(cprj,spaceComm,ierr)
1989 :
1990 : !Arguments ------------------------------------
1991 : !scalars
1992 : integer,intent(in) :: spaceComm
1993 : integer,intent(out) :: ierr
1994 : !arrays
1995 : type(pawcprj_type),intent(inout) :: cprj(:,:)
1996 :
1997 : !Local variables-------------------------------
1998 : !scalars
1999 : integer,parameter :: maxBytes=100*1024*1024 ! 100 MBytes
2000 : integer :: ii,ipck,jj,ncpgr,nlmn,nn,n1dim,n2dim,n2dim1,n2dim2,sizeBytes,step
2001 : logical,parameter :: save_memory=.true.
2002 : !arrays
2003 880558 : real(dp),allocatable :: buffer_cprj(:,:,:)
2004 :
2005 : ! *************************************************************************
2006 :
2007 880558 : if (xmpi_comm_size(spaceComm)<2) return
2008 :
2009 12000 : n1dim=size(cprj,1);n2dim=size(cprj,2)
2010 132000 : nlmn=sum(cprj(:,:)%nlmn)
2011 132000 : ncpgr=maxval(cprj(:,:)%ncpgr)
2012 :
2013 12000 : step=n2dim
2014 : if (save_memory) then
2015 12000 : sizeBytes=2*(1+ncpgr)*nlmn *8
2016 12000 : step=n2dim/max(1,sizeBytes/maxBytes)
2017 12000 : if (step==0) step=1
2018 : end if
2019 :
2020 12000 : do n2dim1=1,n2dim,step
2021 :
2022 12000 : n2dim2=min(n2dim1+step-1,n2dim)
2023 132000 : nlmn=sum(cprj(:,n2dim1:n2dim2)%nlmn)
2024 60000 : LIBPAW_ALLOCATE(buffer_cprj,(2,1+ncpgr,nlmn))
2025 :
2026 3468000 : ipck=0 ; buffer_cprj=zero
2027 24000 : do jj=n2dim1,n2dim2
2028 132000 : do ii=1,n1dim
2029 108000 : nn=cprj(ii,jj)%nlmn
2030 2700000 : buffer_cprj(:,1,ipck+1:ipck+nn)=cprj(ii,jj)%cp(:,1:nn)
2031 108000 : if (cprj(ii,jj)%ncpgr/=0) buffer_cprj(:,2:1+ncpgr,ipck+1:ipck+nn)=cprj(ii,jj)%dcp(:,1:ncpgr,1:nn)
2032 120000 : ipck=ipck+nn
2033 : end do
2034 : end do
2035 :
2036 12000 : call xmpi_sum(buffer_cprj,spaceComm,ierr)
2037 :
2038 12000 : ipck=0
2039 24000 : do jj=n2dim1,n2dim2
2040 132000 : do ii=1,n1dim
2041 108000 : nn=cprj(ii,jj)%nlmn
2042 2700000 : cprj(ii,jj)%cp(:,1:nn)=buffer_cprj(:,1,ipck+1:ipck+nn)
2043 108000 : if (cprj(ii,jj)%ncpgr/=0) cprj(ii,jj)%dcp(:,1:ncpgr,1:nn)=buffer_cprj(:,2:1+ncpgr,ipck+1:ipck+nn)
2044 120000 : ipck=ipck+nn
2045 : end do
2046 : end do
2047 :
2048 24000 : LIBPAW_DEALLOCATE(buffer_cprj)
2049 :
2050 : end do
2051 :
2052 12000 : end subroutine pawcprj_mpi_sum
2053 : !!***
2054 :
2055 : !----------------------------------------------------------------------
2056 :
2057 : !!****f* m_pawcprj/pawcprj_mpi_allgather
2058 : !! NAME
2059 : !! pawcprj_mpi_allgather
2060 : !!
2061 : !! FUNCTION
2062 : !! Perform MPI_ALLGATHER on a pawcprj_type inside a MPI communicator.
2063 : !!
2064 : !! INPUTS
2065 : !! cprj_loc= The cprj on the local proc being all-gathered
2066 : !! natom=Number of atoms (size of first dimension of cprj_loc).
2067 : !! n2dim=Size of the second dimension of cprj_loc.
2068 : !! n2std=Stride of n2 dimension
2069 : !! if n2std=1, cprj_loc(:,1) is on proc 0, cprj_loc(:,2) is on proc 1, cprj_loc(:,3) is on proc 2, etc.
2070 : !! if n2std>1, cprj_loc(:,1:n2std) are on proc 0, cprj_loc(:,n2std+1,2*n2std) are on proc 1, etc.
2071 : !! nlmn(natom)=Number of nlm partial waves for each atom.
2072 : !! ncpgr = number of gradients in cprj_loc
2073 : !! nproc=number of processors being gathered
2074 : !! spaceComm=MPI Communicator.
2075 : !! [rank_ordered]= optional, default=FALSE
2076 : !! TRUE: second dimension of gathered datastructure is rank-ordered
2077 : !! FALSE: second dimension of gathered datastructure is not rank-ordered
2078 : !!
2079 : !! OUTPUT
2080 : !! cprj_gat=the gathered cprjs
2081 : !! ierr=Error status.
2082 : !!
2083 : !! SOURCE
2084 :
2085 30042 : subroutine pawcprj_mpi_allgather(cprj_loc,cprj_gat,natom,n2dim,n2std,nlmn,ncpgr,nproc,spaceComm,ierr,&
2086 : & rank_ordered)
2087 :
2088 : !Arguments ------------------------------------
2089 : !scalars
2090 : integer,intent(in) :: natom,n2dim,n2std,ncpgr,nproc,spaceComm
2091 : integer,intent(out) :: ierr
2092 : logical,optional,intent(in) :: rank_ordered
2093 : !arrays
2094 : integer,intent(in) :: nlmn(natom)
2095 : type(pawcprj_type),intent(in) :: cprj_loc(:,:)
2096 : type(pawcprj_type),intent(inout) :: cprj_gat(:,:)
2097 :
2098 : !Local variables-------------------------------
2099 : !scalars
2100 : integer :: iat,ii,jj,t2dim,tcpgr,tg2dim,n1dim,nn
2101 : integer :: ntotcp,ibuf,ipck,iproc
2102 : logical :: rank_ordered_
2103 : character(len=100) :: msg
2104 : !arrays
2105 30042 : real(dp),allocatable :: buffer_cpgr(:,:,:),buffer_cpgr_all(:,:,:)
2106 :
2107 : ! *************************************************************************
2108 :
2109 30042 : n1dim=0
2110 30042 : t2dim=0
2111 30042 : tg2dim=0
2112 30042 : tcpgr=0
2113 30042 : ierr=0
2114 :
2115 30042 : nn=size(nlmn,dim=1)
2116 30042 : n1dim=size(cprj_loc,dim=1)
2117 30042 : t2dim=size(cprj_loc,dim=2)
2118 30042 : tg2dim=size(cprj_gat,dim=2)
2119 30042 : tcpgr=cprj_loc(1,1)%ncpgr
2120 :
2121 30042 : if (nn/=n1dim) then
2122 0 : msg='size mismatch in natom (pawcprj_mpi_allgather)!'
2123 0 : LIBPAW_BUG(msg)
2124 : end if
2125 30042 : if (t2dim/=n2dim) then
2126 0 : msg='size mismatch in dim=2 (pawcprj_mpi_allgather)!'
2127 0 : LIBPAW_BUG(msg)
2128 : end if
2129 30042 : if (tg2dim/=n2dim*nproc) then
2130 0 : msg='size mismatch in dim=2 (pawcprj_mpi_allgather)!'
2131 0 : LIBPAW_BUG(msg)
2132 : end if
2133 30042 : if (tcpgr/=ncpgr) then
2134 0 : msg='size mismatch in ncpgr (pawcprj_mpi_allgather)!'
2135 0 : LIBPAW_BUG(msg)
2136 : end if
2137 30042 : if (mod(n2dim,n2std)/=0) then
2138 0 : msg='n2std should divide n2dim (pawcprj_mpi_allgather)!'
2139 0 : LIBPAW_BUG(msg)
2140 : end if
2141 :
2142 30042 : rank_ordered_=.false.;if(present(rank_ordered)) rank_ordered_=rank_ordered
2143 :
2144 90126 : ntotcp=n2dim*SUM(nlmn(:))
2145 120168 : LIBPAW_ALLOCATE(buffer_cpgr,(2,1+ncpgr,ntotcp))
2146 120168 : LIBPAW_ALLOCATE(buffer_cpgr_all,(2,1+ncpgr,nproc*ntotcp))
2147 :
2148 : !=== Pack cprj_loc ====
2149 30042 : ipck=0
2150 150250 : do jj=1,n2dim
2151 390666 : do iat=1,natom
2152 240416 : nn=nlmn(iat)
2153 6010400 : buffer_cpgr(:,1,ipck+1:ipck+nn)=cprj_loc(iat,jj)%cp(:,1:nn)
2154 9211168 : if (ncpgr/=0) buffer_cpgr(:,2:1+ncpgr,ipck+1:ipck+nn)=cprj_loc(iat,jj)%dcp(:,:,1:nn)
2155 360624 : ipck=ipck+nn
2156 : end do
2157 : end do
2158 :
2159 : !=== allgather data ===
2160 30042 : call xmpi_allgather(buffer_cpgr,2*(ncpgr+1)*ntotcp,buffer_cpgr_all,spaceComm,ierr)
2161 :
2162 : !=== unpack gathered data into cprj(natom,n2dim*nproc)
2163 : !=== second dimension is rank-ordered if rank_ordered_=true
2164 30042 : ipck=0
2165 93444 : do iproc=1,nproc
2166 231084 : do jj=1,n2dim/n2std
2167 454650 : do ii=1,n2std
2168 253608 : if (rank_ordered_) then
2169 98024 : ibuf=(iproc-1)*n2dim+(jj-1)*n2std+ii
2170 : else
2171 155584 : ibuf=(iproc+(jj-1)*nproc-1)*n2std+ii
2172 : end if
2173 898464 : do iat=1,natom
2174 507216 : nn=nlmn(iat)
2175 12680400 : cprj_gat(iat,ibuf)%cp(:,1:nn)=buffer_cpgr_all(:,1,ipck+1:ipck+nn)
2176 507216 : if (ncpgr/=0) cprj_gat(iat,ibuf)%dcp(:,1:ncpgr,1:nn)=&
2177 20761200 : & buffer_cpgr_all(:,2:1+ncpgr,ipck+1:ipck+nn)
2178 760824 : ipck=ipck+nn
2179 : end do
2180 : end do
2181 : end do
2182 : end do
2183 :
2184 : !=== Clean up ===
2185 30042 : LIBPAW_DEALLOCATE(buffer_cpgr)
2186 30042 : LIBPAW_DEALLOCATE(buffer_cpgr_all)
2187 :
2188 30042 : end subroutine pawcprj_mpi_allgather
2189 : !!***
2190 :
2191 : !----------------------------------------------------------------------
2192 :
2193 : !!****f* m_pawcprj/pawcprj_bcast
2194 : !! NAME
2195 : !! pawcprj_bcast
2196 : !!
2197 : !! FUNCTION
2198 : !! Broadcast a pawcprj_type from master to all nodes inside a MPI communicator.
2199 : !!
2200 : !! INPUTS
2201 : !! natom=Number of atoms (size of the first dimension of Cprj).
2202 : !! n2dim=Size of the second dimension of Cprj.
2203 : !! ncpgr=Number of gradients that have to be cast. It is a bit redundant but, it can be used to
2204 : !! broad cast only the %cp"s without caring about the gradients. Just set it to 0 but be careful!
2205 : !! nlmn(natom)=Number of nlm partial waves for each atom.
2206 : !! master=ID of the sending node in spaceComm.
2207 : !! spaceComm=MPI Communicator.
2208 : !!
2209 : !! OUTPUT
2210 : !! ierr=Error status.
2211 : !! Cprj(natom,n2dim)<pawcprj_type>=The datatype to be transmitted by master and received by the others nodes.
2212 : !!
2213 : !! SOURCE
2214 :
2215 2470 : subroutine pawcprj_bcast(Cprj,natom,n2dim,nlmn,ncpgr,master,spaceComm,ierr)
2216 :
2217 : !Arguments ------------------------------------
2218 : !scalars
2219 : integer,intent(in) :: natom,n2dim,ncpgr,master,spaceComm
2220 : integer,intent(out) :: ierr
2221 : !arrays
2222 : integer,intent(in) :: nlmn(natom)
2223 : type(pawcprj_type),intent(inout) :: Cprj(natom,n2dim)
2224 :
2225 : !Local variables-------------------------------
2226 : !scalars
2227 : integer :: iat,jj,n1dim,nn
2228 : integer :: ntotcp,ipck,rank,nprocs
2229 : character(len=100) :: msg
2230 : !arrays
2231 2470 : real(dp),allocatable :: buffer_cp(:,:),buffer_cpgr(:,:,:)
2232 :
2233 : ! *************************************************************************
2234 :
2235 2470 : ierr=0
2236 2470 : nprocs = xmpi_comm_size(spaceComm)
2237 2470 : if (nprocs==1) return
2238 :
2239 0 : rank = xmpi_comm_rank(spaceComm)
2240 :
2241 0 : nn=size(nlmn,dim=1)
2242 0 : n1dim=size(Cprj,dim=1)
2243 : if (nn/=n1dim) then
2244 : msg='size mismatch in natom (pawcprj_bcast)!'
2245 : LIBPAW_BUG(msg)
2246 : end if
2247 :
2248 0 : ntotcp=n2dim*SUM(nlmn(:))
2249 :
2250 0 : LIBPAW_ALLOCATE(buffer_cp,(2,ntotcp))
2251 0 : if (ncpgr/=0) then
2252 0 : LIBPAW_ALLOCATE(buffer_cpgr,(2,ncpgr,ntotcp))
2253 : end if
2254 :
2255 : !=== Master packs Cprj ===
2256 : !Write a routine to pack/unpack?
2257 0 : if (rank==master) then
2258 : ipck=0
2259 0 : do jj=1,n2dim
2260 0 : do iat=1,natom
2261 0 : nn=nlmn(iat)
2262 0 : buffer_cp(:,ipck+1:ipck+nn)=Cprj(iat,jj)%cp(:,1:nn)
2263 0 : if (ncpgr/=0) buffer_cpgr(:,:,ipck+1:ipck+nn)=Cprj(iat,jj)%dcp(:,:,1:nn)
2264 0 : ipck=ipck+nn
2265 : end do
2266 : end do
2267 : end if
2268 :
2269 : !=== Transmit data ===
2270 0 : call xmpi_bcast(buffer_cp,master,spaceComm,ierr)
2271 0 : if (ncpgr/=0) then
2272 0 : call xmpi_bcast(buffer_cpgr,master,spaceComm,ierr)
2273 : end if
2274 :
2275 : !=== UnPack the received buffer ===
2276 0 : if (rank/=master) then
2277 : ipck=0
2278 0 : do jj=1,n2dim
2279 0 : do iat=1,natom
2280 0 : nn=nlmn(iat)
2281 0 : Cprj(iat,jj)%cp(:,1:nn)=buffer_cp(:,ipck+1:ipck+nn)
2282 0 : if (ncpgr/=0) Cprj(iat,jj)%dcp(:,:,1:nn)=buffer_cpgr(:,:,ipck+1:ipck+nn)
2283 0 : ipck=ipck+nn
2284 : end do
2285 : end do
2286 : end if
2287 :
2288 0 : LIBPAW_DEALLOCATE(buffer_cp)
2289 0 : if (ncpgr/=0) then
2290 0 : LIBPAW_DEALLOCATE(buffer_cpgr)
2291 : end if
2292 :
2293 0 : end subroutine pawcprj_bcast
2294 : !!***
2295 :
2296 : !----------------------------------------------------------------------
2297 :
2298 : !!****f* m_pawcprj/pawcprj_transpose
2299 : !! NAME
2300 : !! pawcprj_transpose
2301 : !!
2302 : !! FUNCTION
2303 : !! Transpose a cprj datastructure FOR A GIVEN (K,SPIN)
2304 : !! in order to change the parallel distribution from atom to band (or the contrary).
2305 : !! At input, cprj is distributed over bands (or atoms); at output, it is distributed over atoms (or bands)
2306 : !!
2307 : !! INPUTS
2308 : !! cprjin(n1indim,n2indim)<pawcprj_type>=the input cprj datastructure
2309 : !! cprj_bandpp=number of bands to be treated simultaneoulsy by a processor
2310 : !! natom=number of atoms in cell
2311 : !! nband=number of bands
2312 : !! nspinor=number of spinorial components
2313 : !! spaceComm=MPI Communicator.
2314 : !!
2315 : !! OUTPUT
2316 : !! cprjout(n1outdim,n2outdim)<pawcprj_type>=the output cprj datastructure with another distribution
2317 : !!
2318 : !! NOTES
2319 : !! On the dimensions:
2320 : !! To transfer cprj from band distribution to atom distribution, dimensions should be:
2321 : !! n1indim =natom n2indim =nband/nproc*nspinor
2322 : !! n1outdim=natom/nproc n2outdim=nband*nspinor
2323 : !! To transfer cprj from atom distribution to band distribution, dimensions should be:
2324 : !! n1indim =natom n2indim =nband/nproc*nspinor
2325 : !! n1outdim=natom/nproc n2outdim=nband*nspinor
2326 : !!
2327 : !! SOURCE
2328 :
2329 0 : subroutine pawcprj_transpose(cprjin,cprjout,cprj_bandpp,natom,nband,nspinor,spaceComm)
2330 :
2331 : !Arguments-------------------------------------
2332 : !scalars
2333 : integer :: cprj_bandpp,natom,nband,nspinor,spaceComm
2334 : !arrays
2335 : type(pawcprj_type),intent(in) :: cprjin(:,:)
2336 : type(pawcprj_type),intent(out) :: cprjout(:,:)
2337 :
2338 : !Local variables-------------------------------
2339 : !scalars
2340 : integer :: bpp,buf_indx
2341 : integer :: iashft,iatom,iatom_max_sd,iatom_max_rc,iatom_1,iatom_2,iatm1_sd,iatm1_rc,iatm2_sd,iatm2_rc
2342 : integer :: ib,iband,iband_1,iband_2,iband_shift,iblock_atom,iblock_band,ibshft
2343 : integer :: ierr,ip,ispinor,me,nba,nbb,nbnp_sd,nbnp_rc,ncpgr,nlmn,np
2344 : integer :: rbufsize,sbufsize,size11,size12,size21,size22,transpose_mode
2345 : character(len=100) :: msg
2346 : !arrays
2347 0 : integer,allocatable :: cprjsz_atom(:),cprjsz_block(:,:)
2348 0 : integer,allocatable,target :: count_atom(:),count_band(:),displ_atom(:),displ_band(:)
2349 0 : integer,pointer :: scount(:),sdispl(:),rcount(:),rdispl(:)
2350 0 : real(dp),allocatable :: rbuf(:),sbuf(:)
2351 :
2352 : ! *************************************************************************
2353 :
2354 : !MPI data
2355 0 : me = xmpi_comm_rank(spaceComm)
2356 0 : np = xmpi_comm_size(spaceComm)
2357 :
2358 : !Nothing to do if nprocs=1
2359 0 : if (np==1) then
2360 0 : call pawcprj_copy(cprjin,cprjout)
2361 0 : return
2362 : end if
2363 :
2364 : !Compute bloc sizes
2365 0 : bpp=cprj_bandpp
2366 0 : nba=natom/np;if (mod(natom,np)/=0) nba=nba+1
2367 0 : nbb=nband/(np*bpp)
2368 :
2369 : !Check sizes, select direction of transposition
2370 0 : transpose_mode=0
2371 0 : size11=size(cprjin,1);size12=size(cprjin,2)
2372 0 : size21=size(cprjout,1);size22=size(cprjout,2)
2373 : if (size11==natom.and.size12==nbb*bpp*nspinor.and.&
2374 0 : & size21==nba.and.size22==nband*nspinor) then
2375 : transpose_mode=1
2376 : else if (size11==nba.and.size12==nband*nspinor.and.&
2377 0 : & size21==natom.and.size22==nbb*bpp*nspinor) then
2378 : else
2379 0 : msg='wrong cprjin/cprjout sizes (pawcprj_transpose)!'
2380 0 : LIBPAW_BUG(msg)
2381 : end if
2382 :
2383 : !Compute size of atom bloc (wr to cprj)
2384 0 : LIBPAW_ALLOCATE(cprjsz_atom,(natom))
2385 0 : LIBPAW_ALLOCATE(cprjsz_block,(np,nba))
2386 0 : cprjsz_atom=0;cprjsz_block=0
2387 0 : if (transpose_mode==1) then
2388 0 : do iatom=1,natom
2389 0 : cprjsz_atom(iatom)=2*cprjin(iatom,1)%nlmn*(1+cprjin(iatom,1)%ncpgr)
2390 : end do
2391 : else
2392 0 : do iblock_atom=1,nba
2393 0 : iatom=(iblock_atom-1)*np+1+me
2394 0 : if (iatom<=natom) cprjsz_atom(iatom)=2*cprjin(iblock_atom,1)%nlmn*(1+cprjin(iblock_atom,1)%ncpgr)
2395 : end do
2396 0 : call xmpi_sum(cprjsz_atom,spaceComm,ierr)
2397 : end if
2398 0 : do iblock_atom=1,nba
2399 0 : iashft=(iblock_atom-1)*np
2400 0 : iatom_1=iashft+1;iatom_2=iashft+np
2401 0 : if (iatom_1>natom) cycle
2402 0 : if (iatom_2>natom) iatom_2=natom
2403 0 : do iatom=iatom_1,iatom_2
2404 0 : cprjsz_block(iatom-iashft,iblock_atom)=cprjsz_atom(iatom)+2 ! +2 for nlmn et ncpgr
2405 : end do
2406 : end do
2407 0 : LIBPAW_DEALLOCATE(cprjsz_atom)
2408 :
2409 : !Allocations for MPI_ALLTOALL
2410 0 : LIBPAW_ALLOCATE(count_atom,(np))
2411 0 : LIBPAW_ALLOCATE(displ_atom,(np))
2412 0 : LIBPAW_ALLOCATE(count_band,(np))
2413 0 : LIBPAW_ALLOCATE(displ_band,(np))
2414 :
2415 : !Loop on blocks of bands
2416 0 : do iblock_band=1,nbb !(note: np divides nband)
2417 0 : ibshft=(iblock_band-1)*np*bpp
2418 0 : iband_1=ibshft+1;iband_2=ibshft+np*bpp
2419 0 : if (iband_1>nband.or.iband_2>nband) cycle ! for security
2420 :
2421 : ! Loop on blocks of atoms
2422 0 : do iblock_atom=1,nba
2423 0 : iashft=(iblock_atom-1)*np
2424 0 : iatom_1=iashft+1;iatom_2=iashft+np
2425 0 : if (iatom_1>natom) cycle
2426 0 : if (iatom_2>natom) iatom_2=natom
2427 :
2428 : ! Computation of displacements and sizes of blocks when data are band-distributed
2429 0 : count_band(1)=cprjsz_block(1,iblock_atom)*nspinor*bpp;displ_band(1)=0
2430 0 : do ip=2,np
2431 0 : count_band(ip)=cprjsz_block(ip,iblock_atom)*nspinor*bpp
2432 0 : displ_band(ip)=displ_band(ip-1)+count_band(ip-1)
2433 : end do
2434 :
2435 : ! Computation of displacements and sizes of blocks when data are atom-distributed
2436 0 : count_atom(1)=cprjsz_block(1+me,iblock_atom)*bpp*nspinor;displ_atom(1)=0
2437 0 : do ip=2,np
2438 0 : count_atom(ip)=count_atom(1)
2439 0 : displ_atom(ip)=displ_atom(ip-1)+count_atom(ip-1)
2440 : end do
2441 :
2442 : ! According to transposition mode, select
2443 : ! - displacements and sizes of blocks
2444 : ! - shifts in arrays
2445 0 : if (transpose_mode==1) then
2446 0 : scount => count_band ; sdispl => displ_band
2447 0 : rcount => count_atom ; rdispl => displ_atom
2448 0 : nbnp_sd=bpp;nbnp_rc=np*bpp
2449 0 : iatm1_sd=iatom_1;iatm2_sd=iatom_2
2450 0 : iatm1_rc=iblock_atom;iatm2_rc=iatm1_rc
2451 0 : iatom_max_sd=iatom_2;iatom_max_rc=iashft+1+me
2452 : else
2453 0 : scount => count_atom ; sdispl => displ_atom
2454 0 : rcount => count_band ; rdispl => displ_band
2455 0 : nbnp_sd=np*bpp;nbnp_rc=bpp
2456 0 : iatm1_sd=iblock_atom;iatm2_sd=iatm1_sd
2457 0 : iatm1_rc=iatom_1;iatm2_rc=iatom_2
2458 0 : iatom_max_sd=iashft+1+me;iatom_max_rc=iatom_2
2459 : end if
2460 :
2461 : ! Allocation of buffers
2462 0 : sbufsize=sdispl(np)+scount(np)
2463 0 : rbufsize=rdispl(np)+rcount(np)
2464 0 : LIBPAW_ALLOCATE(sbuf,(sbufsize))
2465 0 : LIBPAW_ALLOCATE(rbuf,(rbufsize))
2466 :
2467 : ! Coying of input cprj to buffer for sending
2468 0 : buf_indx=0
2469 0 : iband_shift=(iblock_band-1)*nbnp_sd-1
2470 0 : if (iatom_max_sd<=natom) then
2471 0 : do iatom=iatm1_sd,iatm2_sd
2472 0 : do ib=1,nbnp_sd
2473 0 : iband=(iband_shift+ib)*nspinor
2474 0 : do ispinor=1,nspinor
2475 0 : iband=iband+1
2476 0 : nlmn=cprjin(iatom,iband)%nlmn;ncpgr=cprjin(iatom,iband)%ncpgr
2477 0 : sbuf(buf_indx+1)=dble(nlmn) ;buf_indx=buf_indx+1
2478 0 : sbuf(buf_indx+1)=dble(ncpgr);buf_indx=buf_indx+1
2479 0 : sbuf(buf_indx+1:buf_indx+2*nlmn)=reshape(cprjin(iatom,iband)%cp(1:2,1:nlmn),(/2*nlmn/))
2480 0 : buf_indx=buf_indx+2*nlmn
2481 0 : if (ncpgr>0) then
2482 0 : sbuf(buf_indx+1:buf_indx+2*ncpgr*nlmn)=reshape(cprjin(iatom,iband)%dcp(1:2,1:ncpgr,1:nlmn),(/2*ncpgr*nlmn/))
2483 0 : buf_indx=buf_indx+2*ncpgr*nlmn
2484 : end if
2485 : end do
2486 : end do
2487 : end do
2488 : end if
2489 0 : if (buf_indx/=sbufsize) then
2490 0 : msg='wrong buffer size for sending (pawcprj_transpose)!'
2491 0 : LIBPAW_BUG(msg)
2492 : end if
2493 :
2494 : ! Main call to MPI_ALLTOALL
2495 0 : call xmpi_alltoallv(sbuf,scount,sdispl,rbuf,rcount,rdispl,spaceComm,ierr)
2496 :
2497 : ! Retrieving of output cprj for received buffer
2498 0 : buf_indx=0
2499 0 : iband_shift=(iblock_band-1)*nbnp_rc-1
2500 0 : if (iatom_max_rc<=natom) then
2501 0 : do iatom=iatm1_rc,iatm2_rc
2502 0 : do ib=1,nbnp_rc
2503 0 : iband=(iband_shift+ib)*nspinor
2504 0 : do ispinor=1,nspinor
2505 0 : iband=iband+1
2506 0 : nlmn =int(rbuf(buf_indx+1));buf_indx=buf_indx+1
2507 0 : ncpgr=int(rbuf(buf_indx+1));buf_indx=buf_indx+1
2508 0 : cprjout(iatom,iband)%nlmn=nlmn;cprjout(iatom,iband)%ncpgr=ncpgr
2509 0 : cprjout(iatom,iband)%cp(1:2,1:nlmn)=reshape(rbuf(buf_indx+1:buf_indx+2*nlmn),(/2,nlmn/))
2510 0 : buf_indx=buf_indx+2*nlmn
2511 0 : if (ncpgr>0) then
2512 0 : cprjout(iatom,iband)%dcp(1:2,1:ncpgr,1:nlmn)=reshape(rbuf(buf_indx+1:buf_indx+2*nlmn*ncpgr),(/2,ncpgr,nlmn/))
2513 0 : buf_indx=buf_indx+2*nlmn*ncpgr
2514 : end if
2515 : end do
2516 : end do
2517 : end do
2518 : else
2519 0 : cprjout(iatom,iband)%nlmn=0;cprjout(iatom,iband)%ncpgr=0
2520 : end if
2521 0 : if (buf_indx/=rbufsize) then
2522 0 : msg='wrong buffer size for receiving (pawcprj_transpose)!'
2523 0 : LIBPAW_BUG(msg)
2524 : end if
2525 :
2526 : ! Deallocation of buffers
2527 0 : LIBPAW_DEALLOCATE(sbuf)
2528 0 : LIBPAW_DEALLOCATE(rbuf)
2529 :
2530 : ! End of loops
2531 : end do ! do iblock_atom
2532 : end do ! do iblock_atom
2533 :
2534 : !Free memory
2535 0 : LIBPAW_DEALLOCATE(count_atom)
2536 0 : LIBPAW_DEALLOCATE(displ_atom)
2537 0 : LIBPAW_DEALLOCATE(count_band)
2538 0 : LIBPAW_DEALLOCATE(displ_band)
2539 0 : LIBPAW_DEALLOCATE(cprjsz_block)
2540 0 : nullify(scount,rcount,sdispl,rdispl)
2541 :
2542 0 : end subroutine pawcprj_transpose
2543 : !!***
2544 :
2545 : !----------------------------------------------------------------------
2546 :
2547 : !!****f* m_pawcprj/pawcprj_gather_spin
2548 : !! NAME
2549 : !! pawcprj_gather_spin
2550 : !!
2551 : !! FUNCTION
2552 : !!
2553 : !! INPUTS
2554 : !! cprj(:,:)=the input cprj datastructure
2555 : !! n2size=number of cprj datastructures to be gathered (second dim)
2556 : !! nspinor : number of spinorial component (on current proc)
2557 : !! nspinortot : total number of spinorial component
2558 : !!
2559 : !! OUTPUT
2560 : !! cprj_gat(:,:) = the cprj containing all nspinor componants
2561 : !!
2562 : !! NOTES
2563 : !! The cprj has been built like the following:
2564 : !! loop on nsppol
2565 : !! loop on k point
2566 : !! loop over band or block of band
2567 : !! These quantities were build only if treated by the current proc
2568 : !! the inner quantities being nspinor
2569 : !!
2570 : !! SOURCE
2571 464 : subroutine pawcprj_gather_spin(cprj,cprj_gat,natom,n2size,nspinor,nspinortot,&
2572 : & spaceComm_spin,ierr)
2573 :
2574 : !Arguments ------------------------------------
2575 : !scalars
2576 : integer,intent(in) :: natom,nspinor,nspinortot,n2size
2577 : integer,intent(in) :: spaceComm_spin
2578 : integer,intent(out) :: ierr
2579 : !arrays
2580 : type(pawcprj_type),intent(in) :: cprj(:,:)
2581 : type(pawcprj_type),intent(inout) :: cprj_gat(:,:)
2582 :
2583 : !Local variables-------------------------------
2584 : !scalars
2585 : integer :: i1,iatom,ibsp,icpgr,ilmn,isp,ispinor,jj,lmndim,n2dim,n2dim_gat,ncpgr
2586 : character(len=100) :: msg
2587 : !arrays
2588 928 : integer :: nlmn(natom)
2589 464 : real(dp),allocatable :: buffer1(:),buffer2(:)
2590 :
2591 : ! *************************************************************************
2592 :
2593 464 : n2dim =size(cprj,dim=2)
2594 464 : n2dim_gat=size(cprj_gat,dim=2)
2595 464 : if (n2dim_gat/=(nspinortot/nspinor)*n2dim) then
2596 0 : msg='wrong dims (pawcprj_gather_spin)!'
2597 0 : LIBPAW_BUG(msg)
2598 : end if
2599 :
2600 1392 : do iatom=1,natom
2601 1392 : nlmn(iatom)=size(cprj(iatom,1)%cp(1,:))
2602 : end do
2603 464 : ncpgr=cprj(1,1)%ncpgr
2604 1392 : lmndim=2*n2size*sum(nlmn(1:natom))*(1+ncpgr)
2605 1392 : LIBPAW_ALLOCATE(buffer1,(lmndim))
2606 1392 : LIBPAW_ALLOCATE(buffer2,(lmndim*nspinortot))
2607 :
2608 464 : isp=0;ibsp=0
2609 464 : jj=1
2610 8432 : do i1=1,n2size
2611 7968 : isp=isp+1
2612 24368 : do iatom=1,natom
2613 302784 : do ilmn=1,nlmn(iatom)
2614 860544 : buffer1(jj:jj+1)=cprj(iatom,isp)%cp(1:2,ilmn)
2615 302784 : jj=jj+2
2616 : end do
2617 23904 : if (ncpgr>0) then
2618 0 : do ilmn=1,nlmn(iatom)
2619 0 : do icpgr=1,ncpgr
2620 0 : buffer1(jj:jj+1)=cprj(iatom,isp)%dcp(1:2,icpgr,ilmn)
2621 0 : jj=jj+2
2622 : end do
2623 : end do
2624 : end if
2625 : end do
2626 : end do
2627 :
2628 464 : call xmpi_allgather(buffer1,lmndim,buffer2,spaceComm_spin,ierr)
2629 :
2630 464 : jj=1
2631 1392 : do ispinor=1,nspinortot
2632 17328 : do i1 =1,n2size
2633 15936 : ibsp=(i1-1)*nspinortot + ispinor
2634 48736 : do iatom=1,natom
2635 605568 : do ilmn=1,nlmn(iatom)
2636 1721088 : cprj_gat(iatom,ibsp)%cp(1:2,ilmn)=buffer2(jj:jj+1)
2637 605568 : jj=jj+2
2638 : end do
2639 47808 : if (ncpgr>0) then
2640 0 : do ilmn=1,nlmn(iatom)
2641 0 : do icpgr=1,ncpgr
2642 0 : cprj_gat(iatom,ibsp)%dcp(1:2,icpgr,ilmn)=buffer2(jj:jj+1)
2643 0 : jj=jj+2
2644 : end do
2645 : end do
2646 : end if
2647 : end do
2648 : end do
2649 : end do
2650 :
2651 464 : LIBPAW_DEALLOCATE(buffer1)
2652 464 : LIBPAW_DEALLOCATE(buffer2)
2653 :
2654 464 : end subroutine pawcprj_gather_spin
2655 : !!***
2656 :
2657 : !----------------------------------------------------------------------
2658 :
2659 : !!****f* m_pawcprj/pawcprj_getdim
2660 : !! NAME
2661 : !! pawcprj_getdim
2662 : !!
2663 : !! FUNCTION
2664 : !! Helper function returning the number of lmn components in the <p_{lmn}^i|\psi> for the i-th atom.
2665 : !! Used to initialize the dimensioning array that is passed to the pawcprj_alloc routines when the
2666 : !! pawcprj_type structure is allocated and initialized.
2667 : !!
2668 : !! INPUTS
2669 : !! natom=number of atoms in the unit cell
2670 : !! nattyp(ntypat)=number of atoms of each type
2671 : !! ntypat=number of atom types
2672 : !! typat(natom-= type of each atom
2673 : !! Pawtab(ntypat)<pawtab_type>=PAW tabulated starting data.
2674 : !! sort_mode(len=*)=String defining the sorting of the atoms in the Cprj arrays.
2675 : !! Two modes are possible:
2676 : !! -- "O[rdered]", if atoms are sorted by atom type.
2677 : !! -- "R[andom]", if atoms are sorted randomly i.e. according the values of typat specified in the input file.
2678 : !!
2679 : !! OUTPUT
2680 : !! dimcprj(natom)=Number of nlm elements in the <p_{lmn}^i|\psi> matrix elements for i=1,...,natom.
2681 : !!
2682 : !! SOURCE
2683 :
2684 42032 : subroutine pawcprj_getdim(dimcprj,natom,nattyp,ntypat,typat,Pawtab,sort_mode)
2685 :
2686 : !Arguments ------------------------------------
2687 : integer,intent(in) :: natom,ntypat
2688 : character(len=*),intent(in) :: sort_mode
2689 : !arrays
2690 : integer,intent(in) :: nattyp(:),typat(natom)
2691 : integer,intent(inout) :: dimcprj(natom)
2692 : type(Pawtab_type),intent(in) :: Pawtab(ntypat)
2693 :
2694 : !Local variables-------------------------------
2695 : integer :: iatom,itypat
2696 : character(len=500) :: msg
2697 :
2698 : ! *************************************************************************
2699 :
2700 : SELECT CASE (sort_mode(1:1))
2701 :
2702 : CASE ("o","O") ! Ordered by atom-type
2703 :
2704 : iatom=0
2705 91412 : do itypat=1,ntypat
2706 147927 : dimcprj(iatom+1:iatom+nattyp(itypat))=Pawtab(itypat)%lmn_size
2707 91412 : iatom=iatom+nattyp(itypat)
2708 : end do
2709 :
2710 : CASE ("r","R") ! Randomly ordered (typat from input file)
2711 :
2712 12513 : do iatom=1,natom
2713 7931 : itypat=typat(iatom)
2714 12513 : dimcprj(iatom)=Pawtab(itypat)%lmn_size
2715 : end do
2716 :
2717 : CASE DEFAULT
2718 0 : msg='Wrong value for sort_mode: '//TRIM(sort_mode)
2719 42032 : LIBPAW_ERROR(msg)
2720 : END SELECT
2721 :
2722 42032 : end subroutine pawcprj_getdim
2723 : !!***
2724 :
2725 : !----------------------------------------------------------------------
2726 :
2727 : !!****f* m_pawcprj/paw_overlap
2728 : !! NAME
2729 : !! paw_overlap
2730 : !!
2731 : !! FUNCTION
2732 : !! Helper function returning the onsite contribution to the overlap between two states.
2733 : !!
2734 : !! INPUTS
2735 : !! spinor_comm= (optional) communicator over spinorial components
2736 : !! typat(:)=The type of each atom.
2737 : !! Pawtab(ntypat)<type(pawtab_type)>=paw tabulated starting data.
2738 : !! cprj1,cprj2<pawcprj_type>
2739 : !! Projected wave functions <Proj_i|Cnk> with all NL projectors for the left and the right wavefunction,respectively.
2740 : !!
2741 : !! OUTPUT
2742 : !!
2743 : !! SOURCE
2744 :
2745 2310 : function paw_overlap(cprj1,cprj2,typat,pawtab,spinor_comm) result(onsite)
2746 :
2747 : !Arguments ------------------------------------
2748 : !scalars
2749 : integer,intent(in),optional :: spinor_comm
2750 : !arrays
2751 : integer,intent(in) :: typat(:)
2752 : real(dp) :: onsite(2)
2753 : type(pawcprj_type),intent(in) :: cprj1(:,:),cprj2(:,:)
2754 : type(pawtab_type),intent(in) :: pawtab(:)
2755 :
2756 : !Local variables-------------------------------
2757 : !scalars
2758 : integer :: iatom,ilmn,itypat,j0lmn,jlmn,klmn,natom,nspinor,isp
2759 : real(dp) :: sij
2760 : character(len=500) :: msg
2761 : !arrays
2762 :
2763 : ! *************************************************************************
2764 :
2765 2310 : natom=SIZE(typat)
2766 :
2767 2310 : if (SIZE(cprj1,DIM=1)/=SIZE(cprj2,DIM=1) .or. SIZE(cprj1,DIM=1)/=natom) then
2768 0 : write(msg,'(a,3i4)')' Wrong size in typat, cprj1, cprj2 : ',natom,SIZE(cprj1),SIZE(cprj2)
2769 0 : LIBPAW_ERROR(msg)
2770 : end if
2771 :
2772 2310 : nspinor = SIZE(cprj1,DIM=2)
2773 :
2774 6930 : onsite=zero
2775 4620 : do iatom=1,natom
2776 2310 : itypat=typat(iatom)
2777 46200 : do jlmn=1,pawtab(itypat)%lmn_size
2778 41580 : j0lmn=jlmn*(jlmn-1)/2
2779 438900 : do ilmn=1,jlmn
2780 395010 : klmn=j0lmn+ilmn
2781 395010 : sij=pawtab(itypat)%sij(klmn); if (jlmn==ilmn) sij=sij*half
2782 436590 : if (ABS(sij)>tol16) then
2783 124740 : do isp=1,nspinor
2784 :
2785 : onsite(1)=onsite(1) + sij*( &
2786 : & cprj1(iatom,isp)%cp(1,ilmn) * cprj2(iatom,isp)%cp(1,jlmn) &
2787 : & +cprj1(iatom,isp)%cp(2,ilmn) * cprj2(iatom,isp)%cp(2,jlmn) &
2788 : & +cprj1(iatom,isp)%cp(1,jlmn) * cprj2(iatom,isp)%cp(1,ilmn) &
2789 : & +cprj1(iatom,isp)%cp(2,jlmn) * cprj2(iatom,isp)%cp(2,ilmn) &
2790 62370 : & )
2791 :
2792 : onsite(2)=onsite(2) + sij*( &
2793 : & cprj1(iatom,isp)%cp(1,ilmn) * cprj2(iatom,isp)%cp(2,jlmn) &
2794 : & -cprj1(iatom,isp)%cp(2,ilmn) * cprj2(iatom,isp)%cp(1,jlmn) &
2795 : & +cprj1(iatom,isp)%cp(1,jlmn) * cprj2(iatom,isp)%cp(2,ilmn) &
2796 : & -cprj1(iatom,isp)%cp(2,jlmn) * cprj2(iatom,isp)%cp(1,ilmn) &
2797 124740 : & )
2798 : end do
2799 : end if
2800 : end do
2801 : end do
2802 : end do
2803 :
2804 2310 : if (present(spinor_comm)) then
2805 2310 : call xmpi_sum(onsite,spinor_comm,isp)
2806 : end if
2807 :
2808 2310 : end function paw_overlap
2809 : !!***
2810 :
2811 : !----------------------------------------------------------------------
2812 :
2813 : !!****f* m_pawcprj/pawcprj_pack
2814 : !! NAME
2815 : !! pawcprj_pack
2816 : !!
2817 : !! FUNCTION
2818 : !! Pack structured data into a simple buffer
2819 : !!
2820 : !! INPUTS
2821 : !! nlmn(natom)=Number of nlm partial waves for each atom.
2822 : !! ncpgr = number of gradients in cprj_out
2823 : !! cprj= The datatype to be packed.
2824 : !!
2825 : !! OUTPUT
2826 : !! buffer = the data packed, dim : (2, n2dim*sum(nlmn))
2827 : !! [buffer_gr] = if present the gradient data packed, dim : (2, ncpgr, n2dim*sum(nlmn))
2828 : !!
2829 : !! SOURCE
2830 :
2831 600 : subroutine pawcprj_pack(nlmn,cprj,buffer,buffer_gr)
2832 :
2833 : !Arguments ------------------------------------
2834 : !scalars
2835 : !arrays
2836 : integer,intent(in) :: nlmn(:)
2837 : type(pawcprj_type),intent(in) :: cprj(:,:)
2838 : real(dp),intent(out) :: buffer(:,:)
2839 : real(dp),intent(out),optional :: buffer_gr(:,:,:)
2840 :
2841 : !Local variables-------------------------------
2842 : !scalars
2843 : integer :: natom,n2buffer,ncpgr,n2dim
2844 : integer :: iat,jj,n1dim,nn
2845 : integer :: ipck
2846 : character(len=100) :: msg
2847 : !arrays
2848 :
2849 : ! *************************************************************************
2850 :
2851 600 : natom=size(nlmn,dim=1)
2852 600 : n2buffer=size(buffer,dim=2)
2853 600 : n1dim=size(cprj,dim=1)
2854 600 : n2dim=size(cprj,dim=2)
2855 :
2856 600 : if (natom/=n1dim) then
2857 0 : msg='size mismatch in natom (pawcprj_pack)!'
2858 0 : LIBPAW_BUG(msg)
2859 : end if
2860 3600 : if (n2dim*SUM(nlmn)/=n2buffer) then
2861 0 : msg='size mismatch in dim=2 (pawcprj_pack)!'
2862 0 : LIBPAW_BUG(msg)
2863 : end if
2864 600 : ncpgr=0
2865 600 : if (present(buffer_gr)) then
2866 0 : ncpgr=size(buffer_gr,dim=2)
2867 : end if
2868 :
2869 : !=== Pack cprj ====
2870 600 : ipck=0
2871 1200 : do jj=1,n2dim
2872 4200 : do iat=1,natom
2873 3000 : nn=nlmn(iat)
2874 118200 : buffer(:,ipck+1:ipck+nn)=cprj(iat,jj)%cp(:,1:nn)
2875 3000 : if (ncpgr/=0) then
2876 0 : buffer_gr(:,:,ipck+1:ipck+nn)=cprj(iat,jj)%dcp(:,:,1:nn)
2877 : end if
2878 3600 : ipck=ipck+nn
2879 : end do
2880 : end do
2881 :
2882 600 : end subroutine pawcprj_pack
2883 : !!***
2884 :
2885 : !----------------------------------------------------------------------
2886 :
2887 : !!****f* m_pawcprj/pawcprj_unpack
2888 : !! NAME
2889 : !! pawcprj_unpack
2890 : !!
2891 : !! FUNCTION
2892 : !! Unpack structured data from a simple buffer
2893 : !!
2894 : !! INPUTS
2895 : !! nlmn(natom)=Number of nlm partial waves for each atom.
2896 : !! ncpgr = number of gradients in cprj_in
2897 : !! buffer = the data to be unpacked, dim : (2, n2dim*sum(nlmn))
2898 : !! [buffer_gr] = if present the gradient data to be unpacked, dim : (2, ncpgr, n2dim*sum(nlmn))
2899 : !!
2900 : !! OUTPUT
2901 : !! cprj=The datatype unpacked
2902 : !!
2903 : !! SOURCE
2904 :
2905 1508 : subroutine pawcprj_unpack(nlmn,cprj,buffer,buffer_gr)
2906 :
2907 : !Arguments ------------------------------------
2908 : !scalars
2909 : !arrays
2910 : integer,intent(in) :: nlmn(:)
2911 : real(dp),intent(in) :: buffer(:,:)
2912 : real(dp),intent(in),optional :: buffer_gr(:,:,:)
2913 : type(pawcprj_type),intent(inout) :: cprj(:,:)
2914 :
2915 : !Local variables-------------------------------
2916 : !scalars
2917 : integer :: natom,n2buffer,ncpgr,n2dim
2918 : integer :: iat,jj,n1dim,nn
2919 : integer :: ipck
2920 : character(len=100) :: msg
2921 : !arrays
2922 :
2923 : ! *************************************************************************
2924 :
2925 1508 : natom=size(nlmn,dim=1)
2926 1508 : n2buffer=size(buffer,dim=2)
2927 1508 : n1dim=size(cprj,dim=1)
2928 1508 : n2dim=size(cprj,dim=2)
2929 :
2930 1508 : if (natom/=n1dim) then
2931 0 : msg='size mismatch in natom (pawcprj_unpack)!'
2932 0 : LIBPAW_BUG(msg)
2933 : end if
2934 9048 : if (n2dim*SUM(nlmn)/=n2buffer) then
2935 0 : msg='size mismatch in dim=2 (pawcprj_unpack)!'
2936 0 : LIBPAW_BUG(msg)
2937 : end if
2938 1508 : ncpgr=0
2939 1508 : if (present(buffer_gr)) then
2940 0 : ncpgr=size(buffer_gr,dim=2)
2941 : end if
2942 :
2943 : !=== Unpack buffers into cprj ===
2944 1508 : ipck=0
2945 3016 : do jj=1,n2dim
2946 10556 : do iat=1,natom
2947 7540 : nn=nlmn(iat)
2948 297076 : cprj(iat,jj)%cp(:,1:nn)=buffer(:,ipck+1:ipck+nn)
2949 7540 : if (ncpgr/=0) then
2950 0 : cprj(iat,jj)%dcp(:,:,1:nn)=buffer_gr(:,:,ipck+1:ipck+nn)
2951 : end if
2952 9048 : ipck=ipck+nn
2953 : end do
2954 : end do
2955 :
2956 1508 : end subroutine pawcprj_unpack
2957 :
2958 0 : end module m_pawcprj
2959 : !!***
|