Line data Source code
1 : !!****p* ABINIT/fold2Bloch
2 : !! NAME
3 : !! fold2Bloch
4 : !!
5 : !! FUNCTION
6 : !! Main routine for the unfolding of the wavefuntion.
7 : !!
8 : !! COPYRIGHT
9 : !! Copyright (C) 2014-2026 ABINIT group (AB)
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 : !! For the initials of contributors, see ~abinit/doc/developers/contributors.txt .
14 : !!
15 : !! INPUTS
16 : !! (main program)
17 : !!
18 : !! OUTPUT
19 : !! (main program)
20 : !!
21 : !! NOTES
22 : !! folds= Array of folds in X,Y, and Z directions
23 : !!
24 : !! SOURCE
25 :
26 : #if defined HAVE_CONFIG_H
27 : #include "config.h"
28 : #endif
29 :
30 : #include "abi_common.h"
31 :
32 2 : program fold2Bloch
33 :
34 2 : use defs_basis
35 : use m_errors
36 : use m_abicore
37 : use m_wfk
38 : use m_xmpi
39 : use m_nctk
40 : use m_hdr
41 : use m_crystal
42 : use m_ebands
43 : use m_fold2block
44 : use netcdf
45 :
46 : use m_fstrings, only : strcat
47 : use m_io_tools, only : get_unit, iomode_from_fname, open_file, prompt
48 :
49 : implicit none
50 :
51 : !Arguments --------------------------------------------------------------
52 :
53 : !Local variables-------------------------------
54 : !scalars
55 : integer :: ikpt, iband,nspinor,nsppol,mband,nkpt,mcg,csppol, cspinor, nfold, iss, ii
56 : integer :: comm, my_rank, nargs, iomode, ncid, ncerr, fform, kunf_varid, weights_varid, eigunf_varid
57 : integer :: cg_b, count, outfile, outfile1, outfile2, lwcg, hicg, pos
58 : character(fnlen) :: fname, outname,seedname
59 : character(len=500) :: msg
60 8 : type(wfk_t) :: wfk
61 104 : type(crystal_t) :: cryst
62 2 : type(ebands_t) :: ebands
63 : !arrays
64 : integer :: folds(3),fold_matrix(3,3)
65 2 : integer, allocatable :: kg(:,:),nband(:), npwarr(:)
66 4 : real(dp), allocatable :: cg(:,:), eig(:),kpts(:,:), weights(:),coefc(:,:), nkval(:,:)
67 :
68 : !*************************************************************************
69 :
70 : !0) Change communicator for I/O (mandatory!)
71 2 : call abi_io_redirect(new_io_comm=xmpi_world)
72 :
73 2 : call xmpi_init()
74 2 : comm = xmpi_world; my_rank = xmpi_comm_rank(xmpi_world)
75 :
76 : !Initialize memory profiling if it is activated
77 : !if a full abimem.mocc report is desired, set the argument of abimem_init to "2" instead of "0"
78 : !note that abimem.mocc files can easily be multiple GB in size so don't use this option normally
79 : #ifdef HAVE_MEM_PROFILING
80 : call abimem_init(0)
81 : #endif
82 :
83 2 : if (xmpi_comm_size(comm) /= 1) then
84 0 : ABI_ERROR("fold2bloch not programmed for parallel execution.")
85 : end if
86 :
87 2 : nargs = command_argument_count()
88 :
89 2 : if (nargs == 0) then
90 2 : call prompt("Enter WFK file name:", fname)
91 2 : call prompt("Enter x y z integers giving the multiplicity:", folds)
92 : else
93 0 : call getargs(folds, fname) !Process command line arguments
94 : end if
95 : ! Use fold_matrix instead of folds(1:3) to prepare possible generalization.
96 2 : fold_matrix = 0
97 8 : do ii=1,3
98 8 : fold_matrix(ii,ii) = folds(ii)
99 : end do
100 :
101 : ! Test if the netcdf library supports MPI-IO
102 : !call nctk_test_mpiio()
103 :
104 2 : if (nctk_try_fort_or_ncfile(fname, msg) /= 0) then
105 0 : ABI_ERROR(msg)
106 : end if
107 :
108 2 : pos=INDEX(fname, "_")
109 2 : write(seedname,'(a)') fname(1:pos-1)
110 :
111 2 : write(std_out,*) ' '//achar(27)//'[97m ***********************' !print program header in pearl white
112 2 : write(std_out,*) ' ** Fold2Bloch V 1.1 **'
113 2 : write(std_out,*) ' **Build Mar 16, 2015**'
114 2 : write(std_out,*) ' ***********************'//achar(27)//'[0m'
115 :
116 2 : ebands = wfk_read_ebands(fname, xmpi_comm_self)
117 2 : iomode = iomode_from_fname(fname)
118 2 : call wfk%open_read(fname, 0, iomode, get_unit(), comm)
119 :
120 2 : nkpt=wfk%hdr%nkpt
121 2 : ABI_MALLOC(npwarr,(nkpt))
122 2 : ABI_MALLOC(nband,(nkpt))
123 2 : ABI_MALLOC(kpts,(3,nkpt))
124 :
125 2 : nsppol=wfk%hdr%nsppol
126 2 : nspinor=wfk%hdr%nspinor
127 88 : npwarr=wfk%hdr%npwarr
128 340 : kpts=wfk%hdr%kptns
129 88 : nband=wfk%hdr%nband
130 86 : mband=maxval(nband)
131 86 : mcg=maxval(npwarr)*nspinor*mband
132 8 : nfold = product(folds)
133 :
134 2 : cryst = wfk%hdr%get_crystal()
135 :
136 2 : NCF_CHECK(nctk_open_create(ncid, strcat(seedname, "_FOLD2BLOCH.nc"), xmpi_comm_self))
137 2 : fform = fform_from_ext("FOLD2BLOCH.nc")
138 2 : NCF_CHECK(wfk%hdr%ncwrite(ncid, fform, nc_define=.True.))
139 2 : NCF_CHECK(cryst%ncwrite(ncid))
140 2 : NCF_CHECK(ebands%ncwrite(ncid))
141 :
142 : ncerr = nctk_def_dims(ncid, [ &
143 : nctkdim_t("nk_unfolded", nkpt * nfold), &
144 6 : nctkdim_t("nsppol_times_nspinor", wfk%hdr%nsppol * wfk%hdr%nspinor)], defmode=.True.)
145 2 : NCF_CHECK(ncerr)
146 : ncerr = nctk_def_arrays(ncid, [ &
147 : nctkarr_t("fold_matrix", "int", "number_of_reduced_dimensions, number_of_reduced_dimensions"), &
148 : nctkarr_t("reduced_coordinates_of_unfolded_kpoints", "dp", "number_of_reduced_dimensions, nk_unfolded"), &
149 : nctkarr_t("unfolded_eigenvalues", "dp", "max_number_of_states, nk_unfolded, number_of_spins"), &
150 : nctkarr_t("spectral_weights", "dp", "max_number_of_states, nk_unfolded, nsppol_times_nspinor") &
151 10 : ])
152 2 : NCF_CHECK(ncerr)
153 2 : NCF_CHECK(nf90_inq_varid(ncid, "reduced_coordinates_of_unfolded_kpoints", kunf_varid))
154 2 : NCF_CHECK(nf90_inq_varid(ncid, "unfolded_eigenvalues", eigunf_varid))
155 2 : NCF_CHECK(nf90_inq_varid(ncid, "spectral_weights", weights_varid))
156 2 : NCF_CHECK(nctk_set_datamode(ncid))
157 2 : NCF_CHECK(nf90_put_var(ncid, nctk_idname(ncid, "fold_matrix"), fold_matrix))
158 :
159 2 : call cryst%free()
160 2 : call ebands%free()
161 :
162 4 : do csppol=1, nsppol
163 2 : if (nsppol==1) then !Determine spin polarization for output file
164 2 : outname=trim(seedname)//".f2b"
165 0 : elseif ((nsppol==2).and.(csppol==1)) then
166 0 : outname=trim(seedname)//"_UP.f2b"
167 0 : write(std_out,*) " ===================="
168 0 : write(std_out,*) " SPIN POLARIZATION UP"
169 0 : write(std_out,*) " ===================="
170 0 : elseif ((nsppol==2).and.(csppol==2)) then
171 0 : outname=trim(seedname)//"_DOWN.f2b"
172 0 : write(std_out,*) " ======================"
173 0 : write(std_out,*) " SPIN POLARIZATION DOWN"
174 0 : write(std_out,*) " ======================"
175 : end if
176 2 : if (nspinor==2) then
177 : !open output file
178 0 : if (open_file(trim(seedname)//"_SPOR_1.f2b", msg, newunit=outfile1, form="formatted", status="unknown") /= 0) then
179 0 : ABI_ERROR(msg)
180 : end if
181 0 : if (open_file(trim(seedname)//"_SPOR_2.f2b", msg, newunit=outfile2, form="formatted", status="unknown") /= 0) then
182 0 : ABI_ERROR(msg)
183 : end if
184 : else
185 2 : if (open_file(outname, msg, newunit=outfile1,form="formatted", status="unknown") /= 0) then
186 0 : ABI_ERROR(msg)
187 : end if
188 : end if
189 :
190 86 : do ikpt=1, nkpt !For each K point
191 84 : ABI_MALLOC(cg,(2,mcg))
192 84 : ABI_MALLOC(eig,((2*mband)**0*mband))
193 84 : ABI_MALLOC(kg,(3,npwarr(ikpt)))
194 84 : ABI_MALLOC(coefc,(2,nspinor*npwarr(ikpt)))
195 84 : ABI_MALLOC(weights, (nfold))
196 84 : ABI_MALLOC(nkval,(3, nfold))
197 84 : call progress(ikpt,nkpt,kpts(:,ikpt)) !Write progress information
198 :
199 : !Read a block of data
200 252 : call wfk%read_band_block([1, nband(ikpt)], ikpt, csppol, xmpio_single, kg_k=kg, cg_k=cg, eig_k=eig)
201 :
202 : !Determine unfolded K point states
203 84 : call newk(kpts(1,ikpt),kpts(2,ikpt),kpts(3,ikpt),folds(1),folds(2),folds(3),nkval)
204 84 : if (csppol == 1) then
205 420 : NCF_CHECK(nf90_put_var(ncid, kunf_varid, nkval, start=[1, 1 + (ikpt-1) * nfold], count=[3, nfold]))
206 : end if
207 :
208 84 : cg_b=1
209 420 : do iband=1, nband(ikpt) !Foe each Eigenvalue
210 648264 : coefc=cg(:,cg_b:(cg_b+nspinor*npwarr(ikpt)-1)) !Split coefficients per eigen value according to the number of "kg"
211 672 : do cspinor=1,nspinor
212 336 : if (cspinor==1) then
213 336 : outfile=outfile1
214 336 : lwcg=1
215 336 : hicg=npwarr(ikpt)
216 : else
217 : ! Move coefficient span to spinor 2
218 0 : outfile=outfile2
219 0 : lwcg=npwarr(ikpt)+1
220 0 : hicg=npwarr(ikpt)*nspinor
221 : end if
222 336 : call sortc(folds(1),folds(2),folds(3),kg,coefc(:,lwcg:hicg),npwarr(ikpt),weights)
223 : ! Write out results, format: new k states(x, y, and z), eigenvalue, weight
224 2352 : do count=1, nfold
225 2352 : write(outfile,50) nkval(1,count),nkval(2,count),nkval(3,count),eig(iband),weights(count)
226 : 50 format(f11.6, f11.6, f11.6, f11.6, f11.6)
227 : end do
228 336 : iss = csppol; if (nspinor == 2) iss = cspinor
229 : ncerr = nf90_put_var(ncid, weights_varid, weights, start=[iband, 1 + (ikpt-1) * nfold, iss], &
230 3360 : stride=[mband, 1, 1], count=[1, nfold, 1])
231 336 : NCF_CHECK(ncerr)
232 672 : if (cspinor == 1) then
233 2352 : weights = eig(iband) ! Use weights as workspace array.
234 : ncerr = nf90_put_var(ncid, eigunf_varid, weights, start=[iband, 1 + (ikpt-1) * nfold, csppol], &
235 3360 : stride=[mband, 1, 1], count=[1, nfold, 1])
236 : !count=[1, nfold, 1])
237 336 : NCF_CHECK(ncerr)
238 : end if
239 : end do ! cspinor
240 420 : cg_b=cg_b+nspinor*npwarr(ikpt) !shift coefficient pointer for next eigenvalue
241 : end do ! iband
242 :
243 84 : ABI_FREE(cg)
244 84 : ABI_FREE(eig)
245 84 : ABI_FREE(kg)
246 84 : ABI_FREE(coefc)
247 84 : ABI_FREE(weights)
248 86 : ABI_FREE(nkval)
249 : end do
250 4 : if (nspinor==2) then
251 0 : close(outfile1) !close output file
252 0 : close(outfile2)
253 : else
254 2 : close(outfile1)
255 : end if
256 : end do
257 2 : call wfk%close()
258 :
259 2 : ABI_FREE(kpts)
260 2 : ABI_FREE(nband)
261 2 : ABI_FREE(npwarr)
262 :
263 : ! Print summary
264 2 : write(std_out,*) ' '//achar(27)//'[97m Number of K points processed:', nkpt
265 2 : if (nsppol==2) then
266 0 : write(std_out,*) ' Data was written to: ', trim(seedname)//"_UP.f2b", " & ", trim(seedname)//"_DOWN.f2b"
267 : else
268 2 : if (nspinor==2) then
269 0 : write(std_out,*) ' Data was written to: ', trim(seedname)//"_SPOR_1.f2b", " & ", trim(seedname)//"_SPOR_2.f2b"
270 : else
271 2 : write(std_out,*) ' Data was written to: ', trim(seedname)//".f2b"
272 : end if
273 : end if
274 2 : write(std_out,*) ' Data format: KX, KY, KZ, Eigenvalue(Ha), Weight'//achar(27)//'[0m'
275 :
276 2 : NCF_CHECK(nf90_close(ncid))
277 :
278 : !Write information on file about the memory before ending mpi module, if memory profiling is enabled
279 2 : call abinit_doctor("__fold2bloch")
280 :
281 2 : call xmpi_end()
282 :
283 0 : end program fold2Bloch
284 : !!***
|