TABLE OF CONTENTS


ABINIT/m_scf_history [ Modules ]

[ Top ] [ Modules ]

NAME

  m_scf_history

FUNCTION

  This module provides the definition of the scf_history_type used to store
  various arrays obtained from previous SCF cycles (density, positions, wavefunctions ...),
  as needed by the specific SCF algorithm.

COPYRIGHT

 Copyright (C) 2011-2022 ABINIT group (MT)
 This file is distributed under the terms of the
 GNU General Public License, see ~abinit/COPYING
 or http://www.gnu.org/copyleft/gpl.txt .

INPUTS

OUTPUT

SOURCE

22 #if defined HAVE_CONFIG_H
23 #include "config.h"
24 #endif
25 
26 #include "abi_common.h"
27 
28 MODULE m_scf_history
29 
30  use defs_basis
31  use m_abicore
32  use m_dtset
33  use m_errors
34 
35  use defs_abitypes, only : MPI_type
36  use m_pawcprj,  only : pawcprj_type, pawcprj_free
37  use m_pawrhoij, only : pawrhoij_type, pawrhoij_nullify, pawrhoij_free
38 
39  implicit none
40 
41  private
42 
43 ! public procedures.
44  public :: scf_history_init
45  public :: scf_history_free
46  public :: scf_history_nullify

m_scf_history/scf_history_free [ Functions ]

[ Top ] [ m_scf_history ] [ Functions ]

NAME

  scf_history_free

FUNCTION

  Clean and destroy a scf_history datastructure

SIDE EFFECTS

  scf_history(:)=<type(scf_history_type)>=scf_history datastructure

SOURCE

336 subroutine scf_history_free(scf_history)
337 
338 !Arguments ------------------------------------
339 !arrays
340  type(scf_history_type),intent(inout) :: scf_history
341 
342 !Local variables-------------------------------
343 !scalars
344  integer :: jj
345 
346 !************************************************************************
347 
348  !@scf_history_type
349 
350  if (allocated(scf_history%pawrhoij_last)) then
351    call pawrhoij_free(scf_history%pawrhoij_last)
352    ABI_FREE(scf_history%pawrhoij_last)
353  end if
354  if (allocated(scf_history%pawrhoij)) then
355    do jj=1,size(scf_history%pawrhoij,2)
356      call pawrhoij_free(scf_history%pawrhoij(:,jj))
357    end do
358    ABI_FREE(scf_history%pawrhoij)
359  end if
360  if (allocated(scf_history%cprj)) then
361    do jj=1,size(scf_history%cprj,3)
362      call pawcprj_free(scf_history%cprj(:,:,jj))
363    end do
364    ABI_FREE(scf_history%cprj)
365  end if
366 
367  if (allocated(scf_history%hindex))       then
368    ABI_FREE(scf_history%hindex)
369  end if
370  if (allocated(scf_history%deltarhor))    then
371    ABI_FREE(scf_history%deltarhor)
372  end if
373  if (allocated(scf_history%xreddiff))     then
374    ABI_FREE(scf_history%xreddiff)
375  end if
376  if (allocated(scf_history%atmrho_last))  then
377    ABI_FREE(scf_history%atmrho_last)
378  end if
379  if (allocated(scf_history%xred_last))    then
380    ABI_FREE(scf_history%xred_last)
381  end if
382  if (allocated(scf_history%rhor_last))    then
383    ABI_FREE(scf_history%rhor_last)
384  end if
385  if (allocated(scf_history%taur_last))    then
386    ABI_FREE(scf_history%taur_last)
387  end if
388  if (allocated(scf_history%cg))           then
389    ABI_FREE(scf_history%cg)
390  end if
391  if (allocated(scf_history%eigen))           then
392    ABI_FREE(scf_history%eigen)
393  end if
394  if (allocated(scf_history%dotprod_sumdiag_cgcprj_ij))           then
395    ABI_FREE(scf_history%dotprod_sumdiag_cgcprj_ij)
396  end if
397 
398  scf_history%history_size=-1
399  scf_history%usecg=0
400  scf_history%icall=0
401  scf_history%mcprj=0
402  scf_history%mcg=0
403  scf_history%meigen=0
404 
405 end subroutine scf_history_free

