Branch data Line data Source code
1 : : !!****m* ABINIT/m_iowf
2 : : !! NAME
3 : : !! m_iowf
4 : : !!
5 : : !! FUNCTION
6 : : !! Procedures for the IO of the WFK file.
7 : : !!
8 : : !! COPYRIGHT
9 : : !! Copyright (C) 1998-2026 ABINIT group (DCA, XG, GMR, AR, MB, MVer, MG)
10 : : !! This file is distributed under the terms of the
11 : : !! GNU General Public License, see ~abinit/COPYING
12 : : !! or http://www.gnu.org/copyleft/gpl.txt .
13 : : !!
14 : : !! SOURCE
15 : :
16 : : #if defined HAVE_CONFIG_H
17 : : #include "config.h"
18 : : #endif
19 : :
20 : : #include "abi_common.h"
21 : :
22 : : MODULE m_iowf
23 : :
24 : : use defs_basis
25 : : use defs_wvltypes
26 : : use m_abicore
27 : : use m_errors
28 : : use m_dtset
29 : : use m_xmpi
30 : : use m_wffile
31 : : use m_abi_etsf
32 : : use m_nctk
33 : : use m_wfk
34 : : use netcdf
35 : : use m_hdr
36 : : use m_ebands
37 : :
38 : : use m_time, only : cwtime, cwtime_report, timab
39 : : use m_io_tools, only : get_unit, flush_unit, iomode2str
40 : : use m_fstrings, only : endswith, sjoin
41 : : use m_numeric_tools, only : mask2blocks
42 : : use defs_abitypes, only : MPI_type
43 : : use defs_datatypes, only : pseudopotential_type
44 : : use m_cgtools, only : cg_zcopy
45 : : use m_crystal, only : crystal_t
46 : : use m_rwwf, only : rwwf
47 : : use m_mpinfo, only : proc_distrb_cycle
48 : : use m_vkbr, only : calc_vkb
49 : : use m_wvl_rwwf, only : wvl_write
50 : :
51 : : implicit none
52 : :
53 : : private
54 : :
55 : : public :: outwf
56 : : public :: outresid
57 : : public :: prtkbff ! Write KB form factors to WFK in netcdf format.
58 : :
59 : : !!***
60 : :
61 : : CONTAINS !====================================================================================================
62 : : !!***
63 : :
64 : : !!****f* m_iowf/outresid
65 : : !! NAME
66 : : !! outresid
67 : : !!
68 : : !! FUNCTION
69 : : !! - Compute the maximal residual and eventually print it
70 : : !!
71 : : !! INPUTS
72 : : !! dtset <type(dataset_type)>=all input variables for this dataset
73 : : !! kptns(3,nkpt)=k points in terms of recip primitive translations
74 : : !! mband=maximum number of bands
75 : : !! nband=number of bands
76 : : !! nkpt=number of k points
77 : : !! nsppol=1 for unpolarized, 2 for spin-polarized
78 : : !! resid(mband*nkpt*nsppol)=squared residuals for each band and k point
79 : : !! where resid(n,k)=|<C(n,k)|(H-e(n,k))|C(n,k)>|^2 for the ground state
80 : : !!
81 : : !! OUTPUT
82 : : !! (only writing)
83 : : !!
84 : : !! SOURCE
85 : :
86 : 8037 : subroutine outresid(dtset,kptns,mband,nband,nkpt,nsppol,resid)
87 : :
88 : : !Arguments ------------------------------------
89 : : !scalars
90 : : integer,intent(in) :: mband,nkpt,nsppol
91 : : type(dataset_type),intent(in) :: dtset
92 : : !arrays
93 : : integer, intent(in) :: nband(nkpt*nsppol)
94 : : real(dp), intent(in) :: kptns(3,nkpt)
95 : : real(dp), intent(in) :: resid(mband*nkpt*nsppol)
96 : :
97 : : !Local variables-------------------------------
98 : : integer,parameter :: nkpt_max=50
99 : : integer :: band_index,spin,ikpt,ibdkpt,nkpt_eff,nband_k,nband_eff
100 : : integer :: iband, ii
101 : : real(dp) :: resims,residm, residk
102 : : character(len=500) :: msg
103 : :
104 : : !Compute mean square and maximum residual over all bands and k points and spins
105 : : !(disregard k point weights and occupation numbers here)
106 : :
107 : : !Find largest residual over bands, k points, and spins, except for nbdbuf highest bands
108 : : !Already AVAILABLE in hdr ?!
109 : 8037 : ibdkpt=0
110 : 8037 : residm=zero
111 : 8037 : resims=zero
112 : 8037 : band_index=0
113 [ + + ]: 16785 : do spin=1,nsppol
114 [ + + ]: 240079 : do ikpt=1,nkpt
115 : 223294 : nband_k=nband(ikpt+(spin-1)*nkpt)
116 [ + + ]: 223294 : if (dtset%nbdbuf>0) then
117 : 71107 : nband_eff=max(1,nband_k-dtset%nbdbuf)
118 : : else
119 : : nband_eff=nband_k
120 : : end if
121 [ + - - + : 1586239 : residm=max(residm,maxval(resid(ibdkpt+1:ibdkpt+nband_eff)))
- - + + ]
122 [ + + ]: 1362945 : resims=resims + sum(resid(ibdkpt+1:ibdkpt+nband_eff))
123 : 223294 : ibdkpt=ibdkpt+nband_k
124 : 232042 : band_index=band_index + nband_eff
125 : : end do
126 : : end do
127 : 8037 : resims=resims/dble(band_index)
128 : :
129 : 8037 : write(msg,'(a,2p,e12.4,a,e12.4)')' Mean square residual over all n,k,spin= ',resims,'; max=',residm
130 [ + + ]: 24111 : call wrtout([std_out, ab_out], msg)
131 : :
132 : 8037 : band_index=0
133 : 8037 : nkpt_eff=nkpt
134 [ + + ]: 8037 : if( (dtset%prtvol==0 .or. dtset%prtvol==1) .and. nkpt_eff>nkpt_max ) nkpt_eff=nkpt_max
135 : :
136 : : !Loop over spin again
137 [ + + ]: 16785 : do spin=1,nsppol
138 : : ! Give (squared) residuals for all bands at each k
139 [ + + ]: 240079 : do ikpt=1,nkpt
140 : 223294 : nband_k=nband(ikpt+(spin-1)*nkpt)
141 : : ! Will not print all residuals when prtvol=0 or 1
142 [ + + ]: 223294 : if(ikpt<=nkpt_eff)then
143 : : ! Find largest residual over all bands for given k point
144 [ + - - + : 1071750 : residk=maxval(resid(1+band_index:nband_k+band_index))
- - + + ]
145 : 128781 : write(msg,'(1x,3f8.4,3x,i2,1p,e13.5,a)')kptns(1:3,ikpt),spin,residk,' kpt; spin; max resid(k); each band:'
146 [ + + ]: 128781 : if(dtset%prtvol>=2) call wrtout(ab_out, msg)
147 : 128781 : call wrtout(std_out, msg)
148 [ + + ]: 297654 : do ii=0,(nband_k-1)/8
149 [ + - + + ]: 983061 : write(msg,'(1x,1p,8e9.2)')(resid(iband+band_index),iband=1+ii*8,min(nband_k,8+ii*8))
150 [ + + ]: 168873 : if(dtset%prtvol>=2) call wrtout(ab_out, msg)
151 : 297654 : call wrtout(std_out, msg)
152 : : end do
153 [ + + ]: 94513 : else if(ikpt==nkpt_eff+1)then
154 : 759 : write(msg,'(2a)')' outresid : prtvol=0 or 1, do not print more k-points.',ch10
155 [ - + ]: 759 : if(dtset%prtvol>=2) call wrtout(ab_out, msg)
156 : 759 : call wrtout(std_out, msg)
157 : : end if
158 : 232042 : band_index=band_index+nband_k
159 : : end do
160 : : end do
161 : :
162 : 8037 : end subroutine outresid
163 : : !!***
164 : :
165 : : !!****f* m_iowf/outwf
166 : : !! NAME
167 : : !! outwf
168 : : !!
169 : : !! FUNCTION
170 : : !! Conduct output of a "wave-functions" file.
171 : : !! - Compute the maximal residual
172 : : !! - Then open a permanent file wff2 for final output of wf data
173 : : !! - Create a new header for the file.
174 : : !! - Write wave-functions (and energies)
175 : : !!
176 : : !! INPUTS
177 : : !! cg(2,mcg)=wavefunction array (storage if nkpt>1)
178 : : !! dtset <type(dataset_type)>=all input variables for this dataset
179 : : !! psps <type(pseudopotential_type)>=variables related to pseudopotentials
180 : : !! eigen( (2*mband)**response *mband*nkpt*nsppol)= eigenvalues (hartree) for all bands at each k point
181 : : !! filnam= character string giving the root to form the name of the
182 : : !! output WFK or WFQ file if response==0, otherwise it is the filename.
183 : : !! hdr <type(hdr_type)>=the header of wf, den and pot files
184 : : !! kg(3,mpw*mkmem)=reduced planewave coordinates.
185 : : !! kptns(3,nkpt)=k points in terms of recip primitive translations
186 : : !! mband=maximum number of bands
187 : : !! mcg=size of wave-functions array (cg) =mpw*nspinor*mband*mkmem*nsppol
188 : : !! mkmem=Number of k-points treated by this node.
189 : : !! mpi_enreg=information about MPI parallelization
190 : : !! mpw=maximum number of plane waves
191 : : !! natom=number of atoms in unit cell
192 : : !! nband=number of bands
193 : : !! nkpt=number of k points
194 : : !! npwarr(nkpt)=number of planewaves in basis and on boundary for each k
195 : : !! nsppol=1 for unpolarized, 2 for spin-polarized
196 : : !! nstep=desired number of electron iteration steps
197 : : !! occ(mband*nkpt*nsppol)=occupations for all bands at each k point
198 : : !! response: if == 0, GS wavefunctions , if == 1, RF wavefunctions
199 : : !! unwff2=unit for output of wavefunction
200 : : !! wfs <type(wvl_projector_type)>=wavefunctions information for wavelets.
201 : : !!
202 : : !! OUTPUT
203 : : !! (only writing)
204 : : !!
205 : : !! NOTES
206 : : !! * The name of the file wff2 might be the same as that of the file wff1.
207 : : !!
208 : : !! SOURCE
209 : :
210 : 5730 : subroutine outwf(cg,dtset,psps,eigen,filnam,hdr,kg,kptns,mband,mcg,mkmem,&
211 : 5730 : mpi_enreg,mpw,natom,nband,nkpt,npwarr,&
212 : 5730 : nsppol,occ,response,unwff2,&
213 : : wfs,wvl,force_write) ! optional
214 : :
215 : : !Arguments ------------------------------------
216 : : !scalars
217 : : integer,intent(in) :: mband,mcg,mkmem,mpw,natom,nkpt,nsppol,response,unwff2
218 : : !integer,intent(in) :: nstep
219 : : character(len=*),intent(in) :: filnam
220 : : type(MPI_type),intent(in) :: mpi_enreg
221 : : type(dataset_type),intent(in) :: dtset
222 : : type(pseudopotential_type),intent(in) :: psps
223 : : type(hdr_type), intent(inout) :: hdr
224 : : type(wvl_wf_type),intent(in) :: wfs
225 : : type(wvl_internal_type), intent(in) :: wvl
226 : : logical, intent(in), optional :: force_write
227 : : !arrays
228 : : integer, intent(in) :: kg(3,mpw*mkmem),nband(nkpt*nsppol),npwarr(nkpt)
229 : : real(dp), intent(inout) :: cg(2,mcg)
230 : : real(dp), intent(in) :: eigen((2*mband)**response*mband*nkpt*nsppol),kptns(3,nkpt)
231 : : real(dp), intent(in) :: occ(mband*nkpt*nsppol)
232 : :
233 : : !Local variables-------------------------------
234 : : integer :: iomode,action,band_index,fform,formeig,iband,icg !,iat,iproj
235 : : integer :: ierr,ikg,ikpt,spin,master,mcg_disk,me,me0,mtag,my_nspinor
236 : : integer :: nband_k,nmaster,npw_k,option,rdwr,sender,source !npwtot_k,
237 : : integer :: spaceComm,spaceComm_io,spacecomsender,spaceWorld,sread,sskip,tim_rwwf,xfdim2
238 : : #ifdef HAVE_MPI
239 : : integer :: ipwnbd
240 : : #endif
241 : : real(dp) :: cpu,wall,gflops
242 : : logical :: ihave_data,iwrite,iam_master,done,prtwf
243 : : character(len=500) :: msg
244 : : type(wffile_type) :: wff2
245 : : !character(len=fnlen) :: path
246 : : !arrays
247 : 5730 : integer,allocatable :: kg_disk(:,:)
248 : : real(dp) :: tsec(2)
249 : 5730 : real(dp),allocatable :: cg_disk(:,:),eig_k(:),occ_k(:)
250 : :
251 : : ! *************************************************************************
252 : : !For readability of the source file, define a "me" variable also in the sequential case
253 : :
254 : : DBG_ENTER("COLL")
255 : :
256 : : ABI_UNUSED(kptns(1,1))
257 : :
258 : 5730 : xfdim2 = natom+4
259 : : !Init mpi_comm
260 : 5730 : spaceWorld= mpi_enreg%comm_cell
261 : 5730 : spaceComm=spaceWorld
262 : 5730 : spaceComm_io=xmpi_comm_self
263 : :
264 [ + + ]: 5730 : if (mpi_enreg%paral_kgb==1 ) spaceComm_io= mpi_enreg%comm_bandspinorfft
265 : : if (mpi_enreg%paral_kgb==1 ) spaceComm= mpi_enreg%comm_cell
266 : :
267 : : !Paral_kgb=1 and Fortran-I/O is not supported (only for testing purpose)
268 [ + + + + ]: 5730 : if (mpi_enreg%paral_kgb==1.and.dtset%iomode==IO_MODE_FORTRAN) then
269 : 132 : spaceWorld=mpi_enreg%comm_kpt
270 : : write(msg,'(7a)') &
271 : 132 : 'WF file is written using standard Fortran I/O',ch10,&
272 : 132 : 'and Kpt-band-FFT parallelization is active !',ch10,&
273 : 132 : 'This is only allowed for testing purposes.',ch10,&
274 : 264 : 'The produced WF file will be incomplete and not useable.'
275 : 132 : ABI_WARNING(msg)
276 : : end if
277 : :
278 : : !If parallel HF calculation
279 [ - + ]: 5730 : if (mpi_enreg%paral_hf==1 ) spaceComm_io= mpi_enreg%comm_hf
280 [ - + ]: 5730 : if (mpi_enreg%paral_hf==1 ) spaceComm= mpi_enreg%comm_cell
281 : :
282 : : !Paral_hf=1 and Fortran-I/O is not supported (copy from paral_kgb... not tested)
283 [ - + - - ]: 5730 : if (mpi_enreg%paral_hf==1.and.dtset%iomode==IO_MODE_FORTRAN) then
284 : 0 : spaceWorld=mpi_enreg%comm_kpt
285 : : write(msg,'(7a)') &
286 : 0 : 'WF file is written using standard Fortran I/O',ch10,&
287 : 0 : 'and HF parallelization is active !',ch10,&
288 : 0 : 'This is only allowed for testing purposes.',ch10,&
289 : 0 : 'The produced WF file will be incomplete and not useable.'
290 : 0 : ABI_WARNING(msg)
291 : : end if
292 : :
293 : : ! check consistency between dimensions and input hdr.
294 : : !ABI_CHECK(mband == maxval(hdr%nband), "hdr:mband")
295 : : !ABI_CHECK(nkpt == hdr%nkpt, "hdr:nkpt")
296 : : !ABI_CHECK(nsppol == hdr%nsppol, "hdr:nsppol")
297 : : !ABI_CHECK(all(hdr%npwarr == npwarr), "hdr:npwarr")
298 : : !ABI_CHECK(all(hdr%nband == nband), "hdr:nband")
299 : : !ABI_CHECK(maxval(hdr%npwarr) == mpw, "hdr:nband")
300 : :
301 : : !Init me
302 : 5730 : me=mpi_enreg%me_kpt
303 : 5730 : me0=me
304 : : !Define master
305 : 5730 : master=0
306 : :
307 : 5730 : my_nspinor=max(1,dtset%nspinor/mpi_enreg%nproc_spinor)
308 : 5730 : tim_rwwf =0
309 : 5730 : source = master
310 : 5730 : sread = master
311 : 5730 : iam_master=(master==me)
312 : 5730 : iwrite=iam_master
313 : 5730 : sender=-1
314 : :
315 : : !Will write the wavefunction file only when nstep>0
316 : : !MT 07 2015: writing reactivated when nstep=0
317 : : !if (nstep>0 .and. dtset%prtwf/=0) then
318 : : !FB 03/2022: Added an option to force writing (used in RT-TDDFT)
319 : 5730 : prtwf = dtset%prtwf/=0
320 [ + + ]: 5730 : if (present(force_write)) then
321 [ - + ]: 49 : if (force_write) prtwf = .true.
322 : : end if
323 [ + + ]: 5681 : if (prtwf) then
324 : :
325 : : ! Only the master write the file, except if MPI I/O, but the
326 : : ! full wff dataset should be provided to WffOpen in this case
327 : 3238 : iomode=IO_MODE_FORTRAN_MASTER
328 [ + + ]: 3238 : if (dtset%iomode==IO_MODE_MPI) iomode = IO_MODE_MPI
329 [ + + ]: 3238 : if (dtset%iomode==IO_MODE_ETSF) iomode = IO_MODE_ETSF
330 : :
331 [ + + + - ]: 3238 : if (dtset%iomode == IO_MODE_ETSF .and. dtset%usewvl == 0) then
332 : : call cg_ncwrite(filnam,hdr,dtset,response,mpw,mband,nband,nkpt,nsppol,&
333 : 3078 : dtset%nspinor,mcg,mkmem,eigen,occ,cg,npwarr,kg,mpi_enreg,done)
334 : : #ifdef HAVE_NETCDF_DEFAULT
335 [ - + ]: 3078 : ABI_CHECK(done, "cg_ncwrite must handle the output of the WFK file.")
336 : : #endif
337 : :
338 : : ! Write KB form factors. Only master works. G-vectors are read from file to avoid
339 : : ! having to deal with paral_kgb distribution.
340 [ + + + + : 3078 : if (me == master .and. dtset%prtkbff == 1 .and. dtset%iomode == IO_MODE_ETSF .and. dtset%usepaw == 0) then
+ - ]
341 [ - + ]: 1 : ABI_CHECK(done, "cg_ncwrite was not able to generate WFK.nc in parallel. Perhaps hdf5 is not working")
342 : 1 : call prtkbff(filnam, hdr, psps, dtset%prtvol)
343 : : end if
344 : :
345 [ + - ]: 3078 : if (done) return
346 : : ! If cg_ncwrite cannot handle the IO because HDF5 + MPI-IO support is missing, we fallback to Fortran + MPI-IO.
347 : 0 : msg = "Could not produce a netcdf file in parallel (MPI-IO support is missing). Will fallback to MPI-IO with Fortran"
348 : 0 : ABI_WARNING(msg)
349 : 0 : iomode=IO_MODE_MPI
350 : : end if
351 : :
352 : 160 : call cwtime(cpu, wall, gflops, "start")
353 [ + - ]: 160 : call wrtout(std_out, sjoin(ch10,' outwf: writing wavefunctions to:', trim(filnam), "with iomode:", iomode2str(iomode)))
354 : :
355 : : ! Create an ETSF file for the wavefunctions
356 [ - + ]: 160 : if (iomode == IO_MODE_ETSF) then
357 [ # # ]: 0 : ABI_CHECK(xmpi_comm_size(spaceComm) == 1, "Legacy etsf-io code does not support nprocs > 1")
358 : 0 : ABI_ERROR("ETSF_IO is not activated")
359 : : ABI_UNUSED(psps%ntypat)
360 : : end if
361 : :
362 : 160 : call WffOpen(iomode,spaceComm,filnam,ierr,wff2,master,me0,unwff2,spaceComm_io)
363 : : ! Conduct wavefunction output to wff2
364 : :
365 [ - + + - : 480 : ABI_MALLOC(kg_disk,(3,mpw))
+ - + - -
+ - + ]
366 : :
367 : 160 : mcg_disk=mpw*my_nspinor*mband
368 [ - + ]: 160 : formeig=0; if (response==1) formeig=1
369 : :
370 [ - + + - : 480 : ABI_MALLOC(eig_k,( (2*mband)**formeig * mband))
+ - + - -
+ - + ]
371 [ - + + - : 480 : ABI_MALLOC(occ_k,(mband))
+ - + - -
+ - + ]
372 : :
373 : : #ifdef HAVE_MPI
374 : 160 : call xmpi_barrier(spaceComm)
375 : : ! Compute mband and mpw
376 [ - + + - : 480 : ABI_MALLOC_OR_DIE(cg_disk,(2,mcg_disk), ierr)
+ - + - +
- - + -
+ ]
377 : : #endif
378 : :
379 : 160 : band_index=0
380 : 160 : icg=0
381 [ + - ]: 160 : if(mpi_enreg%paralbd==0) tim_rwwf=6
382 [ - + ]: 160 : if(mpi_enreg%paralbd==1) tim_rwwf=12
383 : :
384 : : ! Write header info for new wf file
385 : 160 : rdwr=2
386 [ + - ]: 160 : if (dtset%usewvl==0) then
387 : 160 : fform=2
388 : : else
389 : 0 : fform = 200 ! Use 200 as radical for naming file format used by wavelets.
390 : : end if
391 : :
392 [ + - ]: 160 : if (wff2%iomode < 2) then
393 : 160 : call hdr_io(fform,hdr,rdwr,wff2)
394 : 160 : call WffKg(wff2,1)
395 [ # # # # ]: 0 : else if (wff2%iomode==IO_MODE_ETSF .and. iam_master) then
396 [ # # ]: 0 : NCF_CHECK(hdr%ncwrite(wff2%unwff, fform, nc_define=.True.))
397 : : end if
398 : :
399 [ + + ]: 340 : do spin=1,nsppol
400 : 180 : ikg=0
401 : :
402 [ + + ]: 932 : do ikpt=1,nkpt
403 : 592 : nband_k=nband(ikpt+(spin-1)*nkpt)
404 : 592 : npw_k=npwarr(ikpt)
405 : :
406 : : #ifdef HAVE_MPI
407 [ - + ]: 592 : if (dtset%usewvl == 0) then
408 : 592 : mtag=ikpt+(spin-1)*nkpt
409 : 592 : call xmpi_barrier(spaceWorld)
410 : :
411 : : ! Must transfer the wavefunctions to the master processor
412 : : ! Separate sections for paralbd=1 or other values ; might be merged
413 [ + - ]: 592 : if(mpi_enreg%paralbd==0)then
414 : 592 : nmaster=0
415 [ + + ]: 9648 : source=minval(mpi_enreg%proc_distrb(ikpt,1:nband_k,spin))
416 : 592 : ihave_data=.false.
417 [ + + ]: 592 : if(source==me)ihave_data=.true.
418 : 416 : action=0
419 : : ! I am the master node, and I have the data in cg or cg_disk
420 [ + + ]: 416 : if((iam_master).and.(ihave_data))action=1
421 : : ! I am not the master, and I have the data => send to master
422 [ + + ]: 592 : if((.not.iam_master).and.(ihave_data))action=2
423 : : ! I am the master, and I receive the data
424 [ + + ]: 592 : if((iam_master).and.(.not.ihave_data))action=3
425 : :
426 : : ! I have the data in cg or cg_disk ( MPI_IO case)
427 [ + + ]: 592 : if (iomode==IO_MODE_MPI) then
428 : 60 : action = 0
429 : 60 : sender=-1
430 : 60 : iwrite=.false.
431 [ + + ]: 60 : if (ihave_data)then
432 : 52 : action=1
433 : 52 : iwrite=.true.
434 : 52 : sender=me
435 : : end if
436 : : end if
437 : :
438 : : ! I am the master node, and I have the data in cg or cg_disk
439 : : ! I have the data in cg or cg_disk ( MPI_IO case)
440 [ + + ]: 532 : if(action==1)then
441 : : ! Copy from kg to kg_disk
442 [ + + + + ]: 182300 : kg_disk(:,1:npw_k)=kg(:,1+ikg:npw_k+ikg)
443 : : ! Copy from cg to cg_disk
444 [ + + ]: 874956 : do ipwnbd=1,nband_k*npw_k*my_nspinor
445 : 874604 : cg_disk(1,ipwnbd)=cg(1,ipwnbd+icg)
446 : 874956 : cg_disk(2,ipwnbd)=cg(2,ipwnbd+icg)
447 : : end do
448 : : end if
449 : :
450 : : ! I am not the master, and I have the data => send to master
451 : : ! I am the master, and I receive the data
452 [ + + ]: 592 : if ( action==2.or.action==3) then
453 : : !write(std_out,*)npw_k,nband_k
454 : 128 : call timab(48,1,tsec)
455 [ + + ]: 128 : if(action==2)then
456 : 64 : call xmpi_exch(kg(:,1+ikg:npw_k+ikg),3*npw_k,source,kg_disk,nmaster,spaceWorld,2*mtag+1,ierr)
457 : : call xmpi_exch(cg(:,icg+1:icg+nband_k*npw_k*my_nspinor),2*nband_k*npw_k*my_nspinor, &
458 : 64 : & source,cg_disk,nmaster,spaceWorld,2*mtag+2,ierr)
459 : : else
460 : 64 : call xmpi_exch(kg_disk,3*npw_k,source,kg_disk,nmaster,spaceWorld,2*mtag+1,ierr)
461 : 64 : call xmpi_exch(cg_disk,2*nband_k*npw_k*my_nspinor,source,cg_disk,nmaster,spaceWorld,2*mtag+2,ierr)
462 : : end if
463 : 128 : call timab(48,2,tsec)
464 : : end if
465 : :
466 : :
467 [ # # ]: 0 : else if(mpi_enreg%paralbd==1)then
468 : 0 : nmaster=0
469 : : #ifdef HAVE_MPI_IO
470 : 0 : sender=IO_MODE_FORTRAN_MASTER
471 [ # # ]: 0 : if( iomode==IO_MODE_MPI) then
472 : 0 : nmaster=mpi_enreg%proc_distrb(ikpt,1,spin)
473 : 0 : sender=nmaster
474 : : end if
475 : : #endif
476 : :
477 : : ! Note the loop over bands
478 [ # # ]: 0 : do iband=1,nband_k
479 : :
480 : : ! The message passing related to kg is counted as one band
481 : 0 : action=0
482 : :
483 : : ! I am the master node, and I have the data in cg or cg_disk
484 [ # # # # ]: 0 : if( mpi_enreg%proc_distrb(ikpt,iband,spin)==nmaster .and. me==nmaster) then
485 : : action=1
486 : : ! I am not the master, and I have the data => send to master
487 [ # # # # ]: 0 : elseif( mpi_enreg%proc_distrb(ikpt,iband,spin)==me .and. me/=nmaster ) then
488 : : action = 2
489 : : ! I am the master, and I receive the data
490 [ # # # # ]: 0 : elseif( mpi_enreg%proc_distrb(ikpt,iband,spin)/=me .and. me==nmaster ) then
491 : : action=3
492 : : end if
493 : :
494 : : if(action==1) then
495 : : ! I am the master node, and I have the data in cg or cg_disk
496 : : ! Copy from kg to kg_disk
497 [ # # # # : 0 : if(iband==1)kg_disk(:,1:npw_k)=kg(:,1+ikg:npw_k+ikg)
# # ]
498 : : ! Copy from cg to cg_disk
499 [ # # ]: 0 : do ipwnbd=1,npw_k*my_nspinor
500 : 0 : cg_disk(1,(iband-1)*npw_k*my_nspinor+ipwnbd) = cg(1,(iband-1)*npw_k*my_nspinor+ipwnbd+icg)
501 : 0 : cg_disk(2,(iband-1)*npw_k*my_nspinor+ipwnbd) = cg(2,(iband-1)*npw_k*my_nspinor+ipwnbd+icg)
502 : : end do
503 : : end if ! action=1
504 : :
505 [ # # ]: 0 : if ( action==2.or.action==3) then
506 : : ! action=2 : I am not the master, and I have the data => send to master
507 : : ! action=3 : I am the master, and I receive the data
508 : 0 : call timab(48,1,tsec)
509 [ # # ]: 0 : if ( iband == 1 ) then
510 [ # # ]: 0 : if (action==2) then
511 : : call xmpi_exch(kg(:,1+ikg:npw_k+ikg),3*npw_k,mpi_enreg%proc_distrb(ikpt,iband,spin), &
512 : 0 : & kg_disk,nmaster,spaceWorld,iband*(mtag-1)+1,ierr)
513 : : else
514 : : call xmpi_exch(kg_disk,3*npw_k,mpi_enreg%proc_distrb(ikpt,iband,spin), &
515 : 0 : & kg_disk,nmaster,spaceWorld,iband*(mtag-1)+1,ierr)
516 : : end if
517 : : end if ! iband =1
518 : 0 : ipwnbd=(iband-1)*npw_k*my_nspinor
519 [ # # ]: 0 : if (action==2) then
520 : : call xmpi_exch( cg(:,ipwnbd+icg+1:ipwnbd+icg+npw_k*my_nspinor),2*npw_k*my_nspinor &
521 : : & ,mpi_enreg%proc_distrb(ikpt,iband,spin) &
522 : 0 : & ,cg_disk(:,ipwnbd+1:ipwnbd+npw_k*my_nspinor),nmaster,spaceWorld,iband*(mtag-1)+2,ierr)
523 : : else
524 : : call xmpi_exch( cg_disk(:,ipwnbd+1:ipwnbd+npw_k*my_nspinor),2*npw_k*my_nspinor &
525 : : & ,mpi_enreg%proc_distrb(ikpt,iband,spin) &
526 : 0 : & ,cg_disk(:,ipwnbd+1:ipwnbd+npw_k*my_nspinor),nmaster,spaceWorld,iband*(mtag-1)+2,ierr)
527 : : end if
528 : :
529 : 0 : call timab(48,2,tsec)
530 : : end if ! action=2 or action=3
531 : :
532 [ # # ]: 0 : if(iomode==IO_MODE_MPI) then
533 : : ! I have the data in cg or cg_disk
534 : 0 : iwrite=.false.
535 [ # # ]: 0 : if (nmaster == me) iwrite=.true.
536 : : end if
537 : :
538 : : end do ! End of loop over bands
539 : : end if ! End of paralbd=1
540 : : end if
541 : : #endif
542 : :
543 : : ! Only the master will write to disk the final output wf file.
544 : : ! in MPI_IO case only iwrite will write to disk the final output wf file.
545 [ + + ]: 592 : if(iwrite) then
546 : : ! write(std_out,*) 'outwf : I am master and will write wf file'
547 [ + - ]: 416 : if(formeig==0)then
548 [ + + ]: 7056 : eig_k(1:nband_k)=eigen(1+band_index:nband_k+band_index)
549 [ + + ]: 7056 : occ_k(1:nband_k)=occ(1+band_index:nband_k+band_index)
550 : : else
551 [ # # ]: 0 : eig_k(1:2*nband_k*nband_k)=eigen(1+band_index:2*nband_k*nband_k+band_index)
552 : : end if
553 : 416 : option=2
554 [ - + ]: 416 : if(dtset%prtwf==3)option=5
555 : : ! if (dtset%prtwf == 2 .and. mkmem/=0) option=4
556 : :
557 [ + - ]: 416 : if (dtset%usewvl == 0) then
558 : : #ifdef HAVE_MPI
559 : : call rwwf(cg_disk,eig_k,formeig,0,0,ikpt,spin,kg_disk,mband,mcg_disk,mpi_enreg, &
560 : 416 : & nband_k, nband_k,npw_k,my_nspinor,occ_k,option,1,tim_rwwf,wff2)
561 : :
562 : : #else
563 : : kg_disk(:,1:npw_k)=kg(:,1+ikg:npw_k+ikg)
564 : : call rwwf(cg,eig_k,formeig,0,icg,ikpt,spin,kg_disk,mband,mcg,mpi_enreg,nband_k, &
565 : : & nband_k, npw_k,my_nspinor,occ_k,option,1,tim_rwwf,wff2)
566 : : #endif
567 : : else
568 : 0 : call wvl_write(dtset,eigen,mpi_enreg,option,hdr%rprimd,wff2,wfs,wvl,hdr%xred)
569 : : end if
570 : : end if
571 : :
572 : : ! The wavefunctions for the present k point and spin are written
573 [ + - ]: 592 : if(response==0)band_index=band_index+nband_k
574 [ - + ]: 592 : if(response==1)band_index=band_index+2*nband_k*nband_k
575 : :
576 : 592 : sskip=1
577 : : #ifdef HAVE_MPI
578 [ + - ]: 592 : if (dtset%usewvl == 0) then
579 : 592 : sskip=0
580 [ + + ]: 592 : if(.not.(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,spin,me)))sskip=1
581 : : end if
582 : : #endif
583 : : if(sskip==1)then
584 : 416 : icg=icg+npw_k*my_nspinor*nband_k
585 : 416 : ikg=ikg+npw_k
586 : : end if
587 : :
588 : :
589 : : #ifdef HAVE_MPI_IO
590 : 592 : spacecomsender=spaceComm
591 [ + + ]: 592 : if (mpi_enreg%paral_kgb==1) spacecomsender =mpi_enreg%comm_kpt
592 [ - + ]: 592 : if (mpi_enreg%paral_hf==1) spacecomsender =mpi_enreg%comm_kpt
593 : 772 : call WffOffset(wff2,sender,spacecomsender,ierr)
594 : : #endif
595 : :
596 : : end do ! ikpt
597 : : end do ! spin
598 : :
599 [ - + ]: 160 : ABI_FREE(kg_disk)
600 : : #ifdef HAVE_MPI
601 [ - + ]: 160 : ABI_FREE(cg_disk)
602 : : #endif
603 : :
604 [ - + ]: 160 : ABI_FREE(eig_k)
605 [ - + ]: 160 : ABI_FREE(occ_k)
606 : :
607 : : ! Close the wavefunction file (and do NOT delete it !)
608 : : !if (wff2%iomode /= IO_MODE_NETCDF) then
609 : 160 : call WffClose(wff2,ierr)
610 : : !end if
611 : :
612 : 480 : call cwtime_report(" WFK output", cpu, wall, gflops)
613 : : end if ! End condition of nstep>0
614 : :
615 : : ! Block here because we might need to read the WFK file in the caller.
616 : 2652 : call xmpi_barrier(mpi_enreg%comm_cell)
617 : :
618 : : DBG_EXIT("COLL")
619 : :
620 [ - + - + ]: 2652 : end subroutine outwf
621 : : !!***
622 : :
623 : : !----------------------------------------------------------------------
624 : :
625 : : !!****f* m_iowf/cg_ncwrite
626 : : !! NAME
627 : : !! cg_ncwrite
628 : : !!
629 : : !! FUNCTION
630 : : !! Conduct output of a "wave-functions" file with netcdf
631 : : !!
632 : : !! INPUTS
633 : : !! fname= character string giving the root to form the name of the
634 : : !! output WFK or WFQ file if response==0, otherwise it is the filename.
635 : : !! dtset <type(dataset_type)>=all input variables for this dataset
636 : : !! hdr <type(hdr_type)>=the header of wf, den and pot files
637 : : !! response: if == 0, GS wavefunctions , if == 1, RF wavefunctions
638 : : !! mpw=maximum number of plane waves
639 : : !! mband=maximum number of bands
640 : : !! mcg=size of wave-functions array (cg) =mpw*nspinor*mband*mkmem*nsppol
641 : : !! mkmem=maximum number of k-points treated by this node
642 : : !! nkpt=number of k points
643 : : !! nsppol=1 for unpolarized, 2 for spin-polarized
644 : : !! eigen((2*mband)**response *mband*nkpt*nsppol)= eigenvalues (hartree) for all bands at each k point
645 : : !! occ(mband*nkpt*nsppol)=occupations for all bands at each k point
646 : : !! cg(2,mcg)=wavefunction array (storage if nkpt>1)
647 : : !! npwarr(nkpt)=number of planewaves in basis and on boundary for each k
648 : : !! kg(3,mpw*mkmem)=reduced planewave coordinates.
649 : : !! mpi_enreg=information about MPI parallelization
650 : : !!
651 : : !! OUTPUT
652 : : !! done=.True if cg_ncwrite can handle the output of the WFK file in parallel.
653 : : !!
654 : : !! SOURCE
655 : :
656 : 3078 : subroutine cg_ncwrite(fname,hdr,dtset,response,mpw,mband,nband,nkpt,nsppol,nspinor,mcg,&
657 : 3078 : mkmem,eigen,occ,cg,npwarr,kg,mpi_enreg,done)
658 : :
659 : : !Arguments ------------------------------------
660 : : !scalars
661 : : integer,intent(in) :: response,mband,mcg,mkmem,mpw,nkpt,nsppol,nspinor
662 : : character(len=*),intent(in) :: fname
663 : : logical,intent(out) :: done
664 : : type(dataset_type),intent(in) :: dtset
665 : : type(MPI_type),intent(in) :: mpi_enreg
666 : : type(hdr_type),intent(in) :: hdr
667 : : !arrays
668 : : integer, intent(in) :: nband(nkpt*nsppol),kg(3,mpw*mkmem),npwarr(nkpt)
669 : : real(dp),intent(in) :: cg(2,mcg)
670 : : real(dp),intent(in) :: eigen((2*mband)**response*mband*nkpt*nsppol),occ(mband*nkpt*nsppol)
671 : :
672 : : !Local variables-------------------------------
673 : : !scalars
674 : : integer,parameter :: master=0,fform2=2
675 : : integer :: ii,iomode,icg,iband,ikg,ikpt,spin,me_cell,me_kpt,me_band,me_spinor,my_nspinor,nband_k,npw_k
676 : : integer :: comm_cell,comm_fft,comm_bandfft,mtag,formeig
677 : : integer :: cnt,min_cnt,max_cnt,ierr,action,source,ncid,ncerr,cg_varid,kg_varid !,eig_varid,
678 : : integer :: paral_kgb,npwtot_k !,start_pwblock !,start_cgblock !count_pwblock,
679 : : integer :: ipw,ispinor_index,npwso,npwsotot,npwtot,nspinortot,ikpt_this_proc,ispinor
680 : : integer :: bandpp,nproc_band,nproc_fft,nproc_spinor,me_fft,nproc_cell,nwrites
681 : : integer :: comm_mpiio,nranks,bstart,bcount !nbdblock,nblocks,
682 : : !integer :: band_blocksize,band
683 : : real(dp) :: cpu,wall,gflops
684 : : logical :: ihave_data,iam_master,single_writer,same_layout,use_collective
685 : : character(len=500) :: msg
686 : : character(len=fnlen) :: path
687 [ + + ]: 12312 : type(wfk_t) :: wfk
688 : : type(crystal_t) :: crystal
689 : 3078 : type(ebands_t) :: gs_ebands
690 : : !arrays
691 : 3078 : integer,allocatable :: kg_k(:,:),iter2kscgkg(:,:),ind_cg_mpi_to_seq(:),rank_has_cg(:),ranks_io(:)!,gblock(:,:)
692 : : real(dp) :: tsec(2)
693 : 3078 : real(dp),allocatable :: eigen3d(:,:,:),occ3d(:,:,:),cg_k(:,:),my_cgblock(:,:,:)
694 : :
695 : : ! *************************************************************************
696 : :
697 : : DBG_ENTER("COLL")
698 : 3078 : done = .False.
699 : :
700 : 3078 : path = nctk_ncify(fname)
701 : 3078 : call wrtout(std_out, sjoin(" In cg_ncwrite with path:", path))
702 : :
703 : : ! communicators and ranks
704 : 3078 : comm_cell = mpi_enreg%comm_cell; me_cell = mpi_enreg%me_cell; nproc_cell = mpi_enreg%nproc_cell
705 : 3078 : comm_fft = mpi_enreg%comm_fft; me_fft = mpi_enreg%me_fft; nproc_fft = mpi_enreg%nproc_fft
706 : 3078 : comm_bandfft = mpi_enreg%comm_bandfft
707 : 3078 : me_kpt = mpi_enreg%me_kpt; me_band = mpi_enreg%me_band; me_spinor = mpi_enreg%me_spinor
708 : 3078 : iam_master = (me_kpt == master)
709 : :
710 : 3078 : paral_kgb = dtset%paral_kgb
711 : 3078 : nproc_band = mpi_enreg%nproc_band
712 : 3078 : bandpp = mpi_enreg%bandpp
713 : 3078 : nproc_spinor = mpi_enreg%nproc_spinor
714 : :
715 : : ! FIXME
716 : 3078 : my_nspinor = max(1, nspinor/nproc_spinor)
717 [ + + - + ]: 3078 : if (nspinor == 2 .and. my_nspinor == 1) then
718 : 0 : ABI_ERROR("Spinor parallelization not coded yet")
719 : : end if
720 : :
721 : : ! TODO: Be careful with response == 1 in parallel because the distribution of the cg
722 : : ! can be **VERY** different from cg if nprocs > nkpt * nsppol
723 : : !ABI_CHECK(response==0, "response == 1 not coded")
724 [ + - ]: 3078 : if (size(hdr%nband) == size(nband)) then
725 [ + + + - : 44915 : ABI_CHECK(all(Hdr%nband == nband),"nband")
- + ]
726 : : else
727 : 0 : ABI_ERROR("hdr%nband and nband have different size!")
728 : : end if
729 : :
730 [ + + ]: 3078 : if (xmpi_comm_size(comm_cell) == 1) then
731 [ + + + - : 32799 : ABI_CHECK(all(npwarr == hdr%npwarr), "npwarr != hdr%npwarr")
- + ]
732 : : end if
733 : :
734 : 3078 : crystal = hdr%get_crystal()
735 : :
736 : : ! TODO
737 : : ! Be careful with response == 1.
738 : : ! gs_ebands contains the GS eigenvalues and occupation and will be written if this is a
739 : : ! GS wfk. If we have a DFPT file, eigen stores the GKK matrix element, in this case
740 : : ! we don't write gs_ebands but we define new variables in the netcdf file to store the GKK
741 [ - + + - : 15390 : ABI_MALLOC(occ3d, (mband,nkpt,nsppol))
- + + - -
+ + - + -
+ - - + -
+ ]
742 : 3078 : call unpack_eneocc(nkpt,nsppol,mband,nband,occ,occ3d)
743 : :
744 [ + + ]: 3078 : if (response == 0) then
745 : 3065 : formeig = 0
746 [ - + + - : 12260 : ABI_MALLOC(eigen3d, (mband,nkpt,nsppol))
- + + - -
+ + - - +
- + ]
747 : 3065 : call unpack_eneocc(nkpt,nsppol,mband,nband,eigen,eigen3d)
748 [ + - + - : 398657 : call gs_ebands%from_hdr(hdr, mband, eigen3d); gs_ebands%occ = occ3d
+ - - + -
- + + + +
+ + ]
749 [ - + ]: 3065 : ABI_FREE(eigen3d)
750 : : else
751 : 13 : formeig = 1
752 : : end if
753 : :
754 : : ! same_layout is set to True if the internal representation of the cgs
755 : : ! is compatible with the representation on file.
756 : 3078 : iomode = IO_MODE_ETSF
757 [ + + ]: 3078 : if (response == 0) then
758 [ + + + - : 5914 : same_layout = (paral_kgb == 0 .or. (paral_kgb == 1 .and. all([nproc_fft, nproc_band, nproc_spinor] == 1)))
+ + + - ]
759 : : else
760 : : ! For the time being, these cases are not implemented in the DFPT part.
761 [ - + ]: 13 : ABI_CHECK(nproc_fft==1, "nproc_fft != 1 not coded")
762 [ - + ]: 13 : ABI_CHECK(nproc_band==1, "nproc_band != 1 not coded")
763 [ - + ]: 13 : ABI_CHECK(nproc_spinor==1, "nproc_spinor != 1 not coded")
764 : :
765 : : ! Note: It would be possible to use collective IO if the cg1 are block-distributed
766 : 13 : same_layout = .True.
767 [ + + ]: 26 : spin_loop: do spin=1,nsppol
768 [ + + ]: 1722 : do ikpt=1,nkpt
769 : 1696 : nband_k = nband(ikpt + (spin-1)*nkpt)
770 [ + - ]: 1709 : if (.not. proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,spin,me_kpt)) then
771 [ + + + - : 8736 : if (any(mpi_enreg%proc_distrb(ikpt,:nband_k,spin) /= me_kpt)) then
- + ]
772 : 0 : same_layout = .False.; exit spin_loop
773 : : end if
774 : : end if
775 : : end do
776 : : end do spin_loop
777 : 13 : call xmpi_land(same_layout, comm_cell)
778 : : end if
779 : :
780 : 3078 : call cwtime(cpu, wall, gflops, "start")
781 : :
782 [ + + ]: 3078 : if (same_layout) then
783 : 2862 : single_writer = .True.
784 [ + - ]: 2862 : if (nctk_has_mpiio) single_writer = .False.
785 : : !single_writer = .True.
786 : :
787 : : write(msg,'(5a,l1)')&
788 [ + - ]: 2862 : " same layout --> writing WFK file: ",trim(path),", with iomode: ",trim(iomode2str(iomode)),", single writer: ",single_writer
789 : 2862 : call wrtout(std_out,msg, 'PERS', do_flush=.True.)
790 : :
791 [ + - ]: 2862 : if (.not. single_writer) then
792 : :
793 : : ! master opens the file and write the metadata.
794 [ + + ]: 2862 : if (xmpi_comm_rank(comm_cell) == master) then
795 : : ncerr = nf90_einval
796 : : #ifdef HAVE_NETCDF_MPI
797 : : ncerr = nf90_create(path, cmode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_write), &
798 : 2643 : comm=xmpi_comm_self, info=xmpio_info, ncid=ncid)
799 : : #endif
800 [ - + ]: 2643 : NCF_CHECK_MSG(ncerr, sjoin("create_par: ", path))
801 : :
802 : 2643 : call wfk_ncdef_dims_vars(ncid, hdr, fform2, write_hdr=.True.)
803 [ - + ]: 2643 : NCF_CHECK(crystal%ncwrite(ncid))
804 : :
805 [ + + ]: 2643 : if (response == 0) then
806 : : ! Write Gs bands
807 [ - + ]: 2630 : NCF_CHECK(gs_ebands%ncwrite(ncid))
808 : : else
809 : : ! Write H1 matrix elements and occupancies.
810 : : ! Note that GS eigens are not written here.
811 : 13 : call ncwrite_eigen1_occ(ncid, nband, mband, nkpt, nsppol, eigen, occ3d)
812 : : end if
813 : :
814 [ - + ]: 2643 : NCF_CHECK(nf90_close(ncid))
815 : : end if
816 : :
817 : : ! Compute the table for collective IO.
818 [ - + + - : 208691 : ABI_CALLOC(iter2kscgkg, (4, nkpt*nsppol))
+ - + - -
+ - + + +
+ + ]
819 : 2862 : cnt = 0; icg = 0
820 [ + + ]: 5998 : do spin=1,nsppol
821 : 3136 : ikg = 0
822 [ + + ]: 46019 : do ikpt=1,nkpt
823 : 40021 : nband_k = nband(ikpt + (spin-1)*nkpt)
824 : 40021 : npw_k = npwarr(ikpt)
825 [ + + ]: 43157 : if (.not. proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,spin,me_kpt)) then
826 : : ! FIXME v2[01], v3[6], with 4 procs
827 : : ! v6[7] and v7[68], v7[69] fail but due to a extra line with nprocs
828 : : ! v7[96] with np=4 seems to be more serious but it crashes also in trunk.
829 : : ! v2[88] is slow likely due to outscfv
830 [ + + + - : 322154 : ABI_CHECK(all(mpi_enreg%proc_distrb(ikpt,:nband_k,spin) == me_kpt), "bands are distributed")
- + ]
831 : 34191 : cnt = cnt + 1
832 [ + + ]: 170955 : iter2kscgkg(:,cnt) = [ikpt, spin, icg, ikg]
833 : 34191 : icg = icg + npw_k*my_nspinor*nband_k
834 : 34191 : ikg = ikg + npw_k
835 : : end if
836 : : end do
837 : : end do
838 [ + + ]: 2862 : if (cnt == 0) then
839 : 36 : write(std_out,*)"cnt == 0 for me_cell, me, me_kpt",mpi_enreg%me_cell, mpi_enreg%me, mpi_enreg%me_kpt
840 : : !ABI_CHECK(cnt > 0, "This processor does not have wavefunctions!")
841 : : end if
842 : :
843 : 2862 : call xmpi_min(cnt, min_cnt, comm_cell, ierr)
844 : 2862 : call xmpi_max(cnt, max_cnt, comm_cell, ierr)
845 : :
846 : : ! Handle idle procs, i.e. processors that do not have wavefunctions
847 : : ! This happens if paral_kgb == 0 and nprocs > nkpt * nsppol (Abinit does not stop anymore!)
848 : 2862 : comm_mpiio = comm_cell
849 : :
850 [ + + ]: 2862 : if (min_cnt <= 0) then
851 : 72 : ABI_COMMENT("Will create subcommunicator to exclude idle processors from MPI-IO collective calls")
852 [ - + ]: 72 : ABI_CHECK(paral_kgb == 0, "paral_kgb == 1 with idle processors should never happen")
853 : :
854 : : ! Prepare the call to xmpi_subcomm that will replace comm_mpiio.
855 [ - + + - : 496 : ABI_CALLOC(rank_has_cg, (0:nproc_cell-1))
+ - + - -
+ - + +
+ ]
856 [ + + ]: 160 : do spin=1,nsppol
857 [ + + ]: 258 : do ikpt=1,nkpt
858 : 134 : nband_k = nband(ikpt + (spin-1)*nkpt)
859 [ + + ]: 186 : if (.not. proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,spin,me_kpt)) then
860 : 36 : rank_has_cg(me_kpt) = 1
861 : 36 : exit
862 : : end if
863 : : end do
864 : : end do
865 : :
866 : 72 : call xmpi_sum(rank_has_cg,comm_cell,ierr)
867 [ + + + + ]: 352 : nranks = count(rank_has_cg == 1)
868 [ - + + - : 216 : ABI_MALLOC(ranks_io, (nranks))
+ - + - -
+ - + ]
869 : 72 : cnt = 0
870 [ + + ]: 352 : do ii=0,nproc_cell-1
871 [ + + ]: 352 : if (rank_has_cg(ii) == 1) then
872 : 140 : cnt = cnt + 1
873 : 140 : ranks_io(cnt) = ii
874 : : end if
875 : : end do
876 : : !write(std_out,*)"nranks, ranks_io:", nranks, ranks_io
877 : 72 : comm_mpiio = xmpi_subcomm(comm_cell, nranks, ranks_io)
878 [ + + ]: 72 : if (.not. rank_has_cg(me_kpt) == 1) then
879 : 36 : comm_mpiio = xmpi_comm_null
880 [ - + ]: 36 : ABI_CHECK(rank_has_cg(me_kpt) == 0, "rank_has_cg must be 0 or 1")
881 : : end if
882 [ - + ]: 72 : ABI_FREE(ranks_io)
883 : 144 : ABI_FREE(rank_has_cg)
884 : : end if
885 : :
886 : : ! Open the file in parallel inside comm_mpiio.
887 : 2862 : call xmpi_barrier(comm_cell)
888 [ + + ]: 2862 : if (comm_mpiio == xmpi_comm_null) goto 100
889 : :
890 : : ncerr = nf90_einval
891 : : #ifdef HAVE_NETCDF_MPI
892 : : ncerr = nf90_open(path, mode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_write),&
893 : 2826 : comm=comm_mpiio, info=xmpio_info, ncid=ncid)
894 : : #endif
895 [ - + ]: 2826 : NCF_CHECK_MSG(ncerr, sjoin("open_par: ", path))
896 : :
897 : : ! Use individual IO (default) for the G-vectors [3, mpw, nkpt]
898 [ - + ]: 2826 : NCF_CHECK(nf90_inq_varid(ncid, "reduced_coordinates_of_plane_waves", kg_varid))
899 : 2826 : spin = 1; ikg=0
900 [ + + ]: 41648 : do ikpt=1,nkpt
901 : 38822 : npw_k = npwarr(ikpt)
902 [ + + ]: 38822 : if (proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,1,spin,me_kpt)) cycle
903 [ + + + + ]: 231637 : ncerr = nf90_put_var(ncid, kg_varid, kg(:, 1+ikg:), start=[1,1,ikpt], count=[3,npw_k,1])
904 [ - + ]: 33091 : NCF_CHECK_MSG(ncerr, "put_kg_k")
905 : 41648 : ikg = ikg + npw_k
906 : : end do
907 : :
908 [ - + ]: 2826 : NCF_CHECK(nf90_inq_varid(ncid, "coefficients_of_wavefunctions", cg_varid))
909 : :
910 : 2826 : use_collective = (response == 0) ! or (response == 1 .and. nproc_cell == 1)
911 : :
912 [ + + ]: 2826 : if (use_collective) then
913 : 2813 : call wrtout(std_out,"Using collective IO for the CGs")
914 : : ! Use collective IO for the CGs
915 : : ncerr = nf90_einval
916 : : #ifdef HAVE_NETCDF_MPI
917 : 2813 : ncerr = nf90_var_par_access(ncid, cg_varid, nf90_collective)
918 : : #endif
919 [ - + ]: 2813 : NCF_CHECK(ncerr)
920 : :
921 [ + + ]: 35338 : do cnt=1,max_cnt
922 : 32525 : ikpt = iter2kscgkg(1,cnt)
923 : 32525 : spin = iter2kscgkg(2,cnt)
924 : 32525 : icg = iter2kscgkg(3,cnt)
925 : 32525 : ikg = iter2kscgkg(4,cnt)
926 : :
927 : : ! The array on file has shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
928 [ + + ]: 35338 : if (ikpt /= 0) then
929 : 32495 : nband_k = nband(ikpt + (spin-1)*nkpt)
930 : 32495 : npw_k = npwarr(ikpt)
931 : :
932 : : !write(std_out,*)"About to write ikpt, spin, npw_k, icg: ",ikpt, spin, npw_k, icg
933 : : ncerr = nf90_put_var(ncid, cg_varid, cg(:, 1+icg:), start=[1,1,1,1,ikpt,spin], &
934 [ + + + + ]: 422435 : count=[2,npw_k,nspinor,nband_k,1,1])
935 [ - + ]: 32495 : NCF_CHECK_MSG(ncerr, "writing cg")
936 : : else
937 : : ! This happens when nkpt * nsppol // nprocs != 0
938 : : ! Note that we are using collective MPI-IO hence all processors must call put_var
939 : : ! Here we re-write the ug(0) of the first (k-point, spin) treated by the node.
940 : 30 : ikpt = iter2kscgkg(1,1)
941 : 30 : spin = iter2kscgkg(2,1)
942 [ + + ]: 210 : ncerr = nf90_put_var(ncid, cg_varid, cg, start=[1,1,1,1,ikpt,spin], count=[1,1,1,1,1,1])
943 [ - + ]: 30 : NCF_CHECK_MSG(ncerr, "re-writing cg")
944 : : end if
945 : : end do
946 : :
947 : : else
948 : 13 : call wrtout(std_out, "Using individual IO for the CGs")
949 : : ! Individual IO of the CGs (for debugging purposes)
950 : 13 : icg = 0
951 [ + + ]: 26 : do spin=1,nsppol
952 [ + + ]: 1722 : do ikpt=1,nkpt
953 : 1696 : nband_k = nband(ikpt + (spin-1)*nkpt)
954 : 1696 : npw_k = npwarr(ikpt)
955 [ - + ]: 1696 : if (proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,spin,me_kpt)) cycle
956 [ + + + + ]: 22048 : ncerr = nf90_put_var(ncid, cg_varid, cg(:, 1+icg:), start=[1,1,1,1,ikpt,spin], count=[2,npw_k,nspinor,nband_k,1,1])
957 [ - + ]: 1696 : NCF_CHECK_MSG(ncerr, "put_var")
958 : 1709 : icg = icg+npw_k*my_nspinor*nband_k
959 : : end do
960 : : end do
961 : : end if
962 : :
963 [ - + ]: 2826 : NCF_CHECK(nf90_close(ncid))
964 : :
965 : 2862 : 100 call xmpi_barrier(comm_cell)
966 [ - + ]: 2862 : ABI_FREE(iter2kscgkg)
967 : :
968 : 2862 : done = .True.
969 : 5724 : call cwtime_report(" collective ncwrite", cpu, wall, gflops)
970 : :
971 : : else ! single_writer
972 [ # # ]: 0 : if (nproc_cell > 1) then
973 : 0 : ABI_WARNING("Slow version without MPI-IO support. Processors send data to master...")
974 : : else
975 : 0 : call wrtout(std_out, "Using netcdf library without MPI-IO support")
976 : : end if
977 : :
978 [ # # # # : 0 : ABI_MALLOC(kg_k,(3,mpw))
# # # # #
# # # ]
979 [ # # # # : 0 : ABI_MALLOC_OR_DIE(cg_k,(2,mpw*my_nspinor*mband), ierr)
# # # # #
# # # #
# ]
980 : :
981 [ # # ]: 0 : if (iam_master) then
982 : 0 : call wfk%open_write(hdr,path,formeig,iomode,get_unit(),xmpi_comm_self,write_hdr=.True.)
983 : :
984 [ # # ]: 0 : NCF_CHECK(crystal%ncwrite(wfk%fh))
985 : : !write(std_out,*)"after crystal_ncwrite"
986 : :
987 : : ! Write eigenvalues and occupations (these arrays are not MPI-distributed)
988 [ # # ]: 0 : if (response == 0) then
989 [ # # ]: 0 : NCF_CHECK(gs_ebands%ncwrite(wfk%fh))
990 : : else
991 : 0 : call ncwrite_eigen1_occ(wfk%fh, nband, mband, nkpt, nsppol, eigen, occ3d)
992 : : end if
993 : : end if
994 : :
995 : 0 : icg = 0
996 [ # # ]: 0 : do spin=1,nsppol
997 : 0 : ikg = 0
998 [ # # ]: 0 : do ikpt=1,nkpt
999 : 0 : nband_k = nband(ikpt + (spin-1)*nkpt)
1000 : 0 : npw_k = npwarr(ikpt)
1001 : 0 : mtag = ikpt+(spin-1)*nkpt
1002 : :
1003 : 0 : call xmpi_barrier(comm_cell)
1004 : :
1005 : : ! Transfer the wavefunctions and the g-vectors to the master processor
1006 [ # # ]: 0 : source = minval(mpi_enreg%proc_distrb(ikpt,1:nband_k,spin))
1007 : 0 : ihave_data = (source==me_kpt)
1008 : :
1009 : 0 : action=0
1010 [ # # ]: 0 : if (iam_master .and. ihave_data) action=1 ! I am the master node, and I have the data in cg
1011 [ # # ]: 0 : if (.not.iam_master.and.ihave_data) action=2 ! I am not the master, and I have the data => send to master
1012 [ # # ]: 0 : if (iam_master.and..not.ihave_data) action=3 ! I am the master, and I receive the data
1013 : :
1014 [ # # ]: 0 : if (action==1) then ! Copy from kg and cg
1015 [ # # # # ]: 0 : kg_k(:,1:npw_k) = kg(:,ikg+1:ikg+npw_k)
1016 : 0 : call cg_zcopy(npw_k*my_nspinor*nband_k, cg(1,icg+1), cg_k)
1017 : : end if
1018 : :
1019 : : ! Exchange data
1020 [ # # ]: 0 : if (action==2.or.action==3) then
1021 : 0 : call timab(48,1,tsec)
1022 [ # # ]: 0 : if (action==2) then
1023 : 0 : call xmpi_exch(kg(:,1+ikg:npw_k+ikg),3*npw_k,source,kg_k,master,comm_cell,2*mtag+1,ierr)
1024 : : call xmpi_exch(cg(:,icg+1:icg+nband_k*npw_k*my_nspinor),2*nband_k*npw_k*my_nspinor,&
1025 : 0 : source,cg_k,master,comm_cell,2*mtag+2,ierr)
1026 : : else
1027 : 0 : call xmpi_exch(kg_k,3*npw_k,source,kg_k,master,comm_cell,2*mtag+1,ierr)
1028 : 0 : call xmpi_exch(cg_k,2*nband_k*npw_k*my_nspinor,source,cg_k,master,comm_cell,2*mtag+2,ierr)
1029 : : end if
1030 : 0 : call timab(48,2,tsec)
1031 : : end if
1032 : :
1033 : : ! Master writes this block of bands.
1034 [ # # ]: 0 : if (iam_master) then
1035 [ # # ]: 0 : if (response == 0) then
1036 : : call wfk%write_band_block([1,nband_k],ikpt,spin,xmpio_single,kg_k=kg_k,cg_k=cg_k,&
1037 [ # # ]: 0 : eig_k=gs_ebands%eig(:,ikpt,spin),occ_k=gs_ebands%occ(:,ikpt,spin))
1038 : : !write(std_out,*)"cg_k",cg_k(:,1:2)
1039 : : else
1040 [ # # ]: 0 : call wfk%write_band_block([1,nband_k],ikpt,spin,xmpio_single,kg_k=kg_k,cg_k=cg_k)
1041 : : end if
1042 : : end if
1043 : :
1044 [ # # ]: 0 : if (.not.(proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,spin,me_kpt))) then
1045 : 0 : icg = icg+npw_k*my_nspinor*nband_k
1046 : 0 : ikg = ikg+npw_k
1047 : : end if
1048 : :
1049 : : end do !ikpt
1050 : : end do !spin
1051 : :
1052 [ # # ]: 0 : ABI_FREE(kg_k)
1053 [ # # ]: 0 : ABI_FREE(cg_k)
1054 : :
1055 [ # # ]: 0 : if (iam_master) call wfk%close()
1056 : 0 : call xmpi_barrier(comm_cell)
1057 : :
1058 : 0 : done = .True.
1059 : 0 : call cwtime_report(" individual ncwrite", cpu, wall, gflops)
1060 : : end if
1061 : :
1062 : : else ! not same_layout
1063 : :
1064 [ + - ]: 216 : if (nctk_has_mpiio) then
1065 : : call wrtout(std_out, &
1066 [ + - ]: 216 : sjoin("scattered data. writing WFK file",trim(path),", with iomode: ",iomode2str(iomode)), 'PERS', do_flush=.True.)
1067 : :
1068 : : ! master write the metadata.
1069 [ + + ]: 216 : if (xmpi_comm_rank(comm_cell) == master) then
1070 : : ncerr = nf90_einval
1071 : : #ifdef HAVE_NETCDF_MPI
1072 : : ncerr = nf90_create(path, cmode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_write), &
1073 : 48 : comm=xmpi_comm_self, info=xmpio_info, ncid=ncid)
1074 : : #endif
1075 [ - + ]: 48 : NCF_CHECK_MSG(ncerr, sjoin("create_par:", path))
1076 : :
1077 : 48 : call wfk_ncdef_dims_vars(ncid, hdr, fform2, write_hdr=.True.)
1078 [ - + ]: 48 : NCF_CHECK(crystal%ncwrite(ncid))
1079 : :
1080 : : ! Write eigenvalues and occupations (these arrays are not MPI-distributed)
1081 [ + - ]: 48 : if (response == 0) then
1082 [ - + ]: 48 : NCF_CHECK(gs_ebands%ncwrite(ncid))
1083 : : else
1084 : 0 : call ncwrite_eigen1_occ(ncid, nband, mband, nkpt, nsppol, eigen, occ3d)
1085 : : end if
1086 : :
1087 [ - + ]: 48 : NCF_CHECK(nf90_close(ncid))
1088 : : end if
1089 : :
1090 : : ! Reopen the file inside comm_cell
1091 : 216 : call xmpi_barrier(comm_cell)
1092 : : ncerr = nf90_einval
1093 : : #ifdef HAVE_NETCDF_MPI
1094 : : ncerr = nf90_open(path, mode=ior(ior(nf90_netcdf4, nf90_mpiio), nf90_write), &
1095 : 216 : comm=comm_cell, info=xmpio_info, ncid=ncid)
1096 : : #endif
1097 [ - + ]: 216 : NCF_CHECK_MSG(ncerr, sjoin("create_par:", path))
1098 : :
1099 : : ! Get var ids
1100 [ - + ]: 216 : NCF_CHECK(nf90_inq_varid(ncid, "reduced_coordinates_of_plane_waves", kg_varid))
1101 [ - + ]: 216 : NCF_CHECK(nf90_inq_varid(ncid, "coefficients_of_wavefunctions", cg_varid))
1102 : :
1103 : : ! Write the G-vectors
1104 : 216 : ikg = 0
1105 [ + + ]: 2024 : do ikpt=1,nkpt
1106 : 1808 : npw_k = npwarr(ikpt)
1107 : 1808 : npwtot_k = hdr%npwarr(ikpt)
1108 [ + + ]: 1808 : if (proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,1,1,me_kpt)) cycle
1109 [ - + ]: 1018 : if (me_spinor /= 0) cycle
1110 : : !write(std_out,*)"In G-vector loop ",ikpt,", with me_cell me_kpt me_band, me_spinor ",me_cell,me_kpt,me_band,me_spinor
1111 : : !write(std_out,*)" with ik, npw_k, npwtot_k: ",ikpt, npw_k, npwtot_k
1112 : :
1113 [ - + + - : 3054 : ABI_MALLOC(ind_cg_mpi_to_seq, (npw_k))
+ - + - -
+ - + ]
1114 [ + - ]: 1018 : if (allocated(mpi_enreg%my_kgtab)) then
1115 : 1018 : ikpt_this_proc = mpi_enreg%my_kpttab(ikpt)
1116 [ - + + + ]: 174545 : ind_cg_mpi_to_seq = mpi_enreg%my_kgtab(1:npw_k,ikpt_this_proc)
1117 : : else
1118 [ # # ]: 0 : ABI_CHECK(nproc_fft==1, "nproc_fft !=1 and my_kgtab not allocated")
1119 [ # # # # : 0 : ind_cg_mpi_to_seq(1:npw_k) = [(ipw, ipw=1,npw_k)]
# # ]
1120 : : end if
1121 : :
1122 [ - + + - : 2350134 : ABI_CALLOC(kg_k, (3, npwtot_k))
+ - + - -
+ - + - +
+ + + + ]
1123 [ + + ]: 173527 : do ipw=1,npw_k
1124 [ + + ]: 691054 : kg_k(:, ind_cg_mpi_to_seq(ipw)) = kg(:,ikg+ipw)
1125 : : end do
1126 : 1018 : call xmpi_sum_master(kg_k,master,comm_bandfft,ierr)
1127 [ + + ]: 1018 : if (xmpi_comm_rank(comm_bandfft) == master) then
1128 [ + + + + ]: 2408 : ncerr = nf90_put_var(ncid, kg_varid, kg_k, start=[1,1,ikpt], count=[3,npwtot_k,1])
1129 [ - + ]: 344 : NCF_CHECK_MSG(ncerr, "putting kg_k")
1130 : : else
1131 [ + + ]: 2696 : ncerr = nf90_put_var(ncid, kg_varid, kg_k, start=[1,1,ikpt], count=[0,0,0])
1132 [ - + ]: 674 : NCF_CHECK_MSG(ncerr, "skiping kg_k")
1133 : : end if
1134 [ - + ]: 1018 : ABI_FREE(kg_k)
1135 : :
1136 : : ! gblock contains block-distributed G-vectors inside comm_fft
1137 : : !call kg2seqblocks(npwtot_k,npw_k,kg(:,ikg+1:),ind_cg_mpi_to_seq,comm_fft,start_pwblock,count_pwblock,gblock)
1138 : : !write(std_out,*)"gblock(:,2)",gblock(:,2)
1139 : : !ncerr = nf90_put_var(ncid, kg_varid, gblock, start=[1,start_pwblock,ikpt], count=[3,count_pwblock,1])
1140 : : !NCF_CHECK_MSG(ncerr, "putting kg_k")
1141 : : !ABI_FREE(gblock)
1142 : :
1143 [ - + ]: 1018 : ABI_FREE(ind_cg_mpi_to_seq)
1144 : :
1145 : 2024 : ikg = ikg+npw_k
1146 : : end do
1147 : :
1148 : : ! Write wavefunctions
1149 [ - + ]: 216 : if (response == 1) then
1150 : : ! The cg1 is allocated with size mcg1=mpw1*dtset%nspinor*dtset%mband*mk1mem_rbz*dtset%nsppol (see dfpt_looppert)
1151 : : ! hence bands are not MPI distributed, this is the reason why we have to use the loop over bands
1152 : : ! and the check on nwrites
1153 : 0 : icg = 0
1154 [ # # ]: 0 : do spin=1,nsppol
1155 [ # # ]: 0 : do ikpt=1,nkpt
1156 : 0 : nband_k = nband(ikpt + (spin-1)*nkpt)
1157 : 0 : npw_k = npwarr(ikpt)
1158 : 0 : npwtot_k = hdr%npwarr(ikpt)
1159 : 0 : nwrites = 0
1160 [ # # ]: 0 : do iband=1,nband_k
1161 [ # # ]: 0 : if (mpi_enreg%proc_distrb(ikpt,iband,spin)/=me_kpt) cycle
1162 : : ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
1163 : 0 : nwrites = nwrites + 1
1164 : 0 : ii = 1 + (iband-1)*npw_k*my_nspinor + icg
1165 [ # # # # ]: 0 : ncerr = nf90_put_var(ncid, cg_varid, cg(1:,ii:), start=[1,1,1,iband,ikpt,spin], count=[2,npwtot_k,nspinor,1,1,1])
1166 [ # # ]: 0 : NCF_CHECK_MSG(ncerr, "put_var cg")
1167 : : end do ! iband
1168 [ # # ]: 0 : if (nwrites /= 0) icg = icg + npw_k*my_nspinor*nband_k
1169 : : end do !ikpt
1170 : : end do !spin
1171 : :
1172 : : else
1173 : 216 : icg = 0
1174 [ + + ]: 440 : do spin=1,nsppol
1175 [ + + ]: 2256 : do ikpt=1,nkpt
1176 : 1816 : nband_k = nband(ikpt + (spin-1)*nkpt)
1177 : 1816 : npw_k = npwarr(ikpt)
1178 : 1816 : npwtot_k = hdr%npwarr(ikpt)
1179 : :
1180 [ + + ]: 1816 : if (proc_distrb_cycle(mpi_enreg%proc_distrb,ikpt,1,nband_k,spin,me_kpt)) cycle
1181 : : !write(std_out,*)"In u(g)-vector loop ",ikpt,", with me_cell me_kpt me_band, me_spinor ",me_cell,me_kpt,me_band,me_spinor
1182 : : !write(std_out,*)"nband_k, npw_k, npwtot_k: ",nband_k, npw_k, npwtot_k
1183 : :
1184 : : ! Part taken from writewf.
1185 : : ! ind_cg_mpi_to_seq maps the npw_k*nspinor local coefficients of a band to the (npwtot_k, nspinor)
1186 : : ! block on file: the offset of the second spinor component is the total number of planewaves.
1187 : 1024 : npwtot=npwtot_k; npwso=npw_k*nspinor
1188 : 1024 : npwsotot=npwso
1189 : 1024 : nspinortot = min(2,(1+mpi_enreg%paral_spinor)*nspinor)
1190 : :
1191 [ - + + - : 3072 : ABI_MALLOC(ind_cg_mpi_to_seq, (npwso))
+ - + - -
+ - + ]
1192 [ + - ]: 1024 : if (allocated(mpi_enreg%my_kgtab)) then
1193 : 1024 : ikpt_this_proc=mpi_enreg%my_kpttab(ikpt)
1194 [ + + ]: 2048 : do ispinor=1,nspinor
1195 : 1024 : ispinor_index=ispinor
1196 [ - + ]: 1024 : if (nproc_spinor > 1) ispinor_index = mpi_enreg%me_spinor + 1
1197 : : ind_cg_mpi_to_seq(1+npw_k*(ispinor-1):npw_k*ispinor)=npwtot*(ispinor_index-1) &
1198 [ + + ]: 174946 : & + mpi_enreg%my_kgtab(1:npw_k,ikpt_this_proc)
1199 : : end do
1200 : : else
1201 [ # # ]: 0 : ABI_CHECK(nproc_fft==1, "nproc_fft !=1 and my_kgtab not allocated")
1202 [ # # # # : 0 : ind_cg_mpi_to_seq(1:npwso) = [(ipw, ipw=1,npwso)]
# # ]
1203 : : end if
1204 : :
1205 : : ! TODO: Blocking + collective IO in comm_cell
1206 : : !band_blocksize = nband_k/nbdblock
1207 : : !ibandmin = 1
1208 : : !step=min(ii,MAXBAND, nband_disk)
1209 : : !do iblock=1,nband_disk/step+1
1210 : : ! ibandmax = min(ibandmin+step-1, nband_disk)
1211 : : ! nband_block = ibandmax-ibandmin+1
1212 : : !gbase = icg + npw_k*my_nspinor*nband_k
1213 : :
1214 : : ! Each MPI proc write bcount bands with npwtot_k G-vectors starting from bstart.
1215 : : call cg2seqblocks(npwtot_k,npw_k,nspinor,nband_k,cg(:,icg+1:),ind_cg_mpi_to_seq,comm_bandfft,&
1216 [ - + ]: 1024 : bstart,bcount,my_cgblock)
1217 : :
1218 : : ! The coefficients_of_wavefunctions on file have shape [cplex, mpw, nspinor, mband, nkpt, nsppol]
1219 [ + + + + ]: 13312 : ncerr = nf90_put_var(ncid, cg_varid, my_cgblock, start=[1,1,1,bstart,ikpt,spin], count=[2,npwtot_k,nspinor,bcount,1,1])
1220 [ - + ]: 1024 : NCF_CHECK_MSG(ncerr, "put_var cg")
1221 : :
1222 [ - + ]: 1024 : ABI_FREE(my_cgblock)
1223 [ - + ]: 1024 : ABI_FREE(ind_cg_mpi_to_seq)
1224 : :
1225 : 2040 : icg = icg+npw_k*my_nspinor*nband_k
1226 : : end do !ikpt
1227 : : end do !spin
1228 : : end if
1229 : :
1230 [ - + ]: 216 : NCF_CHECK(nf90_close(ncid))
1231 : 216 : done = .True.
1232 : :
1233 : 216 : call cwtime_report(" scattered ncwrite", cpu, wall, gflops)
1234 : : end if !nctk_has_mpiio
1235 : : end if
1236 : :
1237 : 3078 : call crystal%free()
1238 [ + + ]: 3078 : if (response == 0) call gs_ebands%free()
1239 : :
1240 [ - + ]: 3078 : ABI_FREE(occ3d)
1241 : :
1242 : : DBG_EXIT("COLL")
1243 : :
1244 [ - + - + : 6156 : end subroutine cg_ncwrite
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - -
- - - - -
- - - - -
- - - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + ]
1245 : : !!***
1246 : :
1247 : : !!****f* m_iowf/prtkbff
1248 : : !! NAME
1249 : : !! prtkbff
1250 : : !!
1251 : : !! FUNCTION
1252 : : !! Write KB form factors to WFK in netcdf format.
1253 : : !! Only master works. G-vectors are read from file to avoid
1254 : : !! having to deal with paral_kgb distribution.
1255 : : !!
1256 : : !! INPUTS
1257 : : !!
1258 : : !! OUTPUT
1259 : : !!
1260 : : !! SOURCE
1261 : :
1262 : 2 : subroutine prtkbff(filnam, hdr, psps, prtvol)
1263 : :
1264 : : !Arguments ------------------------------------
1265 : : character(len=*),intent(in) :: filnam
1266 : : type(hdr_type),intent(in) :: hdr
1267 : : type(pseudopotential_type),intent(in) :: psps
1268 : : integer,intent(in) :: prtvol
1269 : :
1270 : : !Local variables-------------------------------
1271 : : !scalars
1272 : : character(len=fnlen) :: path
1273 : : character(len=500) :: msg
1274 : : integer :: ncid, ncerr, kg_varid, mpw_disk, npwk_disk, ikpt, iat, iproj
1275 : 2 : integer,allocatable :: kg_disk(:,:)
1276 : 2 : real(dp),allocatable :: vkb(:,:,:),vkbd(:,:,:),vkbsign(:,:)
1277 [ + + + + : 104 : type(crystal_t) :: crystal
+ + + + +
+ + + + +
+ + + + ]
1278 : :
1279 : : ! *************************************************************************
1280 : :
1281 : 2 : path = nctk_ncify(filnam)
1282 : 2 : call wrtout(std_out, sjoin("Writing KB form factors to:", path))
1283 [ - + ]: 2 : NCF_CHECK(nctk_open_modify(ncid, path, xmpi_comm_self))
1284 [ - + ]: 2 : NCF_CHECK(nf90_inq_varid(ncid, "reduced_coordinates_of_plane_waves", kg_varid))
1285 [ + + ]: 268 : mpw_disk = maxval(hdr%npwarr)
1286 : :
1287 : : ! Dimensions needed by client code to allocate memory when reading.
1288 : : ncerr = nctk_def_dims(ncid, [ &
1289 : : nctkdim_t("mproj", psps%mproj), &
1290 : : nctkdim_t("mpsang", psps%mpsang), &
1291 : : nctkdim_t("mpssoang", psps%mpssoang), &
1292 : : nctkdim_t("lnmax", psps%lnmax), &
1293 : : nctkdim_t("lmnmax", psps%lmnmax) &
1294 [ + + ]: 12 : ])
1295 [ - + ]: 2 : NCF_CHECK(ncerr)
1296 : :
1297 [ + + ]: 4 : ncerr = nctk_def_iscalars(ncid, [character(len=nctk_slen) :: "mpspso"])
1298 [ - + ]: 2 : NCF_CHECK(ncerr)
1299 : :
1300 : : ! Write indlmn table (needed to access vkb arrays)
1301 : : ncerr = nctk_def_arrays(ncid, [ &
1302 : : nctkarr_t("indlmn", "int", "six, lmnmax, number_of_atom_species"), &
1303 : : nctkarr_t("vkbsign", "dp", "lnmax, number_of_atom_species"), &
1304 : : nctkarr_t("vkb", "dp", "max_number_of_coefficients, lnmax, number_of_atom_species, number_of_kpoints"), &
1305 : : nctkarr_t("vkbd", "dp", "max_number_of_coefficients, lnmax, number_of_atom_species, number_of_kpoints") &
1306 [ + + ]: 10 : ], defmode=.True.)
1307 [ - + ]: 2 : NCF_CHECK(ncerr)
1308 : :
1309 : : ! Switch to write mode.
1310 [ - + ]: 2 : NCF_CHECK(nctk_set_datamode(ncid))
1311 [ - + ]: 2 : NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "indlmn"), psps%indlmn))
1312 : :
1313 : : ncerr = nctk_write_iscalars(ncid, [character(len=nctk_slen) :: &
1314 : : "mpspso"], &
1315 [ + + + + ]: 6 : [psps%mpspso])
1316 [ - + ]: 2 : NCF_CHECK(ncerr)
1317 : :
1318 : : ! Calculate KB form factors and derivatives.
1319 : : ! The arrays are allocated with lnmax to support pseudos with more than projector.
1320 : : ! Note that lnmax takes into account lloc hence arrays are in packed form and should be
1321 : : ! accessed with the indices provided by psps%indlmn.
1322 [ - + + - : 8 : ABI_MALLOC(vkbsign, (psps%lnmax, psps%ntypat))
- + + - +
- + - - +
- + ]
1323 [ - + + - : 10 : ABI_MALLOC(vkb, (mpw_disk, psps%lnmax, psps%ntypat))
- + + - -
+ + - + -
+ - - + -
+ ]
1324 [ - + + - : 8 : ABI_MALLOC(vkbd, (mpw_disk, psps%lnmax, psps%ntypat))
- + + - -
+ + - - +
- + ]
1325 [ - + + - : 6 : ABI_MALLOC(kg_disk, (3, mpw_disk))
+ - + - -
+ - + ]
1326 : :
1327 : 2 : crystal = hdr%get_crystal()
1328 : :
1329 : : ! For each k-point: read full G-vector list from file, compute KB data and write to file.
1330 [ + + ]: 268 : do ikpt=1,hdr%nkpt
1331 : 266 : npwk_disk = hdr%npwarr(ikpt)
1332 [ + + + + : 1862 : NCF_CHECK(nf90_get_var(ncid, kg_varid, kg_disk, start=[1, 1, ikpt], count=[3, npwk_disk, 1]))
- + - - -
- ]
1333 [ + + + + : 3518682 : vkb = zero; vkbd = zero
+ + + + +
+ + + ]
1334 : 266 : call calc_vkb(crystal, psps, hdr%kptns(:, ikpt), npwk_disk, mpw_disk, kg_disk, vkbsign, vkb, vkbd)
1335 : :
1336 [ + + ]: 266 : if (ikpt == 1) then
1337 : : ! This for the automatic tests.
1338 [ - + ]: 2 : NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "vkbsign"), vkbsign))
1339 [ + + ]: 2 : if(prtvol >= 2) then
1340 : 1 : write(msg,'(a)') 'prtkbff: writing first and last G-components of the KB form factors'
1341 : 1 : call wrtout(ab_out, msg)
1342 [ + + ]: 3 : do iat=1,psps%ntypat
1343 : 2 : write(msg,'(a10,i5)') 'atype ',iat
1344 : 2 : call wrtout(ab_out, msg)
1345 [ + + ]: 5 : do iproj=1,psps%lnmax
1346 : : write(msg,'(a10,i5,a,a10,e12.4,a,2(a10,2e12.4,a))') &
1347 : 2 : 'projector ', iproj,ch10, &
1348 : 2 : 'vkbsign ', vkbsign(iproj,iat), ch10, &
1349 : 2 : 'vkb ', vkb(1,iproj,iat), vkb(npwk_disk,iproj,iat), ch10, &
1350 : 4 : 'vkbd ', vkbd(1,iproj,iat), vkbd(npwk_disk,iproj,iat), ''
1351 : 4 : call wrtout(ab_out, msg)
1352 : : end do
1353 : : end do
1354 : : end if
1355 : : end if
1356 : :
1357 [ + + - + : 1330 : NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "vkb"), vkb, start=[1, 1, 1, ikpt]))
- - ]
1358 [ + + - + : 1332 : NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "vkbd"), vkbd, start=[1, 1, 1, ikpt]))
- - ]
1359 : : end do
1360 [ - + ]: 2 : NCF_CHECK(nf90_close(ncid))
1361 : :
1362 [ - + ]: 2 : ABI_FREE(kg_disk)
1363 [ - + ]: 2 : ABI_FREE(vkbsign)
1364 [ - + ]: 2 : ABI_FREE(vkb)
1365 [ - + ]: 2 : ABI_FREE(vkbd)
1366 : 2 : call crystal%free()
1367 : :
1368 [ - + - + : 2 : end subroutine prtkbff
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + ]
1369 : : !!***
1370 : :
1371 : : !----------------------------------------------------------------------
1372 : :
1373 : : !!****f* m_iowf/ncwrite_eigen1_occ
1374 : : !! NAME
1375 : : !! ncwrite_eigen1_occ
1376 : : !!
1377 : : !! FUNCTION
1378 : : !! Write the first order DFPT eigenvalues and the occupations.
1379 : : !!
1380 : : !! INPUTS
1381 : : !! ncid=Netcdf file handler.
1382 : : !! nband(nkpt*nsppol)=Number of bands.
1383 : : !! mband=maximum number of bands
1384 : : !! nkpt=Total number of k points
1385 : : !! nsppol=1 for unpolarized, 2 for spin-polarized
1386 : : !! eigen((2*mband**2*nkpt*nsppol)= eigenvalues (hartree) for all bands at each k point
1387 : : !! occ3d(mband*nkpt*nsppol)=occupations for all bands at each k point
1388 : : !! Note that occ3d differes from the occ arrays used in the rest of the code.
1389 : : !! occ3d is an array with constant stride `mband` whereas abinit (for reasons that are not clear to me)
1390 : : !! packs the occupations in a 1d vector with a k-dependent separation (nband_k).
1391 : : !! occ and occ3d differ only if nband is k-dependent but you should never assume this, hence
1392 : : !! remember to *convert* occ into occ3d before calling this routine.
1393 : : !!
1394 : : !! SOURCE
1395 : :
1396 : 13 : subroutine ncwrite_eigen1_occ(ncid, nband, mband, nkpt, nsppol, eigen, occ3d)
1397 : :
1398 : : !Arguments ------------------------------------
1399 : : !scalars
1400 : : integer,intent(in) :: ncid,mband,nkpt,nsppol
1401 : : !arrays
1402 : : integer, intent(in) :: nband(nkpt*nsppol)
1403 : : real(dp),intent(in) :: eigen(2*mband**2*nkpt*nsppol),occ3d(mband,nkpt,nsppol)
1404 : :
1405 : : !Local variables-------------------------------
1406 : : !scalars
1407 : : integer :: idx,spin,ikpt,nband_k,ib2,ib1
1408 : : integer :: ncerr,occ_varid,h1mat_varid,eigens_varid
1409 : : !arrays
1410 : 13 : real(dp),allocatable :: h1mat(:,:,:,:,:), fake_eigens(:,:,:)
1411 : :
1412 : : ! *************************************************************************
1413 : :
1414 : : ! Declare h1 array with abinit conventions
1415 : : ! Cannot use a 3D array since eigen1 are packed and one could have different number of bands
1416 : : ! per k-points. (this is not an official etsf-io variable)!
1417 : : ! elements in eigen are packed in the first positions.
1418 : : ! Remember that eigen are not MPI-distributed so no communication is needed
1419 : 13 : idx=1
1420 [ - + + - : 99451 : ABI_CALLOC(h1mat, (2,mband,mband,nkpt,nsppol))
- + + - -
+ + - - +
+ - + - +
- - + - +
+ + + + +
+ + + +
+ ]
1421 [ + + ]: 26 : do spin=1,nsppol
1422 [ + + ]: 1722 : do ikpt=1,nkpt
1423 : 1696 : nband_k = nband(ikpt + (spin-1)*nkpt)
1424 [ + + ]: 8749 : do ib2=1,nband_k
1425 [ + + ]: 38944 : do ib1=1,nband_k
1426 [ + + ]: 90624 : h1mat(:,ib1,ib2,ikpt,spin) = eigen(idx:idx+1)
1427 : 37248 : idx=idx+2
1428 : : end do
1429 : : end do
1430 : : end do
1431 : : end do
1432 : :
1433 [ - + ]: 13 : NCF_CHECK(nctk_set_defmode(ncid))
1434 : :
1435 : : ncerr = nctk_def_arrays(ncid, nctkarr_t('h1_matrix_elements', "dp", &
1436 : 13 : &"complex, max_number_of_states, max_number_of_states, number_of_kpoints, number_of_spins"))
1437 [ - + ]: 13 : NCF_CHECK(ncerr)
1438 : :
1439 [ - + ]: 13 : NCF_CHECK(nf90_inq_varid(ncid, "occupations", occ_varid))
1440 [ - + ]: 13 : NCF_CHECK(nf90_inq_varid(ncid, "h1_matrix_elements", h1mat_varid))
1441 [ - + ]: 13 : NCF_CHECK(nf90_inq_varid(ncid, "eigenvalues", eigens_varid))
1442 : :
1443 : : ! Write data
1444 [ - + ]: 13 : NCF_CHECK(nctk_set_datamode(ncid))
1445 [ - + ]: 13 : NCF_CHECK_MSG(nf90_put_var(ncid, occ_varid, occ3d), "putting occ3d")
1446 [ - + ]: 13 : NCF_CHECK_MSG(nf90_put_var(ncid, h1mat_varid, h1mat), "putting h1mat")
1447 : :
1448 [ - + ]: 13 : ABI_FREE(h1mat)
1449 : :
1450 : : ! GS eigenvalues are set to zero.
1451 [ - + + - : 8814 : ABI_CALLOC(fake_eigens, (mband,nkpt,nsppol))
- + + - -
+ + - + -
+ - - + -
+ + + + +
+ + ]
1452 [ - + ]: 13 : NCF_CHECK_MSG(nf90_put_var(ncid, eigens_varid, fake_eigens), "putting fake eigens")
1453 [ - + ]: 13 : ABI_FREE(fake_eigens)
1454 : :
1455 : 13 : end subroutine ncwrite_eigen1_occ
1456 : : !!***
1457 : :
1458 : : !----------------------------------------------------------------------
1459 : :
1460 : : !!****f* m_iowf/kg2seqblocks
1461 : : !! NAME
1462 : : !! kg2seqblocks
1463 : : !!
1464 : : !! FUNCTION
1465 : : !!
1466 : : !! INPUTS
1467 : : !! npwtot_k
1468 : : !! npw_k=
1469 : : !! kg(3,npw_k)=reduced planewave coordinates.
1470 : : !! gmpi2seq
1471 : : !! comm_fft=FFT communicator
1472 : : !!
1473 : : !! OUTPUT
1474 : : !! start_pwblock
1475 : : !! count_pwblock
1476 : : !! gblock(:,:)
1477 : : !!
1478 : : !! SOURCE
1479 : :
1480 : : subroutine kg2seqblocks(npwtot_k,npw_k,kg_k,gmpi2seq,comm_fft,start_pwblock,count_pwblock,gblock)
1481 : :
1482 : : !Arguments ------------------------------------
1483 : : !scalars
1484 : : integer,intent(in) :: npwtot_k,npw_k,comm_fft
1485 : : integer,intent(out) :: start_pwblock,count_pwblock
1486 : : !arrays
1487 : : integer,intent(in) :: kg_k(3,npw_k),gmpi2seq(npw_k)
1488 : : integer,allocatable,intent(out) :: gblock(:,:)
1489 : :
1490 : : !Local variables-------------------------------
1491 : : !scalars
1492 : : integer :: me_fft,nproc_fft,rank,ig,igseq,maxnpw,ierr
1493 : : !arrays
1494 : : integer,allocatable :: igstart_rank(:),gbuf(:,:)
1495 : :
1496 : : ! *************************************************************************
1497 : :
1498 : : me_fft = xmpi_comm_rank(comm_fft); nproc_fft = xmpi_comm_size(comm_fft)
1499 : :
1500 : : ! Handle sequential case.
1501 : : if (nproc_fft == 1) then
1502 : : start_pwblock = 1; count_pwblock = npw_k
1503 : : ABI_MALLOC(gblock, (3, npw_k))
1504 : : gblock(:,:) = kg_k
1505 : : return ! DOH
1506 : : end if
1507 : :
1508 : : ABI_MALLOC(igstart_rank, (0:nproc_fft))
1509 : : igstart_rank = [(1 + (npwtot_k/nproc_fft)*rank, rank=0,nproc_fft-1), 1 + npwtot_k]
1510 : :
1511 : : ! Get max dimension for workspace array
1512 : : ! Cannot use npwtot_k / nproc_fft because G-vectors are not equally distributed.
1513 : : call xmpi_max(npw_k, maxnpw, comm_fft, ierr)
1514 : : ABI_MALLOC(gbuf, (3, maxnpw))
1515 : :
1516 : : do rank=0,nproc_fft-1
1517 : : start_pwblock = igstart_rank(rank)
1518 : : count_pwblock = igstart_rank(rank+1) - igstart_rank(rank)
1519 : :
1520 : : gbuf = 0
1521 : : do ig=1,npw_k
1522 : : igseq = gmpi2seq(ig)
1523 : : if (igseq >= igstart_rank(rank) .and. igseq < igstart_rank(rank+1)) then
1524 : : igseq = igseq - start_pwblock + 1
1525 : : !ABI_CHECK(igseq <= maxnpw, "boom")
1526 : : gbuf(:,igseq) = kg_k(:,ig)
1527 : : end if
1528 : : end do
1529 : : call xmpi_sum_master(gbuf,rank,comm_fft,ierr)
1530 : :
1531 : : if (me_fft == rank) then
1532 : : ABI_MALLOC_OR_DIE(gblock, (3, count_pwblock), ierr)
1533 : : gblock = gbuf(:, :count_pwblock)
1534 : : end if
1535 : : end do
1536 : :
1537 : : start_pwblock = igstart_rank(me_fft)
1538 : : count_pwblock = igstart_rank(me_fft+1) - igstart_rank(me_fft)
1539 : :
1540 : : ABI_FREE(gbuf)
1541 : : ABI_FREE(igstart_rank)
1542 : :
1543 : : end subroutine kg2seqblocks
1544 : : !!***
1545 : :
1546 : : !----------------------------------------------------------------------
1547 : :
1548 : : !!****f* m_iowf/cg2seqblocks
1549 : : !! NAME
1550 : : !! cg2seqblocks
1551 : : !!
1552 : : !! FUNCTION
1553 : : !!
1554 : : !! INPUTS
1555 : : !! npwtot_k=total number of planewaves
1556 : : !! npw_k=number of planewaves of this processor
1557 : : !! nspinor=number of spinor components (all of them on this processor)
1558 : : !! cg_k(2,npw_k*nspinor*nband_k)
1559 : : !! gmpi2seq(npw_k*nspinor)=index of each local coefficient in the (npwtot_k, nspinor) block on file
1560 : : !! comm_fft=FFT communicator
1561 : : !!
1562 : : !! OUTPUT
1563 : : !! bstart
1564 : : !! bcount
1565 : : !!
1566 : : !! SOURCE
1567 : :
1568 : 1024 : subroutine cg2seqblocks(npwtot_k,npw_k,nspinor,nband,cg_k,gmpi2seq,comm_bandfft,bstart,bcount,my_cgblock)
1569 : :
1570 : : !Arguments ------------------------------------
1571 : : !scalars
1572 : : integer,intent(in) :: npwtot_k,npw_k,nspinor,nband,comm_bandfft
1573 : : integer,intent(out) :: bstart,bcount
1574 : : !arrays
1575 : : integer,intent(in) :: gmpi2seq(npw_k*nspinor)
1576 : : real(dp),intent(in) :: cg_k(2,npw_k*nspinor,nband)
1577 : : real(dp),allocatable,intent(out) :: my_cgblock(:,:,:)
1578 : :
1579 : : !Local variables-------------------------------
1580 : : !scalars
1581 : : integer :: me,nprocs,rank,ig,ib,igseq,ierr,nbb,nbmax,band
1582 : : !arrays
1583 : 1024 : integer,allocatable :: bstart_rank(:)
1584 : 1024 : real(dp),allocatable :: cgbuf(:,:,:)
1585 : :
1586 : : ! *************************************************************************
1587 : :
1588 : 1024 : me = xmpi_comm_rank(comm_bandfft); nprocs = xmpi_comm_size(comm_bandfft)
1589 : :
1590 : : ! Handle sequential case.
1591 [ - + ]: 1024 : if (nprocs == 1) then
1592 : 0 : bstart = 1; bcount = nband
1593 [ # # # # : 0 : ABI_MALLOC_OR_DIE(my_cgblock, (2, npwtot_k*nspinor, nband), ierr)
# # # # #
# # # # #
# # # # #
# # # ]
1594 [ # # # # : 0 : my_cgblock = cg_k
# # # # #
# # # # #
# # ]
1595 : : return ! DOH
1596 : : end if
1597 : :
1598 [ - + + - : 3072 : ABI_MALLOC(bstart_rank, (0:nprocs))
+ - + - -
+ - + ]
1599 [ - + + + : 14224 : bstart_rank = [(1 + (nband/nprocs)*rank, rank=0,nprocs-1), 1 + nband]
+ + + - -
+ - - +
+ ]
1600 : :
1601 : : ! Allocate MPI buffer (same size on each MPI proc)
1602 : 1024 : nbmax = 0
1603 [ + + ]: 4400 : do rank=0,nprocs-1
1604 : 4400 : nbmax = max(nbmax, bstart_rank(rank+1) - bstart_rank(rank))
1605 : : end do
1606 [ - + + - : 4096 : ABI_MALLOC_OR_DIE(cgbuf, (2, npwtot_k*nspinor, nbmax), ierr)
- + + - +
- + - + -
- + - + ]
1607 : :
1608 : 1024 : nbb = bstart_rank(me+1) - bstart_rank(me)
1609 [ - + + - : 4096 : ABI_MALLOC_OR_DIE(my_cgblock, (2, npwtot_k*nspinor, nbb), ierr)
- + + - +
- + - + -
+ - - + +
- - + ]
1610 : :
1611 : : ! TODO: This should be replaced by gatherv but premature optimization....
1612 [ + + ]: 4400 : do rank=0,nprocs-1
1613 : 3376 : bstart = bstart_rank(rank)
1614 : 3376 : bcount = bstart_rank(rank+1) - bstart_rank(rank)
1615 : :
1616 [ + + ]: 15872 : do band=bstart, bstart+bcount-1
1617 : 12496 : ib = band - bstart + 1
1618 [ + + + + ]: 29910880 : cgbuf(:,:,ib) = zero
1619 [ + + ]: 2872744 : do ig=1,npw_k*nspinor
1620 : 2856872 : igseq = gmpi2seq(ig)
1621 [ + + ]: 8583112 : cgbuf(:,igseq,ib) = cg_k(:,ig,band)
1622 : : end do
1623 : : end do ! band
1624 : :
1625 : 3376 : call xmpi_sum_master(cgbuf,rank,comm_bandfft,ierr)
1626 [ + + + - : 8580060 : if (me == rank) my_cgblock = cgbuf(:,:,:bcount)
+ - + - -
+ - - + +
+ + + + ]
1627 : : end do ! rank
1628 : :
1629 : 1024 : bstart = bstart_rank(me)
1630 : 1024 : bcount = bstart_rank(me+1) - bstart_rank(me)
1631 : :
1632 [ - + ]: 1024 : ABI_FREE(cgbuf)
1633 : 1024 : ABI_FREE(bstart_rank)
1634 : :
1635 : : end subroutine cg2seqblocks
1636 : : !!***
1637 : :
1638 : : !----------------------------------------------------------------------
1639 : :
1640 : : END MODULE m_iowf
1641 : : !!***
|
