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-2018 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

PARENTS

CHILDREN

SOURCE

27 #if defined HAVE_CONFIG_H
28 #include "config.h"
29 #endif
30 
31 #include "abi_common.h"
32 
33 MODULE m_scf_history
34 
35  use defs_basis
36  use defs_abitypes
37  use m_abicore
38  use m_errors
39 
40  use m_pawcprj,  only : pawcprj_type, pawcprj_free
41  use m_pawrhoij, only : pawrhoij_type, pawrhoij_nullify, pawrhoij_free
42 
43  implicit none
44 
45  private
46 
47 ! public procedures.
48  public :: scf_history_init
49  public :: scf_history_free
50  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

PARENTS

      gstateimg,scfcv

CHILDREN

SOURCE

359 subroutine scf_history_free(scf_history)
360 
361 
362 !This section has been created automatically by the script Abilint (TD).
363 !Do not modify the following lines by hand.
364 #undef ABI_FUNC
365 #define ABI_FUNC 'scf_history_free'
366 !End of the abilint section
367 
368  implicit none
369 
370 !Arguments ------------------------------------
371 !arrays
372  type(scf_history_type),intent(inout) :: scf_history
373 
374 !Local variables-------------------------------
375 !scalars
376  integer :: jj
377 
378 !************************************************************************
379 
380  !@scf_history_type
381 
382  if (allocated(scf_history%pawrhoij_last)) then
383    call pawrhoij_free(scf_history%pawrhoij_last)
384    ABI_DATATYPE_DEALLOCATE(scf_history%pawrhoij_last)
385  end if
386  if (allocated(scf_history%pawrhoij)) then
387    do jj=1,size(scf_history%pawrhoij,2)
388      call pawrhoij_free(scf_history%pawrhoij(:,jj))
389    end do
390    ABI_DATATYPE_DEALLOCATE(scf_history%pawrhoij)
391  end if
392  if (allocated(scf_history%cprj)) then
393    do jj=1,size(scf_history%cprj,3)
394      call pawcprj_free(scf_history%cprj(:,:,jj))
395    end do
396    ABI_DATATYPE_DEALLOCATE(scf_history%cprj)
397  end if
398 
399  if (allocated(scf_history%hindex))       then
400    ABI_DEALLOCATE(scf_history%hindex)
401  end if
402  if (allocated(scf_history%deltarhor))    then
403    ABI_DEALLOCATE(scf_history%deltarhor)
404  end if
405  if (allocated(scf_history%xreddiff))     then
406    ABI_DEALLOCATE(scf_history%xreddiff)
407  end if
408  if (allocated(scf_history%atmrho_last))  then
409    ABI_DEALLOCATE(scf_history%atmrho_last)
410  end if
411  if (allocated(scf_history%xred_last))    then
412    ABI_DEALLOCATE(scf_history%xred_last)
413  end if
414  if (allocated(scf_history%rhor_last))    then
415    ABI_DEALLOCATE(scf_history%rhor_last)
416  end if
417  if (allocated(scf_history%taur_last))    then
418    ABI_DEALLOCATE(scf_history%taur_last)
419  end if
420  if (allocated(scf_history%cg))           then
421    ABI_DEALLOCATE(scf_history%cg)
422  end if
423  if (allocated(scf_history%eigen))           then
424    ABI_DEALLOCATE(scf_history%eigen)
425  end if
426  if (allocated(scf_history%dotprod_sumdiag_cgcprj_ij))           then
427    ABI_DEALLOCATE(scf_history%dotprod_sumdiag_cgcprj_ij)
428  end if
429 
430  scf_history%history_size=-1
431  scf_history%usecg=0
432  scf_history%icall=0
433  scf_history%mcprj=0
434  scf_history%mcg=0
435  scf_history%meigen=0
436 
437 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

PARENTS

      gstate,scfcv

CHILDREN

SOURCE