m_scf_history/scf_history_init [ Functions ]

[ Top ] [ m_scf_history ] [ Functions ]

NAME

  scf_history_init

FUNCTION

  Init all scalars and pointers in a scf_history datastructure
  according to scf_history%history_size value which has to be
  defined before calling this routine

INPUTS

  dtset <type(dataset_type)>=all input variables in this dataset
  mpi_enreg=MPI-parallelisation information
  usecg= if ==0 => no handling of wfs (and eigenvalues),
         if==1 => handling of density/potential AND wfs and eigen,
         if==2 => ONLY handling of wfs and eigen

SIDE EFFECTS

  scf_history=<type(scf_history_type)>=scf_history datastructure
    hindex is always allocated
    The density/potential arrays that are possibly allocated are : atmrho_last, deltarhor,
      pawrhoij, pawrhoij_last, rhor_last, taur_last, xreddiff, xred_last.
    The wfs arrays that are possibly allocated are : cg, cprj and eigen

SOURCE

216 subroutine scf_history_init(dtset,mpi_enreg,usecg,scf_history)
217 
218 !Arguments ------------------------------------
219 !scalars
220  integer, intent(in) :: usecg
221  type(dataset_type),intent(in) :: dtset
222  type(MPI_type),intent(in) :: mpi_enreg
223 !arrays
224  type(scf_history_type),intent(inout) :: scf_history
225 
226 !Local variables-------------------------------
227 !scalars
228  integer :: jj,mband_cprj,my_natom,my_nspinor,nfft
229 !arrays
230 
231 !************************************************************************
232 
233  !@scf_history_type
234 
235  if (scf_history%history_size<0) then
236    call scf_history_nullify(scf_history)
237  else
238 
239    scf_history%usecg=usecg
240    scf_history%wfmixalg=dtset%fockoptmix/100
241 
242    nfft=dtset%nfft
243    if (dtset%usepaw==1.and.(dtset%pawecutdg>=1.0000001_dp*dtset%ecut)) nfft=dtset%nfftdg
244    my_natom=mpi_enreg%my_natom
245 
246    if (scf_history%history_size>=0 .and. usecg<2) then
247      ABI_MALLOC(scf_history%rhor_last,(nfft,dtset%nspden))
248      ABI_MALLOC(scf_history%taur_last,(nfft,dtset%nspden*dtset%usekden))
249      ABI_MALLOC(scf_history%xred_last,(3,dtset%natom))
250      ABI_MALLOC(scf_history%pawrhoij_last,(my_natom*dtset%usepaw))
251      if (dtset%usepaw==1) then
252        call pawrhoij_nullify(scf_history%pawrhoij_last)
253      end if
254    end if
255 
256    if (scf_history%history_size>0) then
257 
258      scf_history%natom=dtset%natom
259      scf_history%nfft=nfft
260      scf_history%nspden=dtset%nspden
261      scf_history%beta=zero
262      scf_history%icall=0
263 
264      scf_history%mcg=0
265      scf_history%mcprj=0
266      scf_history%meigen=0
267      if (usecg>0) then
268        my_nspinor=max(1,dtset%nspinor/mpi_enreg%nproc_spinor)
269        scf_history%meigen=dtset%nbandhf*dtset%mkmem*dtset%nsppol
270        scf_history%mcg=dtset%mpw*my_nspinor*scf_history%meigen ! This is for scf_history_wf
271        if(usecg==1)then
272          scf_history%meigen=dtset%mband*dtset%mkmem*dtset%nsppol
273          scf_history%mcg=dtset%mpw*my_nspinor*scf_history%meigen ! This is for scf_history (when extrapwf==1)
274        endif
275        if (dtset%usepaw==1) then
276          mband_cprj=dtset%nbandhf
277          if(usecg==1)mband_cprj=dtset%mband
278          if (dtset%paral_kgb/=0) mband_cprj=mband_cprj/mpi_enreg%nproc_band
279          scf_history%mcprj=my_nspinor*mband_cprj*dtset%mkmem*dtset%nsppol
280        end if
281      end if
282 
283      if (usecg<2) then
284        ABI_MALLOC(scf_history%hindex,(scf_history%history_size))
285        scf_history%alpha=zero
286      else
287        ABI_MALLOC(scf_history%hindex,(2*scf_history%history_size+1))
288        scf_history%alpha=dtset%wfmix
289      endif
290      scf_history%hindex(:)=0
291 
292      if (usecg<2) then
293        ABI_MALLOC(scf_history%deltarhor,(nfft,dtset%nspden,scf_history%history_size))
294        ABI_MALLOC(scf_history%xreddiff,(3,dtset%natom,scf_history%history_size))
295        ABI_MALLOC(scf_history%atmrho_last,(nfft))
296        if (dtset%usepaw==1) then
297          ABI_MALLOC(scf_history%pawrhoij,(my_natom,scf_history%history_size))
298          do jj=1,scf_history%history_size
299            call pawrhoij_nullify(scf_history%pawrhoij(:,jj))
300          end do
301        endif
302      end if
303 
304      if (scf_history%usecg>0) then
305        ABI_MALLOC(scf_history%cg,(2,scf_history%mcg,scf_history%history_size))
306        ABI_MALLOC(scf_history%eigen,(scf_history%meigen,scf_history%history_size))
307 !      Note that the allocation is made even when usepaw==0. Still, scf_history%mcprj=0 ...
308        ABI_MALLOC(scf_history%cprj,(dtset%natom,scf_history%mcprj,scf_history%history_size))
309      end if
310 
311      if (scf_history%usecg==2)then
312 !      This relatively small matrix is always allocated when usecg==1, even if not used
313        ABI_MALLOC(scf_history%dotprod_sumdiag_cgcprj_ij,(2,scf_history%history_size,scf_history%history_size))
314      endif
315 
316    end if
317  end if
318 
319 end subroutine scf_history_init

