TABLE OF CONTENTS


ABINIT/fold2Bloch [ Programs ]

[ Top ] [ Programs ]

NAME

 fold2Bloch

FUNCTION

 Main routine for the unfolding of the wavefuntion.

COPYRIGHT

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

INPUTS

  (main program)

OUTPUT

  (main program)

NOTES

 folds= Array of folds in X,Y, and Z directions

SOURCE

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