225 subroutine scf_history_init(dtset,mpi_enreg,usecg,scf_history)
226 
227 
228 !This section has been created automatically by the script Abilint (TD).
229 !Do not modify the following lines by hand.
230 #undef ABI_FUNC
231 #define ABI_FUNC 'scf_history_init'
232 !End of the abilint section
233 
234  implicit none
235 
236 !Arguments ------------------------------------
237 !scalars
238  integer, intent(in) :: usecg
239  type(dataset_type),intent(in) :: dtset
240  type(MPI_type),intent(in) :: mpi_enreg
241 !arrays
242  type(scf_history_type),intent(inout) :: scf_history
243 
244 !Local variables-------------------------------
245 !scalars
246  integer :: jj,mband_cprj,my_natom,my_nspinor,nfft
247 !arrays
248 
249 !************************************************************************
250 
251  !@scf_history_type
252 
253  if (scf_history%history_size<0) then
254    call scf_history_nullify(scf_history)
255  else
256 
257    scf_history%usecg=usecg
258    scf_history%wfmixalg=dtset%fockoptmix/100
259 
260    nfft=dtset%nfft
261    if (dtset%usepaw==1.and.(dtset%pawecutdg>=1.0000001_dp*dtset%ecut)) nfft=dtset%nfftdg
262    my_natom=mpi_enreg%my_natom
263 
264    if (scf_history%history_size>=0 .and. usecg<2) then
265      ABI_ALLOCATE(scf_history%rhor_last,(nfft,dtset%nspden))
266      ABI_ALLOCATE(scf_history%taur_last,(nfft,dtset%nspden*dtset%usekden))
267      ABI_ALLOCATE(scf_history%xred_last,(3,dtset%natom))
268      ABI_DATATYPE_ALLOCATE(scf_history%pawrhoij_last,(my_natom*dtset%usepaw))
269      if (dtset%usepaw==1) then
270        call pawrhoij_nullify(scf_history%pawrhoij_last)
271      end if
272    end if
273 
274    if (scf_history%history_size>0) then
275 
276      scf_history%natom=dtset%natom
277      scf_history%nfft=nfft
278      scf_history%nspden=dtset%nspden
279      scf_history%beta=zero
280      scf_history%icall=0
281 
282      scf_history%mcg=0
283      scf_history%mcprj=0
284      scf_history%meigen=0
285      if (usecg>0) then
286        my_nspinor=max(1,dtset%nspinor/mpi_enreg%nproc_spinor)
287        scf_history%meigen=dtset%nbandhf*dtset%mkmem*dtset%nsppol
288        scf_history%mcg=dtset%mpw*my_nspinor*scf_history%meigen ! This is for scf_history_wf
289        if(usecg==1)then
290          scf_history%meigen=dtset%mband*dtset%mkmem*dtset%nsppol
291          scf_history%mcg=dtset%mpw*my_nspinor*scf_history%meigen ! This is for scf_history (when extrapwf==1)
292        endif
293        if (dtset%usepaw==1) then
294          mband_cprj=dtset%nbandhf
295          if(usecg==1)mband_cprj=dtset%mband
296          if (dtset%paral_kgb/=0) mband_cprj=mband_cprj/mpi_enreg%nproc_band
297          scf_history%mcprj=my_nspinor*mband_cprj*dtset%mkmem*dtset%nsppol
298        end if
299      end if
300 
301      if (usecg<2) then
302        ABI_ALLOCATE(scf_history%hindex,(scf_history%history_size))
303        scf_history%alpha=zero
304      else
305        ABI_ALLOCATE(scf_history%hindex,(2*scf_history%history_size+1))
306        scf_history%alpha=dtset%wfmix
307      endif
308      scf_history%hindex(:)=0
309 
310      if (usecg<2) then 
311        ABI_ALLOCATE(scf_history%deltarhor,(nfft,dtset%nspden,scf_history%history_size))
312        ABI_ALLOCATE(scf_history%xreddiff,(3,dtset%natom,scf_history%history_size))
313        ABI_ALLOCATE(scf_history%atmrho_last,(nfft))
314        if (dtset%usepaw==1) then
315          ABI_DATATYPE_ALLOCATE(scf_history%pawrhoij,(my_natom,scf_history%history_size))
316          do jj=1,scf_history%history_size
317            call pawrhoij_nullify(scf_history%pawrhoij(:,jj))
318          end do
319        endif
320      end if
321 
322      if (scf_history%usecg>0) then
323        ABI_ALLOCATE(scf_history%cg,(2,scf_history%mcg,scf_history%history_size))
324        ABI_ALLOCATE(scf_history%eigen,(scf_history%meigen,scf_history%history_size))
325 !      Note that the allocation is made even when usepaw==0. Still, scf_history%mcprj=0 ...
326        ABI_DATATYPE_ALLOCATE(scf_history%cprj,(dtset%natom,scf_history%mcprj,scf_history%history_size))
327      end if
328 
329      if (scf_history%usecg==2)then
330 !      This relatively small matrix is always allocated when usecg==1, even if not used
331        ABI_ALLOCATE(scf_history%dotprod_sumdiag_cgcprj_ij,(2,scf_history%history_size,scf_history%history_size))
332      endif
333 
334    end if
335  end if
336 
337 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

PARENTS

      gstateimg,m_scf_history

CHILDREN

SOURCE

459 subroutine scf_history_nullify(scf_history)
460 
461 
462 !This section has been created automatically by the script Abilint (TD).
463 !Do not modify the following lines by hand.
464 #undef ABI_FUNC
465 #define ABI_FUNC 'scf_history_nullify'
466 !End of the abilint section
467 
468  implicit none
469 
470 !Arguments ------------------------------------
471 !arrays
472  type(scf_history_type),intent(inout) :: scf_history
473 !Local variables-------------------------------
474 !scalars
475 
476 !************************************************************************
477 
478  !@scf_history_type
479  scf_history%history_size=-1
480  scf_history%icall=0
481  scf_history%mcprj=0
482  scf_history%mcg=0
483  scf_history%meigen=0
484 
485 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

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