m_scf_history/scf_history_nullify [ Functions ]

[ Top ] [ m_scf_history ] [ Functions ]

NAME

  scf_history_nullify

FUNCTION

  Nullify (set to null) an scf_history datastructure

SIDE EFFECTS

  scf_history(:)=<type(scf_history_type)>=scf_history datastructure

SOURCE

422 subroutine scf_history_nullify(scf_history)
423 
424 !Arguments ------------------------------------
425 !arrays
426  type(scf_history_type),intent(inout) :: scf_history
427 !Local variables-------------------------------
428 !scalars
429 
430 !************************************************************************
431 
432  !@scf_history_type
433  scf_history%history_size=-1
434  scf_history%icall=0
435  scf_history%mcprj=0
436  scf_history%mcg=0
437  scf_history%meigen=0
438 
439 end subroutine scf_history_nullify

m_scf_history/scf_history_type [ Types ]

[ Top ] [ m_scf_history ] [ Types ]

NAME

 scf_history_type

FUNCTION

 This structured datatype contains various arrays obtained from
 previous SCF cycles (density, positions...)

SOURCE

 59  type, public :: scf_history_type
 60 
 61 ! WARNING : if you modify this datatype, please check whether there might be creation/destruction/copy routines,
 62 ! declared in another part of ABINIT, that might need to take into account your modification.
 63 
 64 ! Integer scalar
 65 
 66   integer :: history_size
 67    ! Number of previous SCF cycles stored in history
 68    ! If history_size<0, scf_history is not used
 69    ! If history_size=0, scf_history only contains
 70    !    current values of data (rhor, taur, pawrhoih, xred)
 71    ! If history_size>0, scf_history contains
 72    !    current values of data and also
 73    !    history_size previous values of these data
 74 
 75   integer :: icall
 76    ! Number of call for the routine extraprho or wf_mixing
 77 
 78   integer :: mcg
 79    ! Size of cg array
 80 
 81   integer :: mcprj
 82    ! Size of cprj datsatructure array
 83 
 84   integer :: meigen
 85    ! Size of eigen array
 86 
 87   integer :: natom
 88    ! Number of atoms in cell
 89 
 90   integer :: nfft
 91    ! Size of FFT grid (for density)
 92 
 93   integer :: nspden
 94    ! Number of independant spin components for density
 95 
 96   integer :: usecg
 97    ! usecg=0 if the extrapolation/mixing of the density/potential is active but not the one of the wavefunction
 98    ! usecg=1 if the extrapolation/mixing of the density/potential and wavefunctions is active
 99    ! usecg=2 if the extrapolation/mixing of the wavefunctions is active but not the one of the density/potential
