LCOV - code coverage report
Current view: top level - src/98_main - fold2Bloch.F90 (source / functions) Coverage Total Hit
Test: coverage.info Lines: 81.3 % 139 113
Test Date: 2026-11-06 14:42:53 Functions: 100.0 % 2 2
Branches: 52.2 % 414 216

             Branch data     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                 :             : !!***
        

Generated by: LCOV version 2.3-1