Branch data Line data Source code
1 : : !!****m* ABINIT/m_rttddft_tdks
2 : : !! NAME
3 : : !! m_rttddft_tdks
4 : : !!
5 : : !! FUNCTION
6 : : !! Contains the main object (tdks) to propagate
7 : : !! the time-dependent Kohn-Sham equations in RT-TDDFT
8 : : !!
9 : : !! COPYRIGHT
10 : : !! Copyright (C) 2021-2026 ABINIT group (FB)
11 : : !! This file is distributed under the terms of the
12 : : !! GNU General Public License, see ~abinit/COPYING
13 : : !! or http://www.gnu.org/copyleft/gpl.txt .
14 : : !!
15 : : !! SOURCE
16 : :
17 : : #if defined HAVE_CONFIG_H
18 : : #include "config.h"
19 : : #endif
20 : :
21 : : #include "abi_common.h"
22 : :
23 : : module m_rttddft_tdks
24 : :
25 : : use defs_basis
26 : : use defs_abitypes, only: MPI_type
27 : : use defs_datatypes, only: pseudopotential_type
28 : : use defs_wvltypes, only: wvl_data, nullify_wvl_data
29 : : use libxc_functionals, only: libxc_functionals_get_hybridparams
30 : : use m_bandfft_kpt, only: bandfft_kpt, bandfft_kpt_init1, bandfft_kpt_destroy_array
31 : : use m_cgprj, only: ctocprj
32 : : use m_common, only: setup1
33 : : use m_dtfil, only: datafiles_type
34 : : use m_dtset, only: dataset_type
35 : : use m_ebands, only: ebands_t, unpack_eneocc
36 : : use m_energies, only: energies_type
37 : : use m_errors, only: msg_hndl, assert
38 : : use m_extfpmd, only: extfpmd_type
39 : : use m_gemm_nonlop_projectors, only: init_gemm_nonlop, destroy_gemm_nonlop
40 : : use m_geometry, only: fixsym
41 : : use m_hdr, only: hdr_type
42 : : use m_initylmg, only: initylmg
43 : : use m_invovl, only: init_invovl, destroy_invovl
44 : : use m_io_tools, only: open_file
45 : : use m_inwffil, only: inwffil
46 : : use m_kg, only: kpgio, getph, getcut
47 : : use m_mpinfo, only: proc_distrb_cycle
48 : : use m_occ, only: newocc
49 : : use m_paw_an, only: paw_an_type, paw_an_init, paw_an_free, &
50 : : paw_an_nullify
51 : : use m_pawang, only: pawang_type
52 : : use m_pawcprj, only: pawcprj_type,pawcprj_free,pawcprj_alloc, &
53 : : pawcprj_getdim
54 : : use m_paw_dmft, only: init_sc_dmft,destroy_sc_dmft,paw_dmft_type
55 : : use m_pawfgr, only: pawfgr_type, pawfgr_init, pawfgr_destroy
56 : : use m_pawfgrtab, only: pawfgrtab_type, pawfgrtab_init, pawfgrtab_free
57 : : use m_paw_init, only: pawinit,paw_gencond
58 : : use m_paw_ij, only: paw_ij_type, paw_ij_init, paw_ij_free, paw_ij_nullify
59 : : use m_paw_nhat, only: nhatgrid
60 : : use m_paw_occupancies, only: initrhoij
61 : : use m_pawrad, only: pawrad_type
62 : : use m_pawrhoij, only: pawrhoij_type, pawrhoij_copy, pawrhoij_free
63 : : use m_paw_sphharm, only: setsym_ylm
64 : : use m_pawtab, only: pawtab_type, pawtab_get_lsize
65 : : use m_paw_tools, only: chkpawovlp
66 : : use m_pawxc, only: pawxc_get_usekden
67 : : use m_pspini, only: pspini
68 : : use m_profiling_abi, only: abimem_record
69 : : use m_rttddft_tdef, only: tdef_type
70 : : use m_spacepar, only: setsym
71 : : use m_specialmsg, only: wrtout
72 : : use m_symtk, only: symmetrize_xred
73 : : use m_wffile, only: wffile_type, WffClose
74 : : use m_xmpi, only: xmpi_bcast, xmpi_sum
75 : : use m_drivexc, only: xc_need_kden
76 : :
77 : : implicit none
78 : :
79 : : private
80 : : !!***
81 : :
82 : : !! NAME
83 : : !! tdks_type: Time Dependent Kohn-Sham type
84 : : !! Object containing the TD KS orbitals and all other
85 : : !! important variables required to run RT-TDDFT
86 : : !!
87 : : !! SOURCE
88 : : type,public :: tdks_type
89 : :
90 : : !scalars
91 : : integer :: bantot !total number of bands
92 : : integer :: first_step !start propagation from first_step (for restart)
93 : : integer :: mband_cprj !nb of band per proc (for cprj)
94 : : integer :: mcg !nb of WFs (cg) coeffs
95 : : integer :: mcprj !nb of cprj (projectors applied to WF)
96 : : integer :: nfftf !nb of FFT grid pts (fine grid)
97 : : integer :: nfft !nb of FFT grid pts (coarse grid)
98 : : integer :: nhatgrdim !dimension of nhatgr array
99 : : integer :: ngrvdw !dimension of grvdw array
100 : : integer :: ntime !max nb of time steps
101 : : integer :: current_unit!unit nb of the current density file
102 : : integer :: tdener_unit !unit nb of the energy file
103 : : integer :: tdef_unit !unit nb of the efield file
104 : : integer :: tdrestart_unit !unit nb of the restart file
105 : : integer :: unpaw !paw data tmp file unit
106 : : integer :: usexcnhat !use nhat in the computation of the XC term
107 : : real(dp) :: dt !propagation time step
108 : : real(dp) :: ecore !core energy
109 : : real(dp) :: etot !total energy
110 : : real(dp) :: gsqcut !cut-off on G^2
111 : : real(dp) :: ucvol !primitive cell volume
112 : : real(dp) :: zion !total ionic charge
113 : : logical :: gemm_nonlop_use_gemm !use efficient BLAS call
114 : : !for computing non local potential
115 : : type(energies_type) :: energies !contains various energy values
116 : : type(hdr_type) :: hdr !header: contains various info
117 : : type(paw_dmft_type) :: paw_dmft !paw_dmft object (unused but
118 : : !required by various routines)
119 : : type(pawfgr_type) :: pawfgr !FFT fine grid in PAW sphere
120 : : type(pawang_type),pointer :: pawang => NULL() !angular grid in PAW sphere
121 : : type(tdef_type) :: tdef !Object containing variables related to TD electric field
122 : : type(wvl_data) :: wvl !wavelets ojects (unused but
123 : : !required by various routines)
124 : : character(len=fnlen) :: fname_current!Name of the TDCURRENT file
125 : : character(len=fnlen) :: fname_tdener!Name of the TDENER file
126 : : character(len=fnlen) :: fname_tdef !Name of the TDEFIELD file
127 : : character(len=fnlen) :: fname_wfk0 !Name of the input WFK file containing
128 : : !the intial (t=0) wfs
129 : : !arrays
130 : : integer,allocatable :: atindx(:) !index table of atom ordered by type
131 : : integer,allocatable :: atindx1(:) !nb of the atom for each index in atindx
132 : : integer,allocatable :: dimcprj(:) !Contains dimension for cprj array
133 : : integer,allocatable :: dimcprj_srt(:) !Contains dimension for cprj array ordered by atom type
134 : : integer,allocatable :: indsym(:,:,:) !atom indexing for symmetries
135 : : integer,allocatable :: irrzon(:,:,:) !irreducible Brillouin zone
136 : : integer,allocatable :: kg(:,:) !red. coord. of G vecs
137 : : integer,allocatable :: nattyp(:) !nb of atoms of different types
138 : : integer,allocatable :: npwarr(:) !number of PW at each k-point
139 : : integer,allocatable :: symrec(:,:,:) !sym. operations in recip space
140 : : real(dp) :: gprimd(3,3) !primitive cell vectors in recip space
141 : : real(dp) :: gmet(3,3) !metric tensor in recip space
142 : : real(dp) :: rprimd(3,3) !prim cell vectors in direct space
143 : : real(dp) :: rmet(3,3) !metric tensor in direct space
144 : : real(dp),allocatable :: cg(:,:) !WF coefficients in PW basis <k+G|psi_nk>
145 : : real(dp),allocatable :: cg0(:,:) !Initial WF coefficients in PW basis <k+G|psi_nk>
146 : : real(dp),allocatable :: current(:,:)!Current density
147 : : real(dp),allocatable :: eigen(:) !eigen-energies
148 : : real(dp),allocatable :: eigen0(:) !Initial eigen-energies (at t=0)
149 : : real(dp),allocatable :: grvdw(:,:) !Gradient of the total energy coming
150 : : !from VDW dispersion correction !FB: Needed?
151 : : real(dp),allocatable :: nhat(:,:) !compensation charge density
152 : : real(dp),allocatable :: nhatgr(:,:,:) !gradient of nhat
153 : : real(dp),allocatable :: occ(:) !occupation numbers
154 : : real(dp),allocatable :: occ0(:) !Initial occupation numbers
155 : : real(dp),allocatable :: phnons(:,:,:) !For symmetries (nonsymmorphic translation phases)
156 : : real(dp),allocatable :: ph1d(:,:) !Structure factor phase: exp(2Pi i G.xred)
157 : : !on coarse grid
158 : : real(dp),allocatable :: ph1df(:,:) !Structure factor phase: exp(2Pi i G.xred) for G
159 : : !on fine grid
160 : : real(dp),allocatable :: rhog(:,:) !charge density in recip space
161 : : real(dp),allocatable :: rhor(:,:) !charge density in direct space
162 : : real(dp),allocatable :: taug(:,:) !kin ener density in recip space !FB: Needed?
163 : : real(dp),allocatable :: taur(:,:) !kin ener density in direct space !FB: Needed?
164 : : real(dp),allocatable :: vhartr(:) !Hartree part of the potential
165 : : real(dp),allocatable :: vpsp(:) !PSP part of the potential
166 : : real(dp),allocatable :: vtrial(:,:) !"Trial" potential
167 : : real(dp),allocatable :: vxc(:,:) !XC part of the potential
168 : : real(dp),allocatable :: vxc_hybcomp(:,:) !Hybrid part of the xc potential !FB:Needed?
169 : : real(dp),allocatable :: vxctau(:,:,:) !dV_{XC}/dtau (tau = kin. ener density)
170 : : !for mGGAs !FB: Needed?
171 : : real(dp),allocatable :: xred(:,:,:) !red. coord. of atoms
172 : : real(dp),allocatable :: xccc3d(:) !3D core electron density
173 : : !for XC core correction
174 : : real(dp),allocatable :: xcctau3d(:) !3D core electron kin ener density
175 : : !for XC core correction
176 : : real(dp),allocatable :: ylm(:,:) !real spherical harmonics for each k+G
177 : : real(dp),allocatable :: ylmgr(:,:,:)!real spherical harmonics gradients !FB: Needed?
178 : : type(pawcprj_type), allocatable :: cprj(:,:) !projectors applied on WF <p_lmn|C_nk>
179 : : type(pawcprj_type), allocatable :: cprj0(:,:) !projectors applied on WF <p_lmn|C_nk>
180 : : type(paw_an_type), allocatable :: paw_an(:) !various arrays on angular mesh
181 : : type(pawfgrtab_type),allocatable :: pawfgrtab(:) !PAW atomic data on fine grid
182 : : type(paw_ij_type), allocatable :: paw_ij(:) !various arrays on partial waves (i,j channels)
183 : : type(pawrad_type), pointer :: pawrad(:) => NULL() !radial grid in PAW sphere
184 : : type(pawrhoij_type), pointer :: pawrhoij(:) => NULL() !operator rho_ij= <psi|p_i><p_j|psi>
185 : : type(pawtab_type), pointer :: pawtab(:) => NULL() !tabulated PAW atomic data
186 : :
187 : : contains
188 : :
189 : : procedure :: init => tdks_init
190 : : procedure :: free => tdks_free
191 : :
192 : : end type tdks_type
193 : : !!***
194 : :
195 : : contains
196 : : !!***
197 : :
198 : : !!****f* m_rttddft_tdks/tdks_init
199 : : !!
200 : : !! NAME
201 : : !! tdks_init
202 : : !!
203 : : !! FUNCTION
204 : : !! Initialize the tdks object
205 : : !!
206 : : !! INPUTS
207 : : !! codvsn = code version
208 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
209 : : !! dtset <type(dataset_type)> = all input variables for this dataset
210 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
211 : : !! pawang <type(pawang_type)> = paw angular mesh and related data
212 : : !! pawrad(ntypat*usepaw) <type(pawrad_type)> = paw radial mesh and related data
213 : : !! pawtab(ntypat*usepaw) <type(pawtab_type)> = paw tabulated starting data
214 : : !! psps <type(pseudopotential_type)> = variables related to pseudopotentials
215 : : !!
216 : : !! OUTPUT
217 : : !! tdks <class(tdks_type)> = the tdks object to initialize
218 : : !!
219 : : !! SOURCE
220 : 50 : subroutine tdks_init(tdks ,codvsn, dtfil, dtset, mpi_enreg, pawang, pawrad, pawtab, psps)
221 : :
222 : : !Arguments ------------------------------------
223 : : !scalars
224 : : class(tdks_type), intent(inout) :: tdks
225 : : character(len=8), intent(in) :: codvsn
226 : : type(datafiles_type), intent(in) :: dtfil
227 : : type(dataset_type), intent(inout) :: dtset
228 : : type(MPI_type), intent(inout) :: mpi_enreg
229 : : type(pawang_type), intent(inout),target :: pawang
230 : : type(pseudopotential_type), intent(inout) :: psps
231 : : !arrays
232 : : type(pawrad_type), intent(inout),target :: pawrad(psps%ntypat*psps%usepaw)
233 : : type(pawtab_type), intent(inout),target :: pawtab(psps%ntypat*psps%usepaw)
234 : :
235 : : !Local variables-------------------------------
236 : : !scalars
237 : : integer :: ierr
238 : : integer :: my_natom
239 : : integer :: ncpgr
240 : : integer :: psp_gencond
241 : : real(dp) :: ecut_eff
242 : : character(len=500) :: msg
243 : : character(len=fnlen) :: fname_wfk
244 : : type(extfpmd_type),pointer :: extfpmd => null()
245 : : !arrays
246 : 50 : real(dp),allocatable :: doccde(:)
247 : :
248 : : ! ***********************************************************************
249 : :
250 : 50 : my_natom=mpi_enreg%my_natom
251 : :
252 : : !1) Various initializations & checks (MPI, PW, FFT, PSP, Symmetry ...)
253 : 50 : call first_setup(codvsn,dtfil,dtset,ecut_eff,mpi_enreg,pawrad,pawtab,psps,psp_gencond,tdks)
254 : :
255 : : !2) Deals with restart and setup some basic variables and filenames
256 : 50 : tdks%first_step = 1
257 : 50 : tdks%fname_tdener = dtfil%fnameabo_td_ener
258 : 50 : tdks%fname_wfk0 = dtfil%fnamewffk
259 : 50 : fname_wfk = dtfil%fnamewffk
260 : 50 : tdks%fname_tdef = dtfil%fnameabo_td_ef
261 : 50 : tdks%fname_current = dtfil%fnameabo_td_current
262 [ + + ]: 50 : if (dtset%td_restart > 0) then
263 [ + + ]: 20 : if (mpi_enreg%me == 0) then
264 [ - + ]: 6 : if (open_file('TD_RESTART', msg, newunit=tdks%tdrestart_unit, status='old', form='formatted') /= 0) then
265 : 0 : write(msg,'(a,a,a)') 'Error while trying to open file TD_RESTART needed to restart the calculation.'
266 : 0 : ABI_ERROR(msg)
267 : : end if
268 : 6 : read(tdks%tdrestart_unit,*) tdks%first_step
269 : 6 : tdks%first_step = tdks%first_step + 1
270 : 6 : read(tdks%tdrestart_unit,*) tdks%fname_wfk0
271 : 6 : read(tdks%tdrestart_unit,*) fname_wfk
272 : 6 : read(tdks%tdrestart_unit,*) tdks%fname_tdener
273 [ + + ]: 6 : if (dtset%td_ef_type /= 0) then
274 : 3 : read(tdks%tdrestart_unit,*) tdks%fname_tdef
275 : : end if
276 [ + + ]: 6 : if (dtset%prtcurrent /= 0) then
277 : 3 : read(tdks%tdrestart_unit,*) tdks%fname_current
278 : : end if
279 : : end if
280 : : !Send to all procs
281 : 20 : call xmpi_bcast(tdks%first_step,0,mpi_enreg%comm_world,ierr)
282 : 20 : call xmpi_bcast(tdks%fname_tdener,0,mpi_enreg%comm_world,ierr)
283 : 20 : call xmpi_bcast(tdks%fname_wfk0,0,mpi_enreg%comm_world,ierr)
284 : 20 : call xmpi_bcast(fname_wfk,0,mpi_enreg%comm_world,ierr)
285 : 20 : call xmpi_bcast(tdks%fname_tdef,0,mpi_enreg%comm_world,ierr)
286 : 20 : call xmpi_bcast(tdks%fname_current,0,mpi_enreg%comm_world,ierr)
287 : : else
288 [ + + ]: 30 : if (mpi_enreg%me == 0) then
289 [ - + ]: 16 : if (open_file('TD_RESTART', msg, newunit=tdks%tdrestart_unit, status='replace', form='formatted') /= 0) then
290 : 0 : write(msg,'(a,a,a)') 'Error while trying to open file TD_RESTART.'
291 : 0 : ABI_ERROR(msg)
292 : : end if
293 : : end if
294 : : end if
295 : :
296 : : !3) Reads initial KS orbitals from file (calls inwffil)
297 : 50 : call read_wfk(dtfil,dtset,ecut_eff,fname_wfk,mpi_enreg,tdks)
298 : :
299 : : !4) Init occupation numbers
300 [ - + + - : 150 : ABI_MALLOC(tdks%occ0,(dtset%mband*dtset%nkpt*dtset%nsppol))
+ - + - -
+ - + -
+ ]
301 [ + + ]: 9138 : tdks%occ0(:)=dtset%occ_orig(:,1)
302 : : !calc occupation number with metallic occupation using the previously read WF
303 [ + + ]: 50 : if (dtset%occopt>=3.and.dtset%occopt<=9) then ! allowing for occopt 9
304 [ - + + - : 48 : ABI_MALLOC(doccde,(dtset%mband*dtset%nkpt*dtset%nsppol))
- + - + ]
305 : : call newocc(doccde,tdks%eigen0,tdks%energies%entropy_ks,tdks%energies%e_fermie, &
306 : : & tdks%energies%e_fermih,dtset%ivalence,dtset%spinmagntarget, &
307 : : & dtset%mband,dtset%nband,dtset%nelect,dtset%ne_qFD,dtset%nh_qFD, &
308 : : & dtset%nkpt,dtset%nspinor,dtset%nsppol,tdks%occ0,dtset%occopt, &
309 : 24 : & dtset%prtvol,dtset%tphysel,dtset%tsmear,dtset%wtk,extfpmd=extfpmd)
310 [ - + ]: 24 : ABI_FREE(doccde)
311 : : end if
312 : :
313 : : !5) Some further initialization (Mainly for PAW and allocation of arrays for Hamiltonian and densities)
314 : 50 : call second_setup(dtset,mpi_enreg,pawang,pawrad,pawtab,psps,psp_gencond,tdks)
315 : :
316 : : !6) TD external elec. field perturbation
317 [ + + - + ]: 50 : if (dtset%td_ef_type/=0 .and. psps%useylm/=1) ABI_ERROR("TD Electric field only works with spherical harmonics (useylm=1)")
318 [ + + - + ]: 50 : if (dtset%td_ef_type/=0 .and. psps%usepaw/=1) ABI_ERROR("TD Electric field only works with PAW")
319 : : !Init vector potential and associated constants
320 [ + + - + ]: 50 : if (dtset%td_ef_type/=0 .or. dtset%prtcurrent/=0) then
321 [ - + + - : 63 : ABI_MALLOC(tdks%current,(3,dtset%nsppol))
+ - + - -
+ - + -
+ ]
322 [ + + + + ]: 134 : tdks%current = zero
323 : : end if
324 : : call tdks%tdef%init(dtset%td_ef_type,dtset%td_ef_pol,dtset%td_ef_ezero,dtset%td_ef_tzero, &
325 : 50 : & dtset%td_ef_lambda,dtset%td_ef_tau,dtset%td_ef_induced_vecpot,dtset%nkpt,dtset%kptns)
326 [ + + ]: 50 : if (dtset%td_restart /=0) then
327 : 20 : call tdks%tdef%restart(mpi_enreg, tdks%tdrestart_unit)
328 : : call tdks%tdef%update(dtset,mpi_enreg,(tdks%first_step-1)*dtset%dtele,tdks%rprimd,tdks%gprimd,tdks%kg, &
329 : 20 : & psps%mpsang,tdks%npwarr,tdks%ylm,tdks%ylmgr,tdks%current,update_vecpot_ind=.false.)
330 : : else
331 : : call tdks%tdef%update(dtset,mpi_enreg,(tdks%first_step-1)*dtset%dtele,tdks%rprimd,tdks%gprimd,tdks%kg, &
332 : 30 : & psps%mpsang,tdks%npwarr,tdks%ylm,tdks%ylmgr,tdks%current)
333 : : end if
334 : :
335 : : !7) Keep initial cg and cproj in memory for occupations
336 : : !Keep initial wavefunction in memory
337 [ + + ]: 50 : if (dtset%td_restart == 0) then
338 [ - + + - : 90 : ABI_MALLOC(tdks%cg0,(2,tdks%mcg))
+ - + - -
+ - + -
+ ]
339 [ + + + + ]: 3468894 : tdks%cg0(:,:) = tdks%cg(:,:)
340 : : end if
341 : : !and associated cprojs to compute occupations
342 [ + + ]: 50 : if (psps%usepaw ==1) then
343 : 30 : ncpgr=0
344 [ - + + - : 12216 : ABI_MALLOC(tdks%cprj0,(dtset%natom,tdks%mcprj))
- + + - +
- + - - +
- + - + +
+ + + ]
345 : 30 : call pawcprj_alloc(tdks%cprj0,ncpgr,tdks%dimcprj)
346 : : call ctocprj(tdks%atindx,tdks%cg0,1,tdks%cprj0,tdks%gmet,tdks%gprimd,0,0,0, &
347 : : & dtset%istwfk,tdks%kg,tdks%tdef%kpa,tdks%mcg,tdks%mcprj,dtset%mgfft, &
348 : : & dtset%mkmem,mpi_enreg,psps%mpsang,dtset%mpw,dtset%natom,tdks%nattyp, &
349 : : & dtset%nband,dtset%natom,dtset%ngfft,dtset%nkpt,dtset%nloalg, &
350 : : & tdks%npwarr,dtset%nspinor,dtset%nsppol,dtset%nsppol,psps%ntypat, &
351 : : & dtset%paral_kgb,tdks%ph1d,psps,tdks%rmet,dtset%typat,tdks%ucvol, &
352 : 30 : & tdks%unpaw,tdks%xred,tdks%ylm,tdks%ylmgr)
353 : : end if
354 [ - + + - : 150 : ABI_MALLOC(tdks%occ,(dtset%mband*dtset%nkpt*dtset%nsppol))
+ - + - -
+ - + -
+ ]
355 : :
356 : : !Keep some additional stuff in memory within the tdks object
357 : 50 : tdks%unpaw = dtfil%unpaw
358 : 50 : tdks%dt = dtset%dtele
359 : 50 : tdks%ntime = dtset%ntime
360 : :
361 : 50 : tdks%pawang => pawang
362 : 50 : tdks%pawrad => pawrad
363 : 50 : tdks%pawtab => pawtab
364 : :
365 : 50 : end subroutine tdks_init
366 : : !!***
367 : :
368 : : !!****f* m_rttddft_tdks/tdks_free
369 : : !!
370 : : !! NAME
371 : : !! tdks_free
372 : : !!
373 : : !! FUNCTION
374 : : !! Free all the memory associated with the tdks object
375 : : !!
376 : : !! INPUTS
377 : : !! tdks <class(tdks_type)> = the tdks object to free
378 : : !! dtset <type(dataset_type)> = all input variables for this dataset
379 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
380 : : !! psps <type(pseudopotential_type)> = variables related to pseudopotentials
381 : : !!
382 : : !! OUTPUT
383 : : !!
384 : : !! SOURCE
385 : 49 : subroutine tdks_free(tdks,dtset,mpi_enreg,psps)
386 : :
387 : : !Arguments ------------------------------------
388 : : !scalars
389 : : class(tdks_type), intent(inout) :: tdks
390 : : type(dataset_type), intent(inout) :: dtset
391 : : type(MPI_type), intent(inout) :: mpi_enreg
392 : : type(pseudopotential_type), intent(inout) :: psps
393 : :
394 : : ! ***********************************************************************
395 : :
396 : : !Destroy hidden save variables
397 : 49 : call bandfft_kpt_destroy_array(bandfft_kpt,mpi_enreg)
398 [ + + ]: 49 : if (psps%usepaw ==1) then
399 : 29 : call destroy_invovl(dtset%nkpt,dtset%gpu_option)
400 : : end if
401 [ - + - - ]: 49 : if(tdks%gemm_nonlop_use_gemm .and. dtset%gpu_option==ABI_GPU_DISABLED) then
402 : 0 : call destroy_gemm_nonlop(dtset%gpu_option)
403 : : end if
404 : :
405 : : !Call type destructors
406 : 49 : call destroy_sc_dmft(tdks%paw_dmft)
407 : 49 : call pawfgr_destroy(tdks%pawfgr)
408 : 49 : call tdks%hdr%free()
409 : :
410 : : !Nullify pointers
411 [ + - ]: 49 : if(associated(tdks%pawang)) tdks%pawang => null()
412 [ + - ]: 49 : if(associated(tdks%pawrad)) tdks%pawrad => null()
413 [ + - ]: 49 : if(associated(tdks%pawtab)) tdks%pawtab => null()
414 : :
415 : : !Deallocate allocatables
416 [ + - ]: 49 : ABI_SFREE(tdks%atindx)
417 [ + - ]: 49 : ABI_SFREE(tdks%atindx1)
418 [ + - ]: 49 : ABI_SFREE(tdks%cg)
419 [ + - ]: 49 : ABI_SFREE(tdks%cg0)
420 [ + + ]: 49 : ABI_SFREE(tdks%current)
421 [ + + ]: 49 : ABI_SFREE(tdks%dimcprj)
422 [ + + ]: 49 : ABI_SFREE(tdks%dimcprj_srt)
423 [ + - ]: 49 : ABI_SFREE(tdks%eigen)
424 [ + - ]: 49 : ABI_SFREE(tdks%eigen0)
425 [ + - ]: 49 : ABI_SFREE(tdks%grvdw)
426 [ + - ]: 49 : ABI_SFREE(tdks%indsym)
427 [ + - ]: 49 : ABI_SFREE(tdks%irrzon)
428 [ + - ]: 49 : ABI_SFREE(tdks%kg)
429 [ + - ]: 49 : ABI_SFREE(tdks%nattyp)
430 [ + - ]: 49 : ABI_SFREE(tdks%nhat)
431 [ + - ]: 49 : ABI_SFREE(tdks%nhatgr)
432 [ + - ]: 49 : ABI_SFREE(tdks%npwarr)
433 [ + - ]: 49 : ABI_SFREE(tdks%occ)
434 [ + - ]: 49 : ABI_SFREE(tdks%occ0)
435 [ + - ]: 49 : ABI_SFREE(tdks%ph1d)
436 [ + - ]: 49 : ABI_SFREE(tdks%ph1df)
437 [ + - ]: 49 : ABI_SFREE(tdks%phnons)
438 [ + - ]: 49 : ABI_SFREE(tdks%rhog)
439 [ + - ]: 49 : ABI_SFREE(tdks%rhor)
440 [ + - ]: 49 : ABI_SFREE(tdks%symrec)
441 [ + - ]: 49 : ABI_SFREE(tdks%taug)
442 [ + - ]: 49 : ABI_SFREE(tdks%taur)
443 [ + - ]: 49 : ABI_SFREE(tdks%vhartr)
444 [ + - ]: 49 : ABI_SFREE(tdks%vpsp)
445 [ + - ]: 49 : ABI_SFREE(tdks%vtrial)
446 [ + - ]: 49 : ABI_SFREE(tdks%vxc)
447 [ + - ]: 49 : ABI_SFREE(tdks%vxctau)
448 [ + - ]: 49 : ABI_SFREE(tdks%vxc_hybcomp)
449 [ + - ]: 49 : ABI_SFREE(tdks%xred)
450 [ + - ]: 49 : ABI_SFREE(tdks%xccc3d)
451 [ + + ]: 49 : ABI_SFREE(tdks%xcctau3d)
452 [ + - ]: 49 : ABI_SFREE(tdks%ylm)
453 [ + - ]: 49 : ABI_SFREE(tdks%ylmgr)
454 [ + - ]: 49 : ABI_SFREE(tdks%tdef%kpa)
455 : :
456 [ + + ]: 49 : if(allocated(tdks%cprj)) then
457 : 29 : call pawcprj_free(tdks%cprj)
458 [ + - + + : 8093 : ABI_FREE(tdks%cprj)
- + - + -
+ ]
459 : : end if
460 [ + + ]: 49 : if(allocated(tdks%cprj0)) then
461 : 29 : call pawcprj_free(tdks%cprj0)
462 [ + - + + : 8093 : ABI_FREE(tdks%cprj0)
- + - + -
+ ]
463 : : end if
464 [ + + ]: 49 : if(allocated(tdks%paw_an)) then
465 : 29 : call paw_an_free(tdks%paw_an)
466 [ + - + + : 95 : ABI_FREE(tdks%paw_an)
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + -
+ ]
467 : : end if
468 [ + + ]: 49 : if(allocated(tdks%pawfgrtab)) then
469 : 29 : call pawfgrtab_free(tdks%pawfgrtab)
470 [ + - + + : 95 : ABI_FREE(tdks%pawfgrtab)
- + - + -
+ - + - +
- + - + -
+ - + ]
471 : : end if
472 [ + + ]: 49 : if(allocated(tdks%paw_ij)) then
473 : 29 : call paw_ij_free(tdks%paw_ij)
474 [ + - + + : 95 : ABI_FREE(tdks%paw_ij)
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + ]
475 : : end if
476 [ + - ]: 49 : if(associated(tdks%pawrhoij)) then
477 : 49 : call pawrhoij_free(tdks%pawrhoij)
478 [ - + ]: 49 : ABI_FREE(tdks%pawrhoij)
479 : : end if
480 : :
481 : 49 : end subroutine tdks_free
482 : : !!***
483 : :
484 : : !!****f* m_rttddft_tdks/first_setup
485 : : !!
486 : : !! NAME
487 : : !! first_setup
488 : : !!
489 : : !! FUNCTION
490 : : !! Intialize many important quantities before running RT-TDDFT
491 : : !! (PW, FFT, PSP, Symmetry etc.)
492 : : !!
493 : : !! INPUTS
494 : : !! codvsn = code version
495 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
496 : : !! dtset <type(dataset_type)> = all input variables for this dataset
497 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
498 : : !! pawrad(ntypat*usepaw) <type(pawrad_type)> = paw radial mesh and related data
499 : : !! pawtab(ntypat*usepaw) <type(pawtab_type)> = paw tabulated starting data
500 : : !! psps <type(pseudopotential_type)> = variables related to pseudopotentials
501 : : !! tdks <type(tdks_type)> = the tdks object to initialize
502 : : !!
503 : : !! OUTPUT
504 : : !! psp_gencond <integer> = store conditions for generating psp
505 : : !! ecut_eff <real(dp)> = effective PW cutoff energy
506 : : !!
507 : : !! NOTES
508 : : !! USE OF FFT GRIDS:
509 : : !! =================
510 : : !! In case of PAW:
511 : : !! ---------------
512 : : !! Two FFT grids are used:
513 : : !! - A "coarse" FFT grid (defined by ecut)
514 : : !! for the application of the Hamiltonian on the plane waves basis.
515 : : !! It is defined by nfft, ngfft, mgfft, ...
516 : : !! Hamiltonian, wave-functions, density related to WFs (rhor here), ...
517 : : !! are expressed on this grid.
518 : : !! - A "fine" FFT grid (defined) by ecutdg)
519 : : !! for the computation of the density inside PAW spheres.
520 : : !! It is defined by nfftf, ngfftf, mgfftf, ...
521 : : !! Total density, potentials, ...
522 : : !! are expressed on this grid.
523 : : !! In case of norm-conserving:
524 : : !! ---------------------------
525 : : !! - Only the usual FFT grid (defined by ecut) is used.
526 : : !! It is defined by nfft, ngfft, mgfft, ...
527 : : !! For compatibility reasons, (nfftf,ngfftf,mgfftf)
528 : : !! are set equal to (nfft,ngfft,mgfft) in that case.
529 : : !! In case of wavelets:
530 : : !! --------------------
531 : : !! - Only the usual FFT grid (defined by wvl_crmult) is used.
532 : : !! It is defined by nfft, ngfft, mgfft, ... This is strictly not
533 : : !! an FFT grid since its dimensions are not suited for FFTs. They are
534 : : !! defined by wvl_setngfft().
535 : : !! For compatibility reasons, (nfftf,ngfftf,mgfftf)
536 : : !! are set equal to (nfft,ngfft,mgfft) in that case.
537 : : !!
538 : : !! SOURCE
539 : 50 : subroutine first_setup(codvsn,dtfil,dtset,ecut_eff,mpi_enreg,pawrad,pawtab,psps,psp_gencond,tdks)
540 : :
541 : : !Arguments ------------------------------------
542 : : !scalars
543 : : character(len=8), intent(in) :: codvsn
544 : : integer, intent(out) :: psp_gencond
545 : : real(dp), intent(out) :: ecut_eff
546 : : type(datafiles_type), intent(in) :: dtfil
547 : : type(dataset_type), intent(inout) :: dtset
548 : : type(pseudopotential_type), intent(inout) :: psps
549 : : type(MPI_type), intent(inout) :: mpi_enreg
550 : : type(tdks_type), intent(inout) :: tdks
551 : : !arrays
552 : : type(pawrad_type), intent(inout) :: pawrad(psps%ntypat*psps%usepaw)
553 : : type(pawtab_type), intent(inout) :: pawtab(psps%ntypat*psps%usepaw)
554 : :
555 : : !Local variables-------------------------------
556 : : !scalars
557 : : integer,parameter :: response=0, cplex=1
558 : : integer :: comm_psp
559 : : integer :: gscase
560 : : integer :: iatom, ierr, itypat, indx
561 : : integer :: mgfftf, my_natom
562 : : integer :: npwmin, nfftot
563 : : real(dp) :: gsqcut_eff, gsqcutc_eff
564 : : real(dp) :: ecutdg_eff
565 : 50 : type(ebands_t) :: bstruct
566 : : !arrays
567 : : character(len=500) :: msg
568 : 50 : integer, allocatable :: npwarr_(:)
569 : : integer :: ngfft(18)
570 : : integer :: ngfftf(18)
571 : 100 : integer :: npwtot(dtset%nkpt)
572 : :
573 : : ! ***********************************************************************
574 : :
575 : 50 : my_natom=mpi_enreg%my_natom
576 : :
577 : : !** Init FFT grid(s) sizes (be careful !)
578 : : !See NOTES in the comments at the beginning of this subroutine.
579 : 50 : tdks%nfft = dtset%nfft
580 : : call pawfgr_init(tdks%pawfgr,dtset,mgfftf,tdks%nfftf,ecut_eff,ecutdg_eff, &
581 : 50 : & ngfft,ngfftf)
582 : :
583 : : !** Init to zero different energies
584 : 50 : call tdks%energies%init()
585 : 50 : tdks%ecore = zero
586 : 50 : tdks%etot = zero
587 : :
588 : : !** various additional setup mostly related to fft grids and the box (rprimd, metric..)
589 : : call setup1(dtset%acell_orig,tdks%bantot,dtset,ecutdg_eff,ecut_eff,tdks%gmet, &
590 : : & tdks%gprimd,gsqcut_eff,gsqcutc_eff,ngfftf,ngfft,dtset%nkpt, &
591 : : & dtset%nsppol,response,tdks%rmet,dtset%rprim_orig,tdks%rprimd, &
592 : 50 : & tdks%ucvol,psps%usepaw)
593 : :
594 : : !** Set up the basis sphere of planewaves
595 [ - + + - : 150 : ABI_MALLOC(tdks%npwarr,(dtset%nkpt))
+ - + - -
+ - + -
+ ]
596 [ - + + - : 150 : ABI_MALLOC(tdks%kg,(3,dtset%mpw*dtset%mkmem))
+ - + - -
+ - + -
+ ]
597 : : call kpgio(ecut_eff,dtset%exchn2n3d,tdks%gmet,dtset%istwfk,tdks%kg,dtset%kptns, &
598 : : & dtset%mkmem,dtset%nband,dtset%nkpt,'PERS',mpi_enreg,dtset%mpw, &
599 : 50 : & tdks%npwarr,npwtot,dtset%nsppol)
600 : : call bandfft_kpt_init1(bandfft_kpt,dtset%istwfk,tdks%kg,dtset%mgfft,dtset%mkmem, &
601 : : & mpi_enreg,dtset%mpw,dtset%nband,dtset%nkpt,tdks%npwarr, &
602 : 50 : & dtset%nsppol)
603 : :
604 : : !** Use efficient BLAS calls for computing the non local potential (No GPU yet)
605 [ - + - - ]: 50 : if(dtset%use_gemm_nonlop == 1 .and. dtset%gpu_option==ABI_GPU_DISABLED) then
606 : : ! set global variable
607 : 0 : tdks%gemm_nonlop_use_gemm = .true.
608 : 0 : call init_gemm_nonlop(dtset%gpu_option)
609 : : else
610 : 50 : tdks%gemm_nonlop_use_gemm = .false.
611 : : end if
612 : :
613 : : !** TODO: uncomment when gemm_nonlop can be used on GPU
614 : : ! if(dtset%use_gemm_nonlop == 1 .and. dtset%gpu_option/=ABI_GPU_DISABLED) then
615 : : ! ! set global variable
616 : : ! tdks%gemm_nonlop_use_gemm_gpu = .true.
617 : : ! !call init_gemm_nonlop_gpu(dtset%nkpt)
618 : : ! else
619 : : ! tdks%gemm_nonlop_use_gemm_gpu = .false.
620 : : ! end if
621 : :
622 : : !** Initialize band structure datatype
623 [ - + + - : 150 : ABI_MALLOC(npwarr_,(dtset%nkpt))
+ - + - -
+ - + ]
624 [ + + ]: 958 : npwarr_(:)=tdks%npwarr(:)
625 [ + + ]: 50 : if (dtset%paral_kgb/=0) then
626 : 32 : call xmpi_sum(npwarr_,mpi_enreg%comm_bandfft,ierr)
627 : : end if
628 : 50 : call bstruct%from_dtset(dtset, npwarr_)
629 [ - + ]: 50 : ABI_FREE(npwarr_)
630 : 50 : call unpack_eneocc(dtset%nkpt,dtset%nsppol,bstruct%mband,bstruct%nband,dtset%occ_orig(:,1),bstruct%occ,val=zero)
631 : :
632 : : !** Open and read pseudopotential files
633 : 50 : comm_psp=mpi_enreg%comm_cell
634 : : call pspini(dtset,dtfil,tdks%ecore,psp_gencond,gsqcutc_eff,gsqcut_eff,pawrad, &
635 : 50 : & pawtab,psps,tdks%rprimd,comm_mpi=comm_psp)
636 : :
637 : : !In case of isolated computations, ecore must be set to zero
638 : : !because its contribution is counted in the ewald energy as the ion-ion interaction.
639 [ - + ]: 50 : if (dtset%icoulomb == 1) tdks%ecore = zero
640 : :
641 : : !Include core energy?
642 : 100 : select case(dtset%usepotzero)
643 : : case(0,1)
644 : 50 : tdks%energies%e_corepsp = tdks%ecore / tdks%ucvol
645 : 50 : tdks%energies%e_corepspdc = zero
646 : : case(2)
647 : : ! No need to include the PspCore energy since it is already included in the
648 : : ! local pseudopotential (vpsp)
649 : 0 : tdks%energies%e_corepsp = zero
650 [ + - - ]: 50 : tdks%energies%e_corepspdc = zero
651 : : end select
652 : :
653 : : !** Initialize PAW atomic occupancies
654 [ + + + - : 218 : ABI_MALLOC(tdks%pawrhoij,(my_natom*psps%usepaw))
+ - + + -
+ - + + +
+ - ]
655 [ + + ]: 50 : if (psps%usepaw == 1) then
656 : : call initrhoij(dtset%pawcpxocc,dtset%lexexch,dtset%lpawu,my_natom,dtset%natom, &
657 : : & dtset%nspden,dtset%nspinor,dtset%nsppol,dtset%ntypat, &
658 : : & tdks%pawrhoij,dtset%pawspnorb,pawtab,cplex,dtset%spinat, &
659 : : & dtset%typat,comm_atom=mpi_enreg%comm_atom, &
660 [ + - - - : 30 : & mpi_atmtab=mpi_enreg%my_atmtab)
- - - + -
- - - - -
- - - - -
- - - ]
661 : : end if
662 : :
663 : : !Nullify wvl_data. It is important to do so irregardless of the value of usewvl
664 : : !Only needed here because hdr%init requires a wvl object for the wvl%descr input
665 : 50 : call nullify_wvl_data(tdks%wvl)
666 : :
667 : : !** Initialize header
668 : 50 : gscase=0
669 : : call tdks%hdr%init(bstruct,codvsn,dtset,pawtab,gscase,psps,tdks%wvl%descr,&
670 : 50 : comm_atom=mpi_enreg%comm_atom,mpi_atmtab=mpi_enreg%my_atmtab)
671 : :
672 : : !Clean band structure datatype
673 : 50 : call bstruct%free()
674 : :
675 : : !** PW basis set: test if the problem is ill-defined.
676 [ + + ]: 958 : npwmin=minval(tdks%hdr%npwarr(:))
677 [ - + ]: 50 : if (dtset%mband > npwmin) then
678 : : ! No way we can solve the problem. Abort now!
679 : 0 : write(msg,"(2(a,i0),4a)") "Number of bands nband= ",dtset%mband, &
680 : 0 : & " > number of planewaves npw= ",npwmin,ch10, &
681 : 0 : & "The number of eigenvectors cannot be greater that the size of the Hamiltonian!",&
682 : 0 : & ch10, "Action: decrease nband or, alternatively, increase ecut"
683 [ # # ]: 0 : if (dtset%ionmov/=23) then
684 : 0 : ABI_ERROR(msg)
685 : : else
686 : 0 : ABI_WARNING(msg)
687 : : end if
688 : :
689 [ - + ]: 50 : else if (dtset%mband >= 0.9 * npwmin) then
690 : : ! Warn the user
691 : 0 : write(msg,"(a,i0,a,f6.1,4a)") "Number of bands nband= ",dtset%mband, &
692 : 0 : & " >= 0.9 * maximum number of planewaves= ",0.9*npwmin,ch10,&
693 : 0 : & "This could lead to some instabilities, you might want to decrease nband or increase ecut!", &
694 : 0 : & ch10,"Assume experienced user. Execution will continue."
695 : 0 : ABI_WARNING(msg)
696 : : end if
697 : :
698 : : !** Initialize symmetry
699 : 50 : nfftot=ngfft(1)*ngfft(2)*ngfft(3)
700 [ - + + - : 200 : ABI_MALLOC(tdks%irrzon,(nfftot**(1-1/dtset%nsym),2,(dtset%nspden/dtset%nsppol)-3*(dtset%nspden/4)))
- + + - +
- + - - +
- + - + ]
701 [ - + + - : 200 : ABI_MALLOC(tdks%phnons,(2,nfftot**(1-1/dtset%nsym),(dtset%nspden/dtset%nsppol)-3*(dtset%nspden/4)))
- + + - +
- + - - +
- + - + ]
702 [ - + + - : 200 : ABI_MALLOC(tdks%indsym,(4,dtset%nsym,dtset%natom))
- + + - +
- + - - +
- + - + ]
703 [ - + + - : 150 : ABI_MALLOC(tdks%symrec,(3,3,dtset%nsym))
+ - + - -
+ - + -
+ ]
704 [ + + + + : 1360680 : tdks%irrzon(:,:,:)=0
+ + ]
705 [ + + + + : 2040820 : tdks%phnons(:,:,:)=zero
+ + ]
706 [ + + + + : 104838 : tdks%indsym(:,:,:)=0
+ + ]
707 [ + + + + : 76178 : tdks%symrec(:,:,:)=0
+ + ]
708 : :
709 : : !Do symmetry stuff if nsym>1
710 [ + - ]: 50 : if (dtset%nsym>1) then
711 : : call setsym(tdks%indsym,tdks%irrzon,dtset%iscf,dtset%natom, &
712 : : & nfftot,ngfft,dtset%nspden,dtset%nsppol,dtset%nsym, &
713 : : & tdks%phnons,dtset%symafm,tdks%symrec,dtset%symrel, &
714 : 50 : & dtset%tnons,dtset%typat,dtset%xred_orig)
715 : :
716 : : !Make sure dtset%iatfix does not break symmetry
717 : 50 : call fixsym(dtset%iatfix,tdks%indsym,dtset%natom,dtset%nsym)
718 : : else
719 : : !The symrec array is used by initberry even in case nsym = 1
720 [ # # # # ]: 0 : tdks%symrec(:,:,1) = 0
721 : 0 : tdks%symrec(1,1,1) = 1 ; tdks%symrec(2,2,1) = 1 ; tdks%symrec(3,3,1) = 1
722 : : end if
723 : :
724 : : !** Initialize and eventually symmetrize reduced atomic coordinates
725 [ - + + - : 200 : ABI_MALLOC(tdks%xred,(3,dtset%natom,dtset%nimage))
- + + - +
- + - - +
- + - + ]
726 [ + - + - : 742 : tdks%xred = dtset%xred_orig
+ - - + -
- - - - -
- - + + +
+ + + ]
727 : : !Eventually symmetrize atomic coordinates over space group elements
728 : : call symmetrize_xred(dtset%natom,dtset%nsym,dtset%symrel,dtset%tnons,tdks%xred, &
729 : 50 : & indsym=tdks%indsym)
730 : :
731 : : !** Create the atindx array
732 : : !** index table of atoms, in order for them to be used type after type.
733 [ - + + - : 150 : ABI_MALLOC(tdks%atindx,(dtset%natom))
+ - + - -
+ - + -
+ ]
734 [ - + + - : 100 : ABI_MALLOC(tdks%atindx1,(dtset%natom))
- + - + -
+ ]
735 [ - + + - : 150 : ABI_MALLOC(tdks%nattyp,(psps%ntypat))
+ - + - -
+ - + -
+ ]
736 : 50 : indx=1
737 [ + + ]: 100 : do itypat=1,psps%ntypat
738 : 50 : tdks%nattyp(itypat)=0
739 [ + + ]: 248 : do iatom=1,dtset%natom
740 [ + - ]: 198 : if(dtset%typat(iatom)==itypat)then
741 : 148 : tdks%atindx(iatom)=indx
742 : 148 : tdks%atindx1(indx)=iatom
743 : 148 : indx=indx+1
744 : 148 : tdks%nattyp(itypat)=tdks%nattyp(itypat)+1
745 : : end if
746 : : end do
747 : : end do
748 : :
749 : : !** Calculate zion: the total positive charge acting on the valence electrons
750 : 50 : tdks%zion=zero
751 [ + + ]: 198 : do iatom=1,dtset%natom
752 : 198 : tdks%zion=tdks%zion+psps%ziontypat(dtset%typat(iatom))
753 : : end do
754 : :
755 [ - + - + : 100 : end subroutine first_setup
- + - + -
+ - + - +
- + - + -
+ - + ]
756 : : !!***
757 : :
758 : : !!****f* m_rttddft_tdks/second_setup
759 : : !!
760 : : !! NAME
761 : : !! second_setup
762 : : !!
763 : : !! FUNCTION
764 : : !! Further important initialization required after reading WFK and computing
765 : : !! occupation numbers in paticular related to PAW
766 : : !!
767 : : !! INPUTS
768 : : !! dtset <type(dataset_type)> = all input variables for this dataset
769 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
770 : : !! pawang <type(pawang_type)> = paw angular mesh and related data
771 : : !! pawrad(ntypat*usepaw) <type(pawrad_type)> = paw radial mesh and related data
772 : : !! pawtab(ntypat*usepaw) <type(pawtab_type)> = paw tabulated starting data
773 : : !! psps <type(pseudopotential_type)> = variables related to pseudopotentials
774 : : !! psp_gencond <integer> = store conditions for generating psp
775 : : !! tdks <type(tdks_type)> = the tdks object to initialize
776 : : !!
777 : : !! OUTPUT
778 : : !!
779 : : !! SOURCE
780 : 50 : subroutine second_setup(dtset, mpi_enreg, pawang, pawrad, pawtab, psps, psp_gencond, tdks)
781 : :
782 : : !Arguments ------------------------------------
783 : : !scalars
784 : : integer, intent(in) :: psp_gencond
785 : : type(pawang_type), intent(inout) :: pawang
786 : : type(dataset_type), intent(inout) :: dtset
787 : : type(pseudopotential_type), intent(inout) :: psps
788 : : type(MPI_type), intent(inout) :: mpi_enreg
789 : : !arrays
790 : : type(pawrad_type), intent(inout) :: pawrad(psps%ntypat*psps%usepaw)
791 : : type(pawtab_type), intent(inout) :: pawtab(psps%ntypat*psps%usepaw)
792 : : type(tdks_type), intent(inout) :: tdks
793 : :
794 : : !Local variables-------------------------------
795 : : !scalars
796 : : logical :: call_pawinit
797 : : integer, parameter :: cplex = 1
798 : : integer :: forces_needed
799 : : integer :: gnt_option
800 : : integer :: has_dijhat, has_vhartree, has_dijfock
801 : : integer :: has_dijnd, has_dijU, has_vxctau
802 : : integer :: my_natom, my_nspinor
803 : : integer :: ncpgr
804 : : integer :: optcut, optgr0, optgr1, optgr2, optrad
805 : : integer :: stress_needed
806 : : integer :: ylm_option
807 : : integer :: use_hybcomp, usevxctau
808 : : real(dp) :: boxcut
809 : : real(dp) :: gsqcut_shp
810 : : real(dp) :: hyb_range_fock
811 : : real(dp),parameter :: k0(3)=(/zero,zero,zero/)
812 : : !arrays
813 : 50 : integer,allocatable :: l_size_atm(:)
814 : :
815 : : ! ***********************************************************************
816 : :
817 : 50 : my_natom=mpi_enreg%my_natom
818 : 50 : my_nspinor=max(1,dtset%nspinor/mpi_enreg%nproc_spinor)
819 : :
820 : : !FB: @MT needed?
821 [ + + ]: 50 : if (psps%usepaw==1) then
822 : : call pawrhoij_copy(tdks%hdr%pawrhoij,tdks%pawrhoij,comm_atom=mpi_enreg%comm_atom, &
823 : 30 : & mpi_atmtab=mpi_enreg%my_atmtab)
824 : : end if
825 : :
826 : : !FB: Needed because paw_dmft is required in mkrho
827 : : !PAW related operations
828 : : !Initialize paw_dmft, even if neither dmft not paw are used
829 : 50 : call init_sc_dmft(dtset,psps%mpsang,tdks%paw_dmft,mpi_enreg=mpi_enreg,use_sc_dmft=0)
830 : :
831 : :
832 : : !*** Main PAW initialization ***
833 : :
834 : : !** Setup the Ylm for each k point
835 [ + + ]: 50 : if (psps%useylm==1) then
836 [ - + + - : 120 : ABI_MALLOC(tdks%ylm,(dtset%mpw*dtset%mkmem,psps%mpsang*psps%mpsang*psps%useylm))
- + + - +
- + - - +
- + - + ]
837 [ - + + - : 120 : ABI_MALLOC(tdks%ylmgr,(dtset%mpw*dtset%mkmem,3,psps%mpsang*psps%mpsang*psps%useylm))
- + + - +
- + - - +
- + - + ]
838 : 30 : ylm_option=0
839 : : call initylmg(tdks%gprimd,tdks%kg,dtset%kptns,dtset%mkmem,mpi_enreg,&
840 : : & psps%mpsang,dtset%mpw,dtset%nband,dtset%nkpt,&
841 : 30 : & tdks%npwarr,dtset%nsppol,ylm_option,tdks%rprimd,tdks%ylm,tdks%ylmgr)
842 : : else
843 [ - + - + ]: 20 : ABI_MALLOC(tdks%ylm,(0,0))
844 [ - + - + ]: 20 : ABI_MALLOC(tdks%ylmgr,(0,0,0))
845 : : end if
846 : :
847 : 50 : tdks%mcprj=0;tdks%mband_cprj=0
848 [ + + ]: 50 : if(psps%usepaw==1) then
849 : 30 : gnt_option=1
850 [ + - + - : 30 : if (dtset%pawxcdev==2.or.(dtset%pawxcdev==1.and.dtset%positron/=0)) gnt_option=2
- + ]
851 : :
852 : : !** Test if we have to call pawinit
853 : : ! Some gen-cond have to be added...
854 : 30 : call paw_gencond(dtset,gnt_option,"test",call_pawinit)
855 : :
856 [ + + - + ]: 30 : if (psp_gencond==1.or.call_pawinit) then
857 : 2 : gsqcut_shp=two*abs(dtset%diecut)*dtset%dilatmx**2/pi**2
858 : 2 : hyb_range_fock=zero
859 [ + - ]: 2 : if (dtset%ixc<0) then
860 : 2 : call libxc_functionals_get_hybridparams(hyb_range=hyb_range_fock)
861 : : end if
862 : : call pawinit(dtset%effmass_free,gnt_option,gsqcut_shp,hyb_range_fock, &
863 : : & dtset%pawlcutd,dtset%pawlmix,psps%mpsang,dtset%pawnphi, &
864 : : & dtset%nsym,dtset%pawntheta,pawang,pawrad,dtset%pawspnorb, &
865 : 2 : & pawtab,dtset%pawxcdev,dtset%ixc,dtset%usepotzero)
866 : :
867 : : ! Update internal values
868 : 2 : call paw_gencond(dtset,gnt_option,"save",call_pawinit)
869 : : end if
870 [ + + ]: 60 : psps%n1xccc=maxval(pawtab(1:psps%ntypat)%usetcore)
871 : : call setsym_ylm(tdks%gprimd,pawang%l_max-1,dtset%nsym,dtset%pawprtvol, &
872 : 30 : & tdks%rprimd,tdks%symrec,pawang%zarot)
873 : :
874 : : !** Initialisation of cprj
875 : 30 : tdks%mband_cprj=dtset%mband
876 [ + + ]: 30 : if (dtset%paral_kgb/=0) tdks%mband_cprj=tdks%mband_cprj/mpi_enreg%nproc_band
877 : 30 : tdks%mcprj=my_nspinor*tdks%mband_cprj*dtset%mkmem*dtset%nsppol
878 [ - + + - : 12216 : ABI_MALLOC(tdks%cprj,(dtset%natom,tdks%mcprj))
- + + - +
- + - - +
- + - + +
+ + + +
- ]
879 : 30 : ncpgr=0
880 : : !FB: @MT dimcprj_srt needed?
881 [ - + + - : 90 : ABI_MALLOC(tdks%dimcprj,(dtset%natom))
+ - + - -
+ - + -
+ ]
882 [ - + + - : 60 : ABI_MALLOC(tdks%dimcprj_srt,(dtset%natom))
- + - + -
+ ]
883 : : call pawcprj_getdim(tdks%dimcprj,dtset%natom,tdks%nattyp,dtset%ntypat, &
884 : 30 : & dtset%typat,pawtab,'R')
885 : : call pawcprj_getdim(tdks%dimcprj_srt,dtset%natom,tdks%nattyp,dtset%ntypat, &
886 : 30 : & dtset%typat,pawtab,'O')
887 : : !call pawcprj_alloc(tdks%cprj,ncpgr,dimcprj_srt)
888 : 30 : call pawcprj_alloc(tdks%cprj,ncpgr,tdks%dimcprj)
889 : : !ABI_FREE(dimcprj_srt)
890 : :
891 : : !** Variables/arrays related to the fine FFT grid
892 [ - + + - : 158 : ABI_MALLOC(tdks%pawfgrtab,(my_natom))
+ - + - -
+ - + - +
+ + + - ]
893 [ + - ]: 30 : if (my_natom>0) then
894 : : call pawtab_get_lsize(pawtab,l_size_atm,my_natom,dtset%typat, &
895 : 30 : & mpi_atmtab=mpi_enreg%my_atmtab)
896 : : call pawfgrtab_init(tdks%pawfgrtab,cplex,l_size_atm,dtset%nspden,dtset%typat, &
897 : 30 : & mpi_atmtab=mpi_enreg%my_atmtab,comm_atom=mpi_enreg%comm_atom)
898 [ - + ]: 30 : ABI_FREE(l_size_atm)
899 : : end if
900 [ + + ]: 60 : tdks%usexcnhat=maxval(pawtab(:)%usexcnhat)
901 : :
902 : : !** Variables/arrays related to the PAW spheres
903 [ - + + - : 158 : ABI_MALLOC(tdks%paw_ij,(my_natom))
+ - + - -
+ - + - +
+ + + - ]
904 [ - + + - : 158 : ABI_MALLOC(tdks%paw_an,(my_natom))
+ - + - -
+ - + - +
+ + + - ]
905 : 30 : call paw_an_nullify(tdks%paw_an)
906 : 30 : call paw_ij_nullify(tdks%paw_ij)
907 [ - + ]: 30 : has_dijhat=0; if (dtset%iscf==22) has_dijhat=1
908 [ + - - + ]: 30 : has_vhartree=0; if (dtset%prtvha > 0 .or. dtset%prtvclmb > 0) has_vhartree=1
909 [ + + + + : 302 : has_dijnd=0;if(any(abs(dtset%nucdipmom)>tol8)) has_dijnd=1
+ - - + ]
910 : 30 : has_dijfock=0
911 [ + - ]: 30 : has_dijU=merge(0,1,dtset%usepawu>0) !Be careful on this!
912 : 30 : has_vxctau=pawxc_get_usekden(dtset%ixc)
913 : : call paw_an_init(tdks%paw_an,dtset%natom,dtset%ntypat,0,0,dtset%nspden, &
914 : : & cplex,dtset%pawxcdev,dtset%typat,pawang,pawtab,has_vxc=1, &
915 : : & has_vxctau=has_vxctau,has_vxc_ex=1,has_vhartree=has_vhartree, &
916 : 30 : & comm_atom=mpi_enreg%comm_atom,mpi_atmtab=mpi_enreg%my_atmtab)
917 : : call paw_ij_init(tdks%paw_ij,cplex,dtset%nspinor,dtset%nsppol,dtset%nspden, &
918 : : & dtset%pawspnorb,dtset%natom,dtset%ntypat,dtset%typat,pawtab, &
919 : : & has_dij=1,has_dijfock=has_dijfock,has_dijhartree=1, &
920 : : & has_dijnd=has_dijnd,has_dijso=1,has_dijhat=has_dijhat, &
921 : : & has_dijU=has_dijU,has_pawu_occ=1,has_exexch_pot=1, &
922 : : & nucdipmom=dtset%nucdipmom,comm_atom=mpi_enreg%comm_atom, &
923 : 30 : & mpi_atmtab=mpi_enreg%my_atmtab)
924 : :
925 : : !** Check for non-overlapping spheres
926 : : call chkpawovlp(dtset%natom,psps%ntypat,dtset%pawovlp,pawtab,tdks%rmet, &
927 : 30 : & dtset%typat,tdks%xred)
928 : :
929 : : !** Identify parts of the rectangular grid where the density has to be calculated
930 : 30 : optcut=0;optgr0=dtset%pawstgylm;optgr1=0;optgr2=0;optrad=1-dtset%pawstgylm
931 : 30 : forces_needed=0 !FB TODO Maybe needs to be changed if Ehrenfest?
932 : 30 : stress_needed=0
933 : : if ((forces_needed==1) .or. &
934 [ + + + - : 30 : & (dtset%xclevel==2 .and. dtset%pawnhatxc>0 .and. tdks%usexcnhat>0) .or. &
- + ]
935 : : & (dtset%positron/=0.and.forces_needed==2)) then
936 : 0 : optgr1=dtset%pawstgylm
937 [ # # ]: 0 : if (stress_needed==1) optrad=1; if (dtset%pawprtwf==1) optrad=1
938 : : end if
939 : : call nhatgrid(tdks%atindx1,tdks%gmet,my_natom,dtset%natom, &
940 : : & tdks%nattyp,tdks%pawfgr%ngfft,psps%ntypat,optcut,optgr0,optgr1, &
941 : : & optgr2,optrad,tdks%pawfgrtab,pawtab,tdks%rprimd,dtset%typat, &
942 : : & tdks%ucvol,tdks%xred,comm_atom=mpi_enreg%comm_atom, &
943 : : & mpi_atmtab=mpi_enreg%my_atmtab,comm_fft=mpi_enreg%comm_fft, &
944 : 30 : & distribfft=mpi_enreg%distribfft)
945 : :
946 [ + + ]: 30 : tdks%nhatgrdim=0;if (dtset%xclevel==2) tdks%nhatgrdim=tdks%usexcnhat*dtset%pawnhatxc
947 [ - + ]: 30 : if (tdks%nhatgrdim>0) then
948 [ # # # # : 0 : ABI_MALLOC(tdks%nhatgr,(cplex*tdks%nfftf,dtset%nspden,3*tdks%nhatgrdim))
# # # # #
# # # # #
# # # # #
# # # ]
949 : : else
950 [ - + - + ]: 30 : ABI_MALLOC(tdks%nhatgr,(0,0,0))
951 : : end if
952 : :
953 [ - + + - : 120 : ABI_MALLOC(tdks%nhat,(tdks%nfftf,dtset%nspden*psps%usepaw))
- + + - +
- + - - +
- + - + ]
954 : :
955 : : !Required in the PAW case to compute the inverse of the overlap (invovl) operator
956 : 60 : call init_invovl(dtset%nkpt)
957 : : else
958 [ - + - + ]: 20 : ABI_MALLOC(tdks%nhat,(0,0))
959 [ - + - + ]: 20 : ABI_MALLOC(tdks%nhatgr,(0,0,0))
960 : 20 : tdks%nhatgrdim=0
961 : : end if
962 : :
963 : : !Allocate various required arrays for calculation of the Hamiltonian
964 : : !Potentials
965 [ - + + - : 150 : ABI_MALLOC(tdks%vhartr,(tdks%nfftf))
+ - + - -
+ - + -
+ ]
966 [ + + ]: 989924 : tdks%vhartr=zero
967 [ - + + - : 100 : ABI_MALLOC(tdks%vpsp,(tdks%nfftf))
- + - + -
+ ]
968 [ + + ]: 989924 : tdks%vpsp=zero
969 [ - + + - : 200 : ABI_MALLOC(tdks%vtrial,(tdks%nfftf,dtset%nspden))
- + + - +
- + - - +
- + - + ]
970 [ + + + + ]: 989974 : tdks%vtrial=zero
971 [ - + + - : 150 : ABI_MALLOC(tdks%vxc,(tdks%nfftf,dtset%nspden))
- + + - -
+ - + -
+ ]
972 [ + + + + ]: 989974 : tdks%vxc=zero
973 [ + - ]: 50 : if (psps%n1xccc/=0) then
974 [ - + + - : 100 : ABI_MALLOC(tdks%xccc3d,(tdks%nfftf))
- + - + -
+ ]
975 : : else
976 [ # # # # ]: 0 : ABI_MALLOC(tdks%xccc3d,(0))
977 : : end if
978 [ + + ]: 989924 : tdks%xccc3d=zero
979 [ + + ]: 50 : if (psps%usepaw==1) then
980 [ + - - - : 90 : ABI_MALLOC(tdks%xcctau3d,(tdks%nfftf*dtset%usekden))
+ - - + -
+ - + -
+ ]
981 [ - + ]: 30 : tdks%xcctau3d=zero
982 : : endif
983 : : !For mGGA
984 [ + - ]: 50 : usevxctau=merge(1,0,xc_need_kden(dtset%ixc))
985 [ - + + - : 250 : ABI_MALLOC(tdks%vxctau,(tdks%nfftf,dtset%nspden,4*usevxctau))
- + + - -
+ + - - -
+ - - + -
+ - + -
+ ]
986 [ - + - - : 50 : tdks%vxctau=zero
- - ]
987 : : !For hybrid functionals
988 : 50 : use_hybcomp=0
989 [ - + ]: 50 : if(mod(dtset%fockoptmix,100)==11) use_hybcomp=1
990 [ - + + - : 200 : ABI_MALLOC(tdks%vxc_hybcomp,(tdks%pawfgr%nfft,dtset%nspden*use_hybcomp))
+ - - - +
- - + - +
- + - + ]
991 [ - + - - ]: 50 : tdks%vxc_hybcomp=zero
992 : : !For VDW corrected functionals
993 : 50 : tdks%ngrvdw=0
994 [ - + ]: 50 : if ((dtset%vdw_xc>=5.and.dtset%vdw_xc<=7)) then
995 : 0 : tdks%ngrvdw=dtset%natom
996 : : end if
997 [ + - - - : 150 : ABI_MALLOC(tdks%grvdw,(3,tdks%ngrvdw))
+ - - + -
+ - + -
+ ]
998 [ - + - - ]: 50 : tdks%grvdw=zero
999 : :
1000 : : !Compute large sphere G^2 cut-off (gsqcut) and box / sphere ratio
1001 [ + + ]: 50 : if (psps%usepaw==1) then
1002 : : call getcut(boxcut,dtset%pawecutdg,tdks%gmet,tdks%gsqcut,dtset%iboxcut, &
1003 : 30 : & std_out,k0,tdks%pawfgr%ngfft)
1004 : : else
1005 : : call getcut(boxcut,dtset%ecut,tdks%gmet,tdks%gsqcut,dtset%iboxcut, &
1006 : 20 : & std_out,k0,tdks%pawfgr%ngfft)
1007 : : end if
1008 : :
1009 : : !Compute structure factor phases (exp(2Pi i G.xred)) on coarse and fine grid
1010 [ + - - + : 150 : ABI_MALLOC(tdks%ph1d,(2,3*(2*tdks%pawfgr%mgfftc+1)*dtset%natom))
+ - + - +
- - + - +
- + ]
1011 [ + - - + : 150 : ABI_MALLOC(tdks%ph1df,(2,3*(2*tdks%pawfgr%mgfft+1)*dtset%natom))
+ - + - +
- - + - +
- + ]
1012 : : call getph(tdks%atindx,dtset%natom,tdks%pawfgr%ngfftc(1),tdks%pawfgr%ngfftc(2), &
1013 : 50 : & tdks%pawfgr%ngfftc(3),tdks%ph1d,tdks%xred)
1014 [ + + + - ]: 50 : if (psps%usepaw==1.and.tdks%pawfgr%usefinegrid==1) then
1015 : : call getph(tdks%atindx,dtset%natom,tdks%pawfgr%ngfft(1),tdks%pawfgr%ngfft(2), &
1016 : 30 : & tdks%pawfgr%ngfft(3),tdks%ph1df,tdks%xred)
1017 : : else
1018 [ + + + + ]: 43940 : tdks%ph1df(:,:)=tdks%ph1d(:,:)
1019 : : end if
1020 : :
1021 : : !Allocate memory for density
1022 [ - + + - : 200 : ABI_MALLOC(tdks%rhor,(tdks%nfftf,dtset%nspden))
- + + - +
- + - - +
- + - + ]
1023 [ - + + - : 200 : ABI_MALLOC(tdks%taur,(tdks%nfftf,dtset%nspden*dtset%usekden))
+ - - - +
- - + - +
- + - + ]
1024 [ - + + - : 150 : ABI_MALLOC(tdks%rhog,(2,tdks%nfftf))
+ - + - -
+ - + -
+ ]
1025 [ + - - - : 150 : ABI_MALLOC(tdks%taug,(2,tdks%nfftf*dtset%usekden))
+ - - + -
+ - + -
+ ]
1026 : :
1027 [ - + ]: 50 : end subroutine second_setup
1028 : : !!***
1029 : :
1030 : : !!****f* m_rttddft_tdks/read_wfk
1031 : : !!
1032 : : !! NAME
1033 : : !! read_wfk
1034 : : !!
1035 : : !! FUNCTION
1036 : : !! Reads initial wavefunctions (KS orbitals) in WFK file (call inwffil)
1037 : : !!
1038 : : !! INPUTS
1039 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
1040 : : !! dtset <type(dataset_type)> = all input variables for this dataset
1041 : : !! ecut_eff <real(dp)> = effective PW cutoff energy
1042 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
1043 : : !! tdks <type(tdks_type)> = the tdks object to initialize
1044 : : !!
1045 : : !! OUTPUT
1046 : : !!
1047 : : !! SOURCE
1048 : 50 : subroutine read_wfk(dtfil, dtset, ecut_eff, fname_wfk, mpi_enreg, tdks)
1049 : :
1050 : : !Arguments ------------------------------------
1051 : : !scalars
1052 : : character(len=fnlen), intent(in) :: fname_wfk
1053 : : real(dp), intent(in) :: ecut_eff
1054 : : type(datafiles_type), intent(in) :: dtfil
1055 : : type(dataset_type), intent(inout) :: dtset
1056 : : type(MPI_type), intent(inout) :: mpi_enreg
1057 : : type(tdks_type), intent(inout) :: tdks
1058 : :
1059 : : !Local variables-------------------------------
1060 : : !scalars
1061 : : integer,parameter :: formeig=0
1062 : : integer :: ask_accurate
1063 : : integer :: band
1064 : : integer :: cnt
1065 : : integer :: ierr, ikpt
1066 : : integer :: my_nspinor
1067 : : integer :: optorth
1068 : : integer :: spin
1069 : : type(wffile_type) :: wff1, wffnow
1070 : : !arrays
1071 : : character(len=500) :: msg
1072 : :
1073 : : ! ***********************************************************************
1074 : :
1075 : : !If paral_kgb == 0, it may happen that some processors are idle (no entry in proc_distrb)
1076 : : !but mkmem == nkpt and this can cause integer overflow in mcg or allocation error.
1077 : : !Here we count the number of states treated by the proc. if cnt == 0, mcg is then set to 0.
1078 : 50 : cnt = 0
1079 [ + + ]: 100 : do spin=1,dtset%nsppol
1080 [ + + ]: 1008 : do ikpt=1,dtset%nkpt
1081 [ + + ]: 10046 : do band=1,dtset%nband(ikpt + (spin-1) * dtset%nkpt)
1082 [ + + ]: 9996 : if (.not. proc_distrb_cycle(mpi_enreg%proc_distrb, ikpt, band, band, spin, mpi_enreg%me_kpt)) cnt = cnt + 1
1083 : : end do
1084 : : end do
1085 : : end do
1086 : :
1087 : 50 : my_nspinor=max(1,dtset%nspinor/mpi_enreg%nproc_spinor)
1088 : 50 : tdks%mcg=dtset%mpw*my_nspinor*dtset%mband*dtset%mkmem*dtset%nsppol
1089 [ - + ]: 50 : if (cnt == 0) then
1090 : 0 : tdks%mcg = 0
1091 : 0 : write(msg,"(2(a,i0))")"rank: ",mpi_enreg%me, "does not have wavefunctions to treat. Setting mcg to: ",tdks%mcg
1092 : 0 : ABI_WARNING(msg)
1093 : : end if
1094 : :
1095 [ + - + - : 50 : if (dtset%usewvl == 0 .and. dtset%mpw > 0 .and. cnt /= 0)then
+ - ]
1096 [ - + - + ]: 50 : if (my_nspinor*dtset%mband*dtset%mkmem*dtset%nsppol > floor(real(HUGE(0))/real(dtset%mpw) )) then
1097 : : ierr = 0
1098 : : write (msg,'(9a)')&
1099 : 0 : & "Default integer is not wide enough to store the size of the wavefunction array (mcg).",ch10,&
1100 : 0 : & "This usually happens when paral_kgb == 0 and there are not enough procs to distribute kpts and spins",ch10,&
1101 : 0 : & "Action: if paral_kgb == 0, use nprocs = nkpt * nsppol to reduce the memory per node.",ch10,&
1102 : 0 : & "If tdks does not solve the problem, use paral_kgb 1 with nprocs > nkpt * nsppol and use npfft/npband/npspinor",ch10,&
1103 : 0 : & "to decrease the memory requirements. Consider also OpenMP threads."
1104 : 0 : ABI_ERROR_NOSTOP(msg,ierr)
1105 : : write (msg,'(5(a,i0), 2a)')&
1106 : 0 : & "my_nspinor: ",my_nspinor, ", mpw: ",dtset%mpw, ", mband: ",dtset%mband,&
1107 : 0 : & ", mkmem: ",dtset%mkmem, ", nsppol: ",dtset%nsppol,ch10,&
1108 : 0 : & 'Note: Compiling with large int (int64) requires a full software stack (MPI/FFTW/BLAS...) compiled in int64 mode'
1109 : 0 : ABI_ERROR(msg)
1110 : : end if
1111 : : end if
1112 : :
1113 : : ! Alloc size for wfk and bands
1114 [ - + + - : 150 : ABI_MALLOC_OR_DIE(tdks%cg,(2,tdks%mcg),ierr)
+ - + - +
- + - - +
+ - - + ]
1115 [ - + + - : 150 : ABI_MALLOC(tdks%eigen,(dtset%mband*dtset%nkpt*dtset%nsppol))
+ - + - -
+ - + -
+ ]
1116 [ - + + - : 100 : ABI_MALLOC(tdks%eigen0,(dtset%mband*dtset%nkpt*dtset%nsppol))
- + - + -
+ ]
1117 : :
1118 [ + + ]: 9138 : tdks%eigen(:) = zero
1119 : 50 : ask_accurate=1
1120 : :
1121 : : !Actually read the intial KS orbitals here
1122 [ + + ]: 50 : if (dtset%td_restart /= 1) then
1123 : 30 : write(msg,'(3a)') ch10,'------------------- Reading initial wavefunctions -------------------',ch10
1124 : : else
1125 : 20 : write(msg,'(3a)') ch10,'------------------- Reading wavefunctions for restart ------------------',ch10
1126 : : end if
1127 : 50 : call wrtout(ab_out,msg)
1128 [ + + ]: 50 : if (do_write_log) call wrtout(std_out,msg)
1129 : 50 : wff1%unwff=dtfil%unwff1
1130 : 50 : optorth=0 !No need to orthogonalize the wfk
1131 [ + + + + ]: 650 : tdks%hdr%rprimd=tdks%rprimd
1132 [ + + + + ]: 4724594 : tdks%cg=zero
1133 : : call inwffil(ask_accurate,tdks%cg,dtset,dtset%ecut,ecut_eff,tdks%eigen, &
1134 : : & dtset%exchn2n3d,formeig,tdks%hdr,1,dtset%istwfk,tdks%kg, &
1135 : : & dtset%kptns,dtset%localrdwf,dtset%mband,tdks%mcg,dtset%mkmem, &
1136 : : & mpi_enreg,dtset%mpw,dtset%nband,tdks%pawfgr%ngfft,dtset%nkpt, &
1137 : : & tdks%npwarr,dtset%nsppol,dtset%nsym,dtset%occ_orig,optorth, &
1138 : : & dtset%symafm,dtset%symrel,dtset%tnons,dtfil%unkg,wff1,wffnow, &
1139 : 50 : & dtfil%unwff1,fname_wfk,tdks%wvl)
1140 : :
1141 : : !Close file
1142 : 50 : call WffClose(wff1,ierr)
1143 : :
1144 : : !Keep initial eigenvalues in memory
1145 [ + + ]: 9138 : tdks%eigen0(:) = tdks%eigen(:)
1146 : :
1147 : : !In case of restart also read wfk file containing wave functions at t=0
1148 [ + + + - ]: 50 : if (dtset%td_restart == 1 .and. tdks%fname_wfk0 /= fname_wfk) then
1149 : 20 : write(msg,'(3a)') ch10,'------------------- Reading initial wavefunctions -------------------',ch10
1150 [ - + + - : 60 : ABI_MALLOC_OR_DIE(tdks%cg0,(2,tdks%mcg),ierr)
+ - + - +
- + - - +
+ - - + ]
1151 : 20 : call wrtout(ab_out,msg)
1152 [ + + ]: 20 : if (do_write_log) call wrtout(std_out,msg)
1153 [ + + + + ]: 1255700 : tdks%cg0=zero
1154 : : call inwffil(ask_accurate,tdks%cg0,dtset,dtset%ecut,ecut_eff,tdks%eigen0, &
1155 : : & dtset%exchn2n3d,formeig,tdks%hdr,1,dtset%istwfk,tdks%kg, &
1156 : : & dtset%kptns,dtset%localrdwf,dtset%mband,tdks%mcg,dtset%mkmem, &
1157 : : & mpi_enreg,dtset%mpw,dtset%nband,tdks%pawfgr%ngfft,dtset%nkpt, &
1158 : : & tdks%npwarr,dtset%nsppol,dtset%nsym,dtset%occ_orig,optorth, &
1159 : : & dtset%symafm,dtset%symrel,dtset%tnons,dtfil%unkg,wff1,wffnow, &
1160 : 20 : & dtfil%unwff1,tdks%fname_wfk0,tdks%wvl)
1161 : :
1162 : : !Close file
1163 : 20 : call WffClose(wff1,ierr)
1164 : :
1165 : : end if
1166 : :
1167 : 50 : end subroutine read_wfk
1168 : : !!***
1169 : :
1170 [ # # # # : 0 : end module m_rttddft_tdks
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # # #
# # # ]
1171 : : !!***
|