100 
101   integer :: wfmixalg
102    ! algorithm used to mix the wavefunctions (in case usecg=2)
103 
104   real(dp) :: alpha
105    ! alpha mixing coefficient for the prediction of density and wavefunctions
106    ! In the case of wavefunction simple mixing, contain wfmix factor
107 
108   real(dp) :: beta
109    ! beta mixing coefficient for the prediction of density and wavefunctions
110 
111 ! Integer arrays
112 
113   integer,allocatable :: hindex(:)
114    ! Indexes of SCF cycles in the history
115    !
116    ! For the density-based schemes (with or without wavefunctions) :
117    ! hindex(history_size)
118    ! hindex(1) is the newest SCF cycle
119    ! hindex(history_size) is the oldest SCF cycle
120    !
121    ! For wavefunction-based schemes (outer loop of a double loop SCF):
122    ! hindex(2*history_size+1)
123    ! The odd indices refer to the out wavefunction,
124    ! the even indices refer to the in wavefunction (not all such wavefunctions being stored, though).
125    ! hindex(1:2) is the newest SCF cycle, hindex(3:4) is the SCF cycle before the newest one ... In case of an
126    ! algorithm based on a biorthogonal ensemble of wavefunctions, the reference is stored in hindex(2*history_size+1)
127    ! When the index points to a location beyond history_size, the corresponding wavefunction set must be reconstructed
128    ! from the existing wavefunctions sets (to be implemented)
129 
130 ! Real (real(dp)) arrays
131 
132    real(dp),allocatable :: cg(:,:,:)
133     ! cg(2,mcg,history_size)
134     ! wavefunction coefficients needed for each SCF cycle of history
135     ! Might also contain the wf residuals
136 
137    real(dp),allocatable :: deltarhor(:,:,:)
138     ! deltarhor(nfft,nspden,history_size)
139     ! Difference between electronic density (in real space)
140     ! and sum of atomic densities at the end of each SCF cycle of history
141 
142    real(dp),allocatable :: eigen(:,:)
143     ! eigen(meigen,history_size)
144     ! eigenvalues for each SCF cycle of history
145 
146    real(dp),allocatable :: atmrho_last(:)
147     ! atmrho_last(nfft)
148     ! Sum of atomic densities at the end of the LAST SCF cycle
149 
150    real(dp),allocatable :: rhor_last(:,:)
151     ! rhor_last(nfft,nspden)
152     ! Last computed electronic density (in real space)
153 
154    real(dp),allocatable :: taur_last(:,:)
155     ! taur_last(nfft,nspden*usekden)
156     ! Last computed kinetic energy density (in real space)
157 
158    real(dp),allocatable :: xreddiff(:,:,:)
159     ! xreddiff(3,natom,history_size)
160     ! Difference of reduced coordinates of atoms between a
161     ! SCF cycle and the previous
162 
163    real(dp),allocatable :: xred_last(:,:)
164     ! xred_last(3,natom)
165     ! Last computed atomic positions (reduced coordinates)
166 
167    real(dp),allocatable :: dotprod_sumdiag_cgcprj_ij(:,:,:)
168     ! dotprod_sumdiag_cgcprj_mn(2,history_size,history_size)
169     ! Container for the scalar products between aligned sets of wavefunctions or their residuals
170     ! S_ij=Sum_nk <wf_nk(set i)|wf_nk(set j)> possibly with some weighting factor that might depend on nk.
171 
172 ! Structured datatypes arrays
173 
174   type(pawrhoij_type), allocatable :: pawrhoij(:,:)
175     ! pawrhoij(natom,history_size)
176     ! PAW only: occupancies matrix at the end of each SCF cycle of history
177 
178   type(pawrhoij_type), allocatable :: pawrhoij_last(:)
179     ! pawrhoij_last(natom)
180     ! PAW only: last computed occupancies matrix
181 
182   type(pawcprj_type),allocatable :: cprj(:,:,:)
183     !cprj(natom,nspinor*mband*mkmem*nsppol,history_size)
184 
185  end type scf_history_type