Branch data Line data Source code
1 : : !!****m* ABINIT/m_rttddft_output
2 : : !! NAME
3 : : !! m_rttddft_ouptut
4 : : !!
5 : : !! FUNCTION
6 : : !! Manages most output of RT-TDDFT runs
7 : : !!
8 : : !! COPYRIGHT
9 : : !! Copyright (C) 2021-2026 ABINIT group (FB)
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 : : !!
14 : : !! SOURCE
15 : :
16 : : #if defined HAVE_CONFIG_H
17 : : #include "config.h"
18 : : #endif
19 : :
20 : : #include "abi_common.h"
21 : :
22 : : module m_rttddft_output
23 : :
24 : : use netcdf
25 : :
26 : : use defs_basis
27 : : use defs_abitypes, only: MPI_type
28 : : use defs_datatypes, only: pseudopotential_type
29 : : use m_common, only: prteigrs
30 : : use m_crystal, only: crystal_t
31 : : use m_dtfil, only: datafiles_type
32 : : use m_dtset, only: dataset_type
33 : : use m_ebands, only: ebands_t
34 : : use m_epjdos, only: epjdos_t
35 : : use m_errors, only: msg_hndl, assert
36 : : use m_ioarr, only: fftdatar_write
37 : : use m_io_tools, only: open_file, close_unit
38 : : use m_iowf, only: outwf
39 : : use m_mpinfo, only: iwrite_fftdatar
40 : : use m_paral_atom, only: get_my_atmtab, free_my_atmtab
41 : : use m_profiling_abi, only: abimem_record
42 : : use m_rttddft_tdks, only: tdks_type
43 : : use m_specialmsg, only: wrtout
44 : : use m_xmpi, only: xmpi_comm_rank
45 : :
46 : : implicit none
47 : :
48 : : private
49 : : !!***
50 : :
51 : : public :: rttddft_output
52 : : !!***
53 : :
54 : : contains
55 : : !!***
56 : :
57 : : !!****f* m_rttddft_output/rttddft_output
58 : : !!
59 : : !! NAME
60 : : !! rttddft_output
61 : : !!
62 : : !! FUNCTION
63 : : !! Main output subroutine
64 : : !!
65 : : !! INPUTS
66 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
67 : : !! dtset <type(dataset_type)> = all input variables for this dataset
68 : : !! istep <integer> = step number
69 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
70 : : !! psps <type(pseudopotential_type)> = variables related to pseudopotentials
71 : : !! tdks <type(tdks_type)> = the tdks object to initialize
72 : : !!
73 : : !! OUTPUT
74 : : !!
75 : : !! SOURCE
76 : 1725 : subroutine rttddft_output(dtfil, dtset, istep, mpi_enreg, psps, tdks)
77 : :
78 : : !Arguments ------------------------------------
79 : : !scalars
80 : : integer, intent(in) :: istep
81 : : type(datafiles_type), intent(inout) :: dtfil
82 : : type(dataset_type), intent(inout) :: dtset
83 : : type(MPI_type), intent(inout) :: mpi_enreg
84 : : type(pseudopotential_type), intent(inout) :: psps
85 : : type(tdks_type), intent(inout) :: tdks
86 : : !arrays
87 : :
88 : : !Local variables-------------------------------
89 : : !scalars
90 : : character(len=500) :: msg
91 : : character(len=100) :: fmt
92 : : character(len=20) :: access
93 : : integer :: i
94 : : !arrays
95 : :
96 : : ! *************************************************************************
97 : :
98 : : !** Special case of first step
99 [ + + ]: 1725 : if (istep == tdks%first_step) then
100 : 50 : access = "sequential"
101 [ + + ]: 50 : if (dtset%td_restart>0) access = "append"
102 : :
103 : : ! Open energy file and writes header if needed
104 [ + + ]: 50 : if (mpi_enreg%me == 0) then
105 [ - + ]: 22 : if (open_file(tdks%fname_tdener,msg,newunit=tdks%tdener_unit,status='unknown',form='formatted',access=access) /= 0) ABI_ERROR(msg)
106 : : end if
107 [ + + ]: 50 : if (dtset%td_restart==0) then
108 : 30 : write(msg,'(a)') "# RT-TDDFT -- Energy file. All quantities are in Hartree atomic units."
109 : 30 : call wrtout(tdks%tdener_unit,msg)
110 : : write(msg,'(a)') "# step time E_total E_kinetic E_hartree E_xc E_ewald &
111 : 30 : & E_corepsp E_localpsp E_nonlocalpsp E_paw E_entropy E_vdw"
112 : 30 : call wrtout(tdks%tdener_unit,msg)
113 : : end if
114 : :
115 : : ! Open electric field file and writes header if needed
116 [ + + ]: 50 : if (dtset%td_ef_type /= 0) then
117 [ + + ]: 21 : if (mpi_enreg%me == 0) then
118 [ - + ]: 7 : if (open_file(tdks%fname_tdef,msg,newunit=tdks%tdef_unit,status='unknown',form='formatted',access=access) /= 0) ABI_ERROR(msg)
119 : : end if
120 [ + + ]: 21 : if (dtset%td_restart==0) then
121 : 11 : write(msg,'(a)') "# RT-TDDFT -- Electric field file. All quantities are in Hartree atomic units."
122 : 11 : call wrtout(tdks%tdef_unit,msg)
123 : 11 : write(msg,'(a)') "# step time E_x E_y E_z A_x A_y A_z A_ext_x A_ext_y A_ext_z A_ind_x A_ind_y A_ind_z"
124 : 11 : call wrtout(tdks%tdef_unit,msg)
125 : : end if
126 : : end if
127 : :
128 : : ! Open current file and writes header if needed
129 [ + + ]: 50 : if (dtset%prtcurrent /= 0) then
130 [ + + ]: 21 : if (mpi_enreg%me == 0) then
131 [ - + ]: 7 : if (open_file(tdks%fname_current,msg,newunit=tdks%current_unit,status='unknown',form='formatted',access=access) /= 0) ABI_ERROR(msg)
132 : : end if
133 [ + + ]: 21 : if (dtset%td_restart==0) then
134 : 11 : write(msg,'(a)') "# RT-TDDFT -- Current density file. All quantities are in Hartree atomic units."
135 : 11 : call wrtout(tdks%current_unit,msg)
136 : 11 : write(msg,'(a)') "# step time J_x J_y J_z"
137 : 11 : call wrtout(tdks%current_unit,msg)
138 : : end if
139 : : end if
140 : : end if
141 : :
142 : : !** Writes some info in main output file
143 : 1725 : write(msg,'(a,a,f14.6,a)') ch10,'Total energy = ', tdks%etot,' Ha'
144 : 1725 : call wrtout(ab_out,msg)
145 [ + + ]: 1725 : if (do_write_log) call wrtout(std_out,msg)
146 : :
147 : 1725 : write(msg,'(a,f18.10,a)') 'Integrated density (ie. total nb of electrons) = ', &
148 [ + + ]: 29326665 : & sum(tdks%rhor(:,1))*tdks%ucvol/tdks%nfftf, ch10
149 : 1725 : call wrtout(ab_out,msg)
150 [ + + ]: 1725 : if (do_write_log) call wrtout(std_out,msg)
151 : :
152 : : !** Writes in energy file
153 : 1725 : write(msg,'(i0,1X,f15.5,11(f14.8,1X))') istep-1, (istep-1)*tdks%dt, tdks%etot, tdks%energies%e_kinetic, &
154 : 1725 : & tdks%energies%e_hartree, tdks%energies%e_xc, tdks%energies%e_ewald, &
155 : 1725 : & tdks%energies%e_corepsp, tdks%energies%e_localpsp, tdks%energies%e_nlpsp_vfock, &
156 : 3450 : & tdks%energies%paw%epaw, tdks%energies%e_entropy, tdks%energies%e_vdw_dftd
157 : 1725 : call wrtout(tdks%tdener_unit,msg)
158 : :
159 : : !** Writes TD elec. field and associated vector potential if needed
160 [ + + ]: 1725 : if (dtset%td_ef_type /= 0) then
161 : 1100 : write(msg,'(i0,1X,f15.5,1X,12(f14.8,1X))') istep, istep*tdks%dt, tdks%tdef%efield(:), tdks%tdef%vecpot(:), &
162 : 2200 : & tdks%tdef%vecpot_ext(:), tdks%tdef%vecpot_ind(:,1)
163 : 1100 : call wrtout(tdks%tdef_unit,msg)
164 : : end if
165 : :
166 : : !** Writes TD current density if needed
167 [ + + ]: 1725 : if (dtset%prtcurrent /= 0) then
168 [ + - ]: 1100 : if (dtset%nsppol == 1) then
169 : 1100 : fmt = '(i0,1X,f15.5,1X,3(f14.8,1X))'
170 : : else
171 : 0 : fmt = '(i0,1X,f15.5,1X,3(f14.8,1X),3(f14.8,1X))'
172 : : end if
173 [ + - + + ]: 2200 : write(msg,fmt) istep, istep*tdks%dt, (tdks%current(:,i),i=1,dtset%nsppol)
174 : 1100 : call wrtout(tdks%current_unit,msg)
175 : : end if
176 : :
177 : : !** Writes additional optional properties
178 : : !Update header, with evolving variables
179 : : call tdks%hdr%update(tdks%bantot,tdks%etot,tdks%energies%e_fermie,tdks%energies%e_fermih, &
180 : : & tdks%hdr%residm,tdks%rprimd,tdks%occ0,tdks%pawrhoij, &
181 : : & tdks%xred,dtset%amu_orig,comm_atom=mpi_enreg%comm_atom, &
182 : 1725 : & mpi_atmtab=mpi_enreg%my_atmtab)
183 : :
184 : : !Computed at actual step
185 [ + - ]: 1725 : if (mod(istep,dtset%td_prtstr) == 0) then
186 : 1725 : call prt_den(dtfil,dtset,istep,mpi_enreg,psps,tdks)
187 [ - + ]: 1725 : if (dtset%prtwf > 0) then
188 : 0 : call prt_wfk(dtfil,dtset,istep,mpi_enreg,psps,tdks)
189 : 0 : call prt_restart(dtfil,dtset,istep,mpi_enreg,tdks)
190 : : end if
191 : : end if
192 : :
193 : : !Computed at previous step
194 [ + - ]: 1725 : if (mod(istep-1,dtset%td_prtstr) == 0) then
195 : 1725 : call prt_eig(dtfil,dtset,istep-1,mpi_enreg,tdks)
196 : 1725 : call prt_occ(dtfil,dtset,istep-1,mpi_enreg,tdks)
197 : 1725 : call prt_dos(dtfil,dtset,istep-1,mpi_enreg,psps,tdks)
198 : : end if
199 : :
200 : : !** Special case of last step
201 [ + + ]: 1725 : if (istep == tdks%first_step+tdks%ntime-1) then
202 [ + - + - ]: 49 : if (mod(istep,dtset%td_prtstr) /= 0 .or. dtset%prtwf <= 0) then
203 : 49 : call prt_wfk(dtfil,dtset,istep,mpi_enreg,psps,tdks,force_write=.TRUE.)
204 : 49 : call prt_restart(dtfil,dtset,istep,mpi_enreg,tdks)
205 : : end if
206 [ + + ]: 49 : if (mpi_enreg%me == 0) then
207 : : !close all files
208 [ - + ]: 21 : if (close_unit(tdks%tdener_unit,msg) /= 0) ABI_ERROR(msg)
209 [ - + ]: 21 : if (close_unit(tdks%tdrestart_unit,msg) /= 0) ABI_ERROR(msg)
210 [ + + ]: 21 : if (dtset%td_ef_type /= 0) then
211 [ - + ]: 7 : if (close_unit(tdks%tdef_unit,msg) /= 0) ABI_ERROR(msg)
212 : : end if
213 [ + + ]: 21 : if (dtset%prtcurrent /= 0) then
214 [ - + ]: 7 : if (close_unit(tdks%current_unit,msg) /= 0) ABI_ERROR(msg)
215 : : end if
216 : : end if
217 : : end if
218 : :
219 : 1725 : end subroutine rttddft_output
220 : : !!***
221 : :
222 : : !!****f* m_rttddft_output/prt_eig
223 : : !!
224 : : !! NAME
225 : : !! prt_eig
226 : : !!
227 : : !! FUNCTION
228 : : !! Outputs eigenvalues
229 : : !!
230 : : !! INPUTS
231 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
232 : : !! dtset <type(dataset_type)> = all input variables for this dataset
233 : : !! istep <integer> = step number
234 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
235 : : !! tdks <type(tdks_type)> = the tdks object to initialize
236 : : !!
237 : : !! OUTPUT
238 : : !!
239 : : !! SIDE EFFECTS
240 : : !!
241 : : !! SOURCE
242 : 1725 : subroutine prt_eig(dtfil, dtset, istep, mpi_enreg, tdks)
243 : :
244 : : !Arguments ------------------------------------
245 : : !scalars
246 : : integer, intent(in) :: istep
247 : : type(datafiles_type), intent(inout) :: dtfil
248 : : type(dataset_type), intent(inout) :: dtset
249 : : type(MPI_type), intent(inout) :: mpi_enreg
250 : : type(tdks_type), intent(inout) :: tdks
251 : : !arrays
252 : :
253 : : !Local variables-------------------------------
254 : : !scalars
255 : : integer,parameter :: enunit=0, option=3
256 : : integer :: me
257 : : integer :: spacecomm
258 : : real(dp) :: vxcavg_dum
259 : : character(len=fnlen) :: fname
260 : : character(len=24) :: step_nb
261 : : !arrays
262 : 3450 : real(dp) :: resid(dtset%mband*dtset%nkpt*dtset%nsppol)
263 : :
264 : : ! *************************************************************************
265 : :
266 : 1725 : spacecomm = mpi_enreg%comm_cell
267 : 1725 : me = xmpi_comm_rank(spacecomm)
268 : :
269 : : !to avoid some useless printing in the output file when both prteig and prtvol are set to 0
270 [ + - - + ]: 1725 : if (dtset%prteig /= 0 .or. dtset%prtvol /= 0) then
271 : :
272 : 0 : write(step_nb,*) istep
273 [ # # # # : 0 : fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_EIG'
# # ]
274 [ # # ]: 0 : resid = zero
275 : 0 : vxcavg_dum=zero
276 : :
277 [ # # ]: 0 : if(me==0)then
278 : : call prteigrs(tdks%eigen,enunit,tdks%energies%e_fermie,tdks%energies%e_fermih, &
279 : : & fname,ab_out,dtset%iscf,dtset%kptns,dtset%kptopt,dtset%mband, &
280 : : & dtset%nband,dtset%nbdbuf,dtset%nkpt,0,dtset%nsppol,tdks%occ0, &
281 : : & dtset%occopt,option,dtset%prteig,dtset%prtvol,resid,dtset%tolwfr, &
282 : 0 : & vxcavg_dum,dtset%wtk)
283 : : end if
284 : : end if
285 : :
286 : 1725 : end subroutine prt_eig
287 : : !!***
288 : :
289 : : !!****f* m_rttddft_output/prt_occ
290 : : !!
291 : : !! NAME
292 : : !! prt_occ
293 : : !!
294 : : !! FUNCTION
295 : : !! Outputs occupation numbers
296 : : !!
297 : : !! INPUTS
298 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
299 : : !! dtset <type(dataset_type)> = all input variables for this dataset
300 : : !! istep <integer> = step number
301 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
302 : : !! tdks <type(tdks_type)> = the tdks object to initialize
303 : : !!
304 : : !! OUTPUT
305 : : !!
306 : : !! SOURCE
307 : 1725 : subroutine prt_occ(dtfil, dtset, istep, mpi_enreg, tdks)
308 : :
309 : : !Arguments ------------------------------------
310 : : !scalars
311 : : integer, intent(in) :: istep
312 : : type(datafiles_type), intent(inout) :: dtfil
313 : : type(dataset_type), intent(inout) :: dtset
314 : : type(MPI_type), intent(inout) :: mpi_enreg
315 : : type(tdks_type), intent(inout) :: tdks
316 : : !arrays
317 : :
318 : : !Local variables-------------------------------
319 : : !scalars
320 : : integer :: band_index
321 : : integer :: iband, ii, ikpt, isppol
322 : : integer :: me
323 : : integer :: nkpt
324 : : integer :: nband_k, nsppol
325 : : integer :: temp_unit
326 : : !arrays
327 : : character(len=fnlen) :: fname
328 : : character(len=4) :: ibnd_fmt, ikpt_fmt
329 : : character(len=500) :: msg
330 : : character(len=24) :: step_nb
331 : :
332 : :
333 : : ! *************************************************************************
334 : :
335 [ - + ]: 1725 : if (dtset%prtocc > 0) then
336 : 0 : me = xmpi_comm_rank(mpi_enreg%comm_cell)
337 : :
338 : 0 : write(step_nb,*) istep
339 [ # # # # : 0 : fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_OCC'
# # ]
340 : :
341 [ # # ]: 0 : if (open_file(fname, msg, newunit=temp_unit, status='unknown', form='formatted') /= 0) then
342 : 0 : ABI_ERROR(msg)
343 : : end if
344 : :
345 : 0 : nkpt = dtset%nkpt
346 : 0 : nsppol = dtset%nsppol
347 : :
348 [ # # ]: 0 : if(me==0)then
349 : : band_index=0
350 [ # # ]: 0 : do isppol=1,nsppol
351 : :
352 [ # # ]: 0 : if(nsppol==2)then
353 [ # # ]: 0 : if(isppol==1)write(msg, '(2a)' ) ch10,' SPIN UP channel '
354 [ # # ]: 0 : if(isppol==2)write(msg, '(2a)' ) ch10,' SPIN DOWN channel '
355 : 0 : call wrtout(temp_unit,msg)
356 : : end if
357 [ # # # # ]: 0 : ikpt_fmt="i4" ; if(nkpt>=10000)ikpt_fmt="i6" ; if(nkpt>=1000000)ikpt_fmt="i9"
358 [ # # ]: 0 : if (nsppol==2.and.isppol==1) then
359 : : write(msg, '(a,'//ikpt_fmt//',2x,a)' ) &
360 : 0 : 'Occupation numbers for nkpt=',nkpt,'k points, SPIN UP:'
361 [ # # ]: 0 : else if (nsppol==2.and.isppol==2) then
362 : : write(msg, '(a,'//ikpt_fmt//',2x,a)' ) &
363 : 0 : 'Occupation numbers for nkpt=',nkpt,'k points, SPIN DOWN:'
364 : : else
365 : : write(msg, '(a,'//ikpt_fmt//',2x,a)' ) &
366 : 0 : 'Occupation numbers for nkpt=',nkpt,'k points:'
367 : : end if
368 : 0 : call wrtout(temp_unit,msg)
369 [ # # ]: 0 : do ikpt=1,nkpt
370 : 0 : nband_k=dtset%nband(ikpt+(isppol-1)*nkpt)
371 [ # # # # ]: 0 : ikpt_fmt="i4" ; if(nkpt>=10000)ikpt_fmt="i6" ; if(nkpt>=1000000)ikpt_fmt="i9"
372 [ # # # # ]: 0 : ibnd_fmt="i3" ; if(nband_k>=1000)ibnd_fmt="i6" ; if(nband_k>=1000000)ibnd_fmt="i9"
373 : : write(msg, '(a,'//ikpt_fmt//',a,'//ibnd_fmt//',a,f9.5,a,3f8.4,a)' ) &
374 : 0 : & ' kpt#',ikpt,', nband=',nband_k,', wtk=',dtset%wtk(ikpt)+tol10,', kpt=',&
375 [ # # ]: 0 : & dtset%kptns(1:3,ikpt)+tol10,' (reduced coord)'
376 : 0 : call wrtout(temp_unit,msg)
377 [ # # ]: 0 : do ii=0,(nband_k-1)/6
378 [ # # # # ]: 0 : write(msg, '(1p,6e12.4)')(tdks%occ(iband+band_index),iband=1+6*ii,min(6+6*ii,nband_k))
379 : 0 : call wrtout(temp_unit,msg)
380 : : end do
381 : 0 : band_index=band_index+nband_k
382 : : end do
383 : : end do
384 : : end if
385 : :
386 : 0 : close(temp_unit)
387 : : end if
388 : :
389 : 1725 : end subroutine prt_occ
390 : : !!***
391 : :
392 : : !!****f* m_rttddft_output/prt_den
393 : : !!
394 : : !! NAME
395 : : !! prt_den
396 : : !!
397 : : !! FUNCTION
398 : : !! Outputs the electronic density
399 : : !!
400 : : !! INPUTS
401 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
402 : : !! dtset <type(dataset_type)> = all input variables for this dataset
403 : : !! istep <integer> = step number
404 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
405 : : !! psps <type(pseudopotential_type)> = variables related to pseudopotentials
406 : : !! tdks <type(tdks_type)> = the tdks object to initialize
407 : : !!
408 : : !! OUTPUT
409 : : !!
410 : : !! SOURCE
411 : 1725 : subroutine prt_den(dtfil, dtset, istep, mpi_enreg, psps, tdks)
412 : :
413 : : !Arguments ------------------------------------
414 : : !scalars
415 : : integer, intent(in) :: istep
416 : : type(datafiles_type), intent(inout) :: dtfil
417 : : type(dataset_type), intent(inout) :: dtset
418 : : type(MPI_type), intent(inout) :: mpi_enreg
419 : : type(pseudopotential_type), intent(inout) :: psps
420 : : type(tdks_type), intent(inout) :: tdks
421 : : !arrays
422 : :
423 : : !Local variables-------------------------------
424 : : !scalars
425 : : integer,parameter :: cplex1=1
426 : : integer :: bantot
427 : : integer :: iatom
428 : : integer :: spacecomm
429 : : integer :: my_comm_atom, my_natom
430 : : integer :: me
431 : : integer :: natom
432 : : ! integer :: ncid
433 : : integer :: timrev
434 : : character(len=fnlen) :: fname
435 : : character(len=24) :: step_nb
436 : : logical :: paral_atom
437 : : logical :: remove_inv
438 : : logical :: my_atmtab_allocated
439 [ + + + + : 108675 : type(crystal_t) :: crystal
+ + + + +
+ + + + +
+ + + + +
+ ]
440 : 1725 : type(ebands_t) :: ebands
441 : : !arrays
442 : 1725 : integer, pointer :: my_atmtab(:)
443 : 1725 : real(dp), allocatable :: doccde(:)
444 : :
445 : : ! *************************************************************************
446 : :
447 : 1725 : if (dtset%prtden /= 0) then
448 : 0 : spacecomm = mpi_enreg%comm_cell
449 : 0 : me = xmpi_comm_rank(spacecomm)
450 : :
451 : 0 : natom = dtset%natom
452 : 0 : my_natom = mpi_enreg%my_natom
453 : 0 : paral_atom=(my_natom/=natom)
454 : 0 : my_comm_atom = mpi_enreg%comm_atom
455 : 0 : nullify(my_atmtab)
456 [ # # ]: 0 : if (paral_atom) then
457 : 0 : call get_my_atmtab(mpi_enreg%comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
458 : : else
459 [ # # # # : 0 : ABI_MALLOC(my_atmtab, (natom))
# # # # #
# # # ]
460 [ # # # # : 0 : my_atmtab = (/ (iatom, iatom=1, natom) /)
# # ]
461 : 0 : my_atmtab_allocated = .true.
462 : : end if
463 : :
464 : : !FB: Maybe this should be moved out of that subroutine if needed in other outputs than densities
465 : 0 : remove_inv=.false.
466 [ # # ]: 0 : timrev = 2; if (any(dtset%kptopt == [3, 4])) timrev= 1
467 : : call crystal%init(dtset%amu_orig(:,1),dtset%spgroup,natom,dtset%npsp,psps%ntypat, &
468 : : dtset%nsym,tdks%rprimd,dtset%typat,tdks%xred,dtset%ziontypat,dtset%znucl,timrev,&
469 : : dtset%nspden==2.and.dtset%nsppol==1,remove_inv,tdks%hdr%title,&
470 [ # # # # ]: 0 : dtset%symrel,dtset%tnons,dtset%symafm)
471 : : !Electron band energies.
472 : 0 : bantot= dtset%mband*dtset%nkpt*dtset%nsppol
473 [ # # # # : 0 : ABI_CALLOC(doccde, (bantot))
# # # # #
# # # #
# ]
474 : : call ebands%init(bantot, dtset%nelect,dtset%ne_qFD,dtset%nh_qFD,dtset%ivalence, &
475 : : doccde,tdks%eigen,dtset%istwfk,dtset%kptns,dtset%nband,dtset%nkpt,tdks%npwarr,dtset%nsppol, &
476 : : dtset%nspinor,dtset%tphysel,dtset%tsmear,dtset%occopt,tdks%occ0,dtset%wtk,&
477 : : dtset%cellcharge(1),dtset%kptopt,dtset%kptrlatt_orig,dtset%nshiftk_orig,dtset%shiftk_orig, &
478 : 0 : dtset%kptrlatt,dtset%nshiftk,dtset%shiftk)
479 [ # # ]: 0 : ABI_FREE(doccde)
480 : :
481 : 0 : write(step_nb,*) istep
482 : :
483 : : !** Outputs the density
484 : : !Warnings :
485 : : !- core charge is excluded from the charge density;
486 : : !- the potential is the INPUT vtrial.
487 [ # # ]: 0 : if (iwrite_fftdatar(mpi_enreg)) then
488 [ # # # # : 0 : fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_DEN'
# # ]
489 : : call fftdatar_write("density",fname,dtset%iomode,tdks%hdr,crystal,tdks%pawfgr%ngfft, &
490 : 0 : & cplex1,tdks%pawfgr%nfft,dtset%nspden,tdks%rhor,mpi_enreg,ebands=ebands)
491 : : end if
492 : :
493 : 0 : call crystal%free()
494 : 0 : call ebands%free()
495 [ # # ]: 0 : ABI_FREE(my_atmtab)
496 : : end if
497 : :
498 [ - + - + : 1725 : end subroutine prt_den
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ ]
499 : : !!***
500 : :
501 : : !!****f* m_rttddft_output/prt_dos
502 : : !!
503 : : !! NAME
504 : : !! prt_dos
505 : : !!
506 : : !! FUNCTION
507 : : !! Computes and outputs the electronic DOS
508 : : !!
509 : : !! INPUTS
510 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
511 : : !! dtset <type(dataset_type)> = all input variables for this dataset
512 : : !! istep <integer> = step number
513 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
514 : : !! psps <type(pseudopotential_type)> = variables related to pseudopotentials
515 : : !! tdks <type(tdks_type)> = the tdks object to initialize
516 : : !!
517 : : !! OUTPUT
518 : : !!
519 : : !! SOURCE
520 : 1725 : subroutine prt_dos(dtfil, dtset, istep, mpi_enreg, psps, tdks)
521 : :
522 : : !Arguments ------------------------------------
523 : : !scalars
524 : : integer, intent(in) :: istep
525 : : type(datafiles_type), intent(inout) :: dtfil
526 : : type(dataset_type), intent(inout) :: dtset
527 : : type(MPI_type), intent(inout) :: mpi_enreg
528 : : type(pseudopotential_type), intent(inout) :: psps
529 : : type(tdks_type), intent(inout) :: tdks
530 : : !arrays
531 : :
532 : : !Local variables-------------------------------
533 : : !scalars
534 : : integer,parameter :: master=0
535 : : integer :: bantot
536 : : integer :: collect
537 : : integer :: iatom
538 : : integer :: spacecomm
539 : : integer :: my_comm_atom, my_natom
540 : : integer :: me
541 : : integer :: natom
542 : : ! integer :: ncid
543 : : integer :: timrev
544 : : character(len=fnlen) :: fname
545 : : character(len=24) :: step_nb
546 : : logical :: paral_atom
547 : : logical :: remove_inv
548 : : logical :: my_atmtab_allocated
549 [ + + + + : 117300 : type(crystal_t) :: crystal
+ + + + +
+ + + + +
+ + + + +
+ ]
550 : 1725 : type(epjdos_t) :: dos
551 : 1725 : type(ebands_t) :: ebands
552 : : !arrays
553 : 1725 : integer, pointer :: my_atmtab(:)
554 : 1725 : real(dp), allocatable :: doccde(:)
555 : :
556 : : ! *************************************************************************
557 : :
558 : 1725 : spacecomm = mpi_enreg%comm_cell
559 : 3450 : me = xmpi_comm_rank(spacecomm)
560 : :
561 : : !FB: @MT - Is this needed?
562 : 1725 : natom = dtset%natom
563 : 1725 : my_natom = mpi_enreg%my_natom
564 : 1725 : paral_atom=(my_natom/=natom)
565 : 1725 : my_comm_atom = mpi_enreg%comm_atom
566 : 1725 : nullify(my_atmtab)
567 [ - + ]: 1725 : if (paral_atom) then
568 : 0 : call get_my_atmtab(mpi_enreg%comm_atom,my_atmtab,my_atmtab_allocated,paral_atom,natom,my_natom_ref=my_natom)
569 : : else
570 [ - + + - : 5175 : ABI_MALLOC(my_atmtab, (natom))
+ - + - -
+ - + ]
571 [ + - + + : 11310 : my_atmtab = (/ (iatom, iatom=1, natom) /)
+ + ]
572 : 1725 : my_atmtab_allocated = .true.
573 : : end if
574 : :
575 : 1725 : remove_inv=.false.
576 [ + + ]: 1725 : timrev = 2; if (any(dtset%kptopt == [3, 4])) timrev= 1
577 : : call crystal%init(dtset%amu_orig(:,1),dtset%spgroup,natom,dtset%npsp,psps%ntypat, &
578 : : dtset%nsym,tdks%rprimd,dtset%typat,tdks%xred,dtset%ziontypat,dtset%znucl,timrev,&
579 : : dtset%nspden==2.and.dtset%nsppol==1,remove_inv,tdks%hdr%title,&
580 [ - + - - ]: 3450 : dtset%symrel,dtset%tnons,dtset%symafm)
581 : : !Electron band energies.
582 : 1725 : bantot= dtset%mband*dtset%nkpt*dtset%nsppol
583 [ - + + - : 398935 : ABI_CALLOC(doccde, (bantot))
+ - + - -
+ - + +
+ ]
584 : : call ebands%init(bantot,dtset%nelect,dtset%ne_qFD,dtset%nh_qFD,dtset%ivalence, &
585 : : doccde,tdks%eigen,dtset%istwfk,dtset%kptns,dtset%nband,dtset%nkpt,tdks%npwarr,dtset%nsppol, &
586 : : dtset%nspinor,dtset%tphysel,dtset%tsmear,dtset%occopt,tdks%occ0,dtset%wtk,&
587 : : dtset%cellcharge(1),dtset%kptopt,dtset%kptrlatt_orig,dtset%nshiftk_orig,dtset%shiftk_orig, &
588 : 1725 : dtset%kptrlatt,dtset%nshiftk,dtset%shiftk)
589 [ - + ]: 1725 : ABI_FREE(doccde)
590 : :
591 : 1725 : write(step_nb,*) istep
592 : :
593 : : !** Generate DOS using the tetrahedron method or using Gaussians
594 [ + - - + ]: 1725 : if (dtset%prtdos>=2.or.dtset%pawfatbnd>0) then
595 [ # # # # : 0 : call dos%init(dtset, psps, tdks%pawtab)
# # # # ]
596 : :
597 [ # # # # ]: 0 : if (dos%partial_dos_flag>=1 .or. dos%fatbands_flag==1)then
598 : : ! Generate fractions for partial DOSs if needed partial_dos 1,2,3,4 give different decompositions
599 : 0 : collect = 1 !; if (psps%usepaw==1 .and. dos%partial_dos_flag /= 2) collect = 0
600 [ # # # # : 0 : if ((psps%usepaw==0.or.dtset%pawprtdos/=2) .and. dos%partial_dos_flag>=1) then
# # ]
601 : 0 : call dos%partial_dos_fractions(crystal,dtset,tdks%eigen,tdks%occ0,tdks%npwarr,tdks%kg,tdks%cg,tdks%mcg,collect,mpi_enreg)
602 : : end if
603 : :
604 [ # # # # ]: 0 : if (psps%usepaw==1 .and. dos%partial_dos_flag /= 2) then
605 : : ! TODO: update partial_dos_fractions_paw for extra atoms - no PAW contribution normally, but check bounds and so on.
606 [ # # # # : 0 : call dos%partial_dos_fractions_paw(tdks%cprj,tdks%dimcprj,dtset,tdks%mcprj,dtset%mkmem,mpi_enreg,tdks%pawrad,tdks%pawtab)
# # # # #
# # # # #
# # ]
607 : : end if
608 : : else
609 [ # # # # : 0 : dos%fractions(:,:,:,1)=one
# # ]
610 : : end if
611 : :
612 : : !Here, print out fatbands for the k-points given in file appended _FATBANDS
613 [ # # # # : 0 : if (me == master .and. dtset%pawfatbnd>0 .and. dos%fatbands_flag==1) then
# # ]
614 [ # # # # : 0 : fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_FATBANDS'
# # ]
615 [ # # # # : 0 : call dos%prtfatbands(dtset,ebands,fname,dtset%pawfatbnd,tdks%pawtab)
# # # # ]
616 : : end if
617 : :
618 : : !Here, computation and output of DOS and partial DOS _DOS
619 [ # # # # ]: 0 : if (dos%fatbands_flag == 0 .and. dos%prtdos /= 4) then
620 [ # # # # : 0 : fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_DOS'
# # ]
621 : 0 : call dos%calcnwrite(dtset,crystal,ebands,fname,spacecomm)
622 : : end if
623 : : end if
624 : :
625 : 1725 : call dos%free()
626 : 1725 : call crystal%free()
627 : 1725 : call ebands%free()
628 [ - + ]: 1725 : ABI_FREE(my_atmtab)
629 : :
630 [ - + - + : 3450 : end subroutine prt_dos
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + -
+ - + - +
- + - + ]
631 : : !!***
632 : :
633 : : !!****f* m_rttddft_output/prt_wfk
634 : : !!
635 : : !! NAME
636 : : !! prt_wfk
637 : : !!
638 : : !! FUNCTION
639 : : !! Outputs wavefunctions in WFK file
640 : : !!
641 : : !! INPUTS
642 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
643 : : !! dtset <type(dataset_type)> = all input variables for this dataset
644 : : !! istep <integer> = step number
645 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
646 : : !! psps <type(pseudopotential_type)> = variables related to pseudopotentials
647 : : !! tdks <type(tdks_type)> = the tdks object to initialize
648 : : !! force_write <logical> = force the writing of WFK (useful for last step) - optional
649 : : !!
650 : : !! OUTPUT
651 : : !!
652 : : !! SOURCE
653 : 49 : subroutine prt_wfk(dtfil, dtset, istep, mpi_enreg, psps, tdks, force_write)
654 : :
655 : : !Arguments ------------------------------------
656 : : !scalars
657 : : integer, intent(in) :: istep
658 : : type(datafiles_type), intent(inout) :: dtfil
659 : : type(dataset_type), intent(inout) :: dtset
660 : : type(MPI_type), intent(inout) :: mpi_enreg
661 : : type(pseudopotential_type), intent(inout) :: psps
662 : : type(tdks_type), intent(inout) :: tdks
663 : : logical, optional, intent(in) :: force_write
664 : : !arrays
665 : :
666 : : !Local variables-------------------------------
667 : : !scalars
668 : : integer,parameter :: response=0
669 : : character(len=fnlen) :: fname
670 : : character(len=24) :: step_nb
671 : : logical :: lforce_write = .FALSE.
672 : : !arrays
673 : :
674 : : ! *************************************************************************
675 : :
676 : 49 : write(step_nb,*) istep
677 [ + - + - : 49 : fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_WFK'
+ - ]
678 : :
679 [ + - ]: 49 : if (present(force_write)) then
680 [ + - ]: 49 : if (force_write) lforce_write = .TRUE.
681 : : end if
682 : :
683 : : !Use initial eigenvalues to ensure that we get the same occupation upon restart
684 [ + - ]: 49 : if (lforce_write) then
685 : : call outwf(tdks%cg,dtset,psps,tdks%eigen0,fname,tdks%hdr,tdks%kg,dtset%kptns, &
686 : : & dtset%mband,tdks%mcg,dtset%mkmem,mpi_enreg,dtset%mpw,dtset%natom, &
687 : : & dtset%nband,dtset%nkpt,tdks%npwarr,dtset%nsppol,tdks%occ0,response, &
688 : 49 : & dtfil%unwff2,tdks%wvl%wfs,tdks%wvl%descr, force_write=.TRUE.)
689 : : else
690 : : call outwf(tdks%cg,dtset,psps,tdks%eigen0,fname,tdks%hdr,tdks%kg,dtset%kptns, &
691 : : & dtset%mband,tdks%mcg,dtset%mkmem,mpi_enreg,dtset%mpw,dtset%natom, &
692 : : & dtset%nband,dtset%nkpt,tdks%npwarr,dtset%nsppol,tdks%occ0,response, &
693 : 0 : & dtfil%unwff2,tdks%wvl%wfs,tdks%wvl%descr)
694 : : end if
695 : :
696 : 49 : end subroutine prt_wfk
697 : : !!***
698 : :
699 : : !!****f* m_rttddft_output/prt_restart
700 : : !!
701 : : !! NAME
702 : : !! prt_restart
703 : : !!
704 : : !! FUNCTION
705 : : !! Print restart file
706 : : !!
707 : : !! INPUTS
708 : : !! dtfil <type datafiles_type> = infos about file names, file unit numbers
709 : : !! istep <integer> = step number
710 : : !! mpi_enreg <MPI_type> = MPI-parallelisation information
711 : : !! tdks <type(tdks_type)> = the tdks object to initialize
712 : : !!
713 : : !! OUTPUT
714 : : !!
715 : : !! SOURCE
716 : 49 : subroutine prt_restart(dtfil, dtset, istep, mpi_enreg, tdks)
717 : :
718 : : !Arguments ------------------------------------
719 : : !scalars
720 : : integer, intent(in) :: istep
721 : : type(datafiles_type), intent(inout) :: dtfil
722 : : type(dataset_type), intent(inout) :: dtset
723 : : type(MPI_type), intent(inout) :: mpi_enreg
724 : : type(tdks_type), intent(inout) :: tdks
725 : : !arrays
726 : :
727 : : !Local variables-------------------------------
728 : : !scalars
729 : : character(len=500) :: msg
730 : : character(len=fnlen) :: fname
731 : : character(len=24 ) :: step_nb
732 : : !arrays
733 : :
734 : : ! *************************************************************************
735 : :
736 : 49 : write(step_nb,*) istep
737 [ + + ]: 49 : if (mpi_enreg%me == 0) rewind(tdks%tdrestart_unit)
738 : 49 : write(msg,'(a)') step_nb
739 : 49 : call wrtout(tdks%tdrestart_unit,msg)
740 : 49 : write(msg,'(a)') trim(tdks%fname_wfk0)
741 : 49 : call wrtout(tdks%tdrestart_unit,msg)
742 [ + - + - : 49 : fname = trim(dtfil%filnam_ds(4))//'_'//trim(adjustl(step_nb))//'_WFK'
+ - ]
743 : 49 : write(msg,'(a)') trim(fname)
744 : 49 : call wrtout(tdks%tdrestart_unit,msg)
745 : 49 : write(msg,'(a)') trim(tdks%fname_tdener)
746 : 49 : call wrtout(tdks%tdrestart_unit,msg)
747 [ + + ]: 49 : if (dtset%td_ef_type /= 0) then
748 : 21 : write(msg,'(a)') trim(tdks%fname_tdef)
749 : 21 : call wrtout(tdks%tdrestart_unit,msg)
750 : : end if
751 [ + + ]: 49 : if (dtset%prtcurrent /= 0) then
752 : 21 : write(msg,'(a)') trim(tdks%fname_current)
753 : 21 : call wrtout(tdks%tdrestart_unit,msg)
754 : : end if
755 [ + + ]: 49 : if (dtset%td_ef_induced_vecpot /= 0) then
756 : 2 : write(msg,*) tdks%tdef%vecpot_ind(:,2)
757 : 2 : call wrtout(tdks%tdrestart_unit,msg)
758 : 2 : write(msg,*) tdks%tdef%vecpot_ind(:,1)
759 : 2 : call wrtout(tdks%tdrestart_unit,msg)
760 : : end if
761 : :
762 : 49 : end subroutine prt_restart
763 : : !!***
764 : :
765 : : end module m_rttddft_output
766 : : !!***
|
