Octopus
td.F90
Go to the documentation of this file.
1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch
2!!
3!! This program is free software; you can redistribute it and/or modify
4!! it under the terms of the GNU General Public License as published by
5!! the Free Software Foundation; either version 2, or (at your option)
6!! any later version.
7!!
8!! This program is distributed in the hope that it will be useful,
9!! but WITHOUT ANY WARRANTY; without even the implied warranty of
10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11!! GNU General Public License for more details.
12!!
13!! You should have received a copy of the GNU General Public License
14!! along with this program; if not, write to the Free Software
15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
16!! 02110-1301, USA.
17!!
18
19#include "global.h"
20
21module td_oct_m
27 use debug_oct_m
32 use epot_oct_m
34 use forces_oct_m
36 use global_oct_m
37 use grid_oct_m
40 use io_oct_m
42 use ions_oct_m
43 use kick_oct_m
44 use, intrinsic :: iso_fortran_env
45 use lasers_oct_m
46 use lda_u_oct_m
49 use loct_oct_m
51 use mesh_oct_m
53 use mpi_oct_m
56 use output_oct_m
58 use parser_oct_m
59 use pes_oct_m
69 use scf_oct_m
71 use space_oct_m
75 use stress_oct_m
77 use types_oct_m
78 use unit_oct_m
80 use v_ks_oct_m
83 use xc_oct_m
84
85 implicit none
86
87 private
88 public :: &
89 td_t, &
90 td_run, &
92 td_init, &
94 td_end, &
95 td_end_run, &
98 td_dump, &
106
108 integer, parameter, public :: &
109 EHRENFEST = 1, &
110 bo = 2
111
112 type td_t
113 private
114 type(propagator_base_t), public :: tr
115 type(scf_t), public :: scf
116 type(ion_dynamics_t), public :: ions_dyn
117 real(real64), public :: dt
118 integer, public :: max_iter
119 integer, public :: iter
120 logical, public :: recalculate_gs
121
122 type(pes_t), public :: pesv
123
124 integer, public :: dynamics
125 integer, public :: energy_update_iter
126 real(real64) :: scissor
127
128 logical :: freeze_occ
129 logical :: freeze_u
130 integer :: freeze_orbitals
131
132 logical :: from_scratch = .false.
133
134 type(td_write_t), public :: write_handler
135 type(restart_t) :: restart_load
136 type(restart_t) :: restart_dump
137 end type td_t
138
139
140contains
141
142 subroutine td_run_init()
143
144 push_sub(td_run_init)
145
146 call calc_mode_par%set_parallelization(p_strategy_states, default = .true.)
147
148 pop_sub(td_run_init)
149 end subroutine td_run_init
150
151 ! ---------------------------------------------------------
152
153 subroutine td_init(td, namespace, space, gr, ions, st, ks, hm, ext_partners, outp)
154 type(td_t), intent(inout) :: td
155 type(namespace_t), intent(in) :: namespace
156 class(space_t), intent(in) :: space
157 type(grid_t), intent(in) :: gr
158 type(ions_t), intent(inout) :: ions
159 type(states_elec_t), intent(in) :: st
160 type(v_ks_t), intent(in) :: ks
161 type(hamiltonian_elec_t), intent(in) :: hm
162 type(partner_list_t), intent(in) :: ext_partners
163 type(output_t), intent(in) :: outp
164
165 integer :: default
166 real(real64) :: propagation_time
167 type(lasers_t), pointer :: lasers
168
169 push_sub(td_init)
170
171 if (hm%pcm%run_pcm) call messages_experimental("PCM for CalculationMode = td", namespace=namespace)
172
173 call ion_dynamics_init(td%ions_dyn, namespace, ions)
174
175 if (ion_dynamics_ions_move(td%ions_dyn)) then
176 if (hm%kpoints%use_symmetries) then
177 message(1) = "KPoints symmetries cannot be used with moving ions."
178 message(2) = "Please set KPointsSymmetries = no."
179 call messages_fatal(2, namespace=namespace)
180 end if
181 if (st%symmetrize_density) then
182 message(1) = "Symmetrization of the density cannot be used with moving ions."
183 message(2) = "Please set SymmetrizeDensity = no."
184 call messages_fatal(2, namespace=namespace)
185 end if
186 end if
187
188 td%iter = 0
189
190 !%Variable TDTimeStep
191 !%Type float
192 !%Section Time-Dependent::Propagation
193 !%Description
194 !% The time-step for the time propagation. For most propagators you
195 !% want to use the largest value that is possible without the
196 !% evolution becoming unstable.
197 !%
198 !% While prior versions of Octopus used to have a default time step, however,
199 !% this now needs to be systematically defined.
200 !%End
202 call parse_variable(namespace, 'TDTimeStep', -m_one, td%dt, unit = units_inp%time)
203
204 if (td%dt <= m_zero) then
205 write(message(1),'(a)') 'A positive value for TDTimeStep must be defined in the input file.'
206 call messages_fatal(1, namespace=namespace)
207 end if
209 call messages_print_var_value('TDTimeStep', td%dt, unit = units_out%time, namespace=namespace)
212 if (parse_is_defined(namespace, 'TDMaxSteps') .and. parse_is_defined(namespace, 'TDPropagationTime')) then
213 call messages_write('You cannot set TDMaxSteps and TDPropagationTime at the same time')
214 call messages_fatal(namespace=namespace)
215 end if
216
217 !%Variable TDPropagationTime
218 !%Type float
219 !%Section Time-Dependent::Propagation
220 !%Description
221 !% The length of the time propagation. You cannot set this variable
222 !% at the same time as <tt>TDMaxSteps</tt>. By default this variable will
223 !% not be used.
224 !%
225 !% The units for this variable are <math>\hbar</math>/Hartree (or <math>\hbar</math>/eV if you
226 !% selected <tt>ev_angstrom</tt> as input units). The approximate conversions to
227 !% femtoseconds are 1 fs = 41.34 <math>\hbar</math>/Hartree = 1.52 <math>\hbar</math>/eV.
228 !%End
229 call parse_variable(namespace, 'TDPropagationTime', -1.0_real64, propagation_time, unit = units_inp%time)
230
231 call messages_obsolete_variable(namespace, 'TDMaximumIter', 'TDMaxSteps')
232
233 !%Variable TDMaxSteps
234 !%Type integer
235 !%Default 1500
236 !%Section Time-Dependent::Propagation
237 !%Description
238 !% Number of time-propagation steps that will be performed. You
239 !% cannot use this variable together with <tt>TDPropagationTime</tt>.
240 !%End
241 default = 1500
242 if (propagation_time > m_zero) default = nint(propagation_time/td%dt)
243 call parse_variable(namespace, 'TDMaxSteps', default, td%max_iter)
244
245 if (propagation_time <= m_zero) propagation_time = td%dt*td%max_iter
247 call messages_print_var_value('TDPropagationTime', propagation_time, unit = units_out%time, namespace=namespace)
248 call messages_print_var_value('TDMaxSteps', td%max_iter, namespace=namespace)
249
250 if (td%max_iter < 1) then
251 write(message(1), '(a,i6,a)') "Input: '", td%max_iter, "' is not a valid value for TDMaxSteps."
252 message(2) = '(TDMaxSteps <= 1)'
253 call messages_fatal(2, namespace=namespace)
254 end if
255
256 td%iter = 0
257
258 td%dt = td%dt
259
260 lasers => list_get_lasers(ext_partners)
261
262 ! now the photoelectron stuff
263 call pes_init(td%pesv, namespace, space, gr, gr%box, st, outp%restart_write_interval, hm%kpoints, &
264 hm%abs_boundaries, ext_partners, td%max_iter, td%dt)
265
266 !%Variable TDDynamics
267 !%Type integer
268 !%Default ehrenfest
269 !%Section Time-Dependent::Propagation
270 !%Description
271 !% Type of dynamics to follow during a time propagation.
272 !% For BO, you must set <tt>MoveIons = yes</tt>.
273 !%Option ehrenfest 1
274 !% Ehrenfest dynamics.
275 !%Option bo 2
276 !% Born-Oppenheimer (Experimental).
277 !%End
278
279 call parse_variable(namespace, 'TDDynamics', ehrenfest, td%dynamics)
280 if (.not. varinfo_valid_option('TDDynamics', td%dynamics)) call messages_input_error(namespace, 'TDDynamics')
281 call messages_print_var_option('TDDynamics', td%dynamics, namespace=namespace)
282 if (td%dynamics .ne. ehrenfest) then
283 if (.not. ion_dynamics_ions_move(td%ions_dyn)) call messages_input_error(namespace, 'TDDynamics')
284 end if
285
286 !%Variable RecalculateGSDuringEvolution
287 !%Type logical
288 !%Default no
289 !%Section Time-Dependent::Propagation
290 !%Description
291 !% In order to calculate some information about the system along the
292 !% evolution (e.g. projection onto the ground-state KS determinant,
293 !% projection of the TDKS spin-orbitals onto the ground-state KS
294 !% spin-orbitals), the ground-state KS orbitals are needed. If the
295 !% ionic potential changes -- that is, the ions move -- one may want
296 !% to recalculate the ground state. You may do this by setting this
297 !% variable.
298 !%
299 !% The recalculation is not done every time step, but only every
300 !% <tt>RestartWriteInterval</tt> time steps.
301 !%End
302 call parse_variable(namespace, 'RecalculateGSDuringEvolution', .false., td%recalculate_gs)
303 if (hm%lda_u_level /= dft_u_none .and. td%recalculate_gs) then
304 call messages_not_implemented("DFT+U with RecalculateGSDuringEvolution=yes", namespace=namespace)
305 end if
306
307 !%Variable TDScissor
308 !%Type float
309 !%Default 0.0
310 !%Section Time-Dependent
311 !%Description
312 !% (experimental) If set, a scissor operator will be applied in the
313 !% Hamiltonian, shifting the excitation energies by the amount
314 !% specified. By default, it is not applied.
315 !%End
316 call parse_variable(namespace, 'TDScissor', m_zero, td%scissor)
317 td%scissor = units_to_atomic(units_inp%energy, td%scissor)
318 call messages_print_var_value('TDScissor', td%scissor, namespace=namespace)
319
320 call propagator_elec_init(gr, namespace, st, td%tr, hm%ks_pot, ion_dynamics_ions_move(td%ions_dyn) .or. &
321 list_has_gauge_field(ext_partners), family_is_mgga_with_exc(ks%xc))
322
323 if (associated(lasers) .and. mpi_grp_is_root(mpi_world)) then
324 call messages_print_with_emphasis(msg="Time-dependent external fields", namespace=namespace)
325 call laser_write_info(lasers%lasers, dt=td%dt, max_iter=td%max_iter, namespace=namespace)
326 call messages_print_with_emphasis(namespace=namespace)
327 end if
328
329 !%Variable TDEnergyUpdateIter
330 !%Type integer
331 !%Section Time-Dependent::Propagation
332 !%Description
333 !% This variable controls after how many iterations Octopus
334 !% updates the total energy during a time-propagation run. For
335 !% iterations where the energy is not updated, the last calculated
336 !% value is reported. If you set this variable to 1, the energy
337 !% will be calculated in each step.
338 !%End
339
340 default = 10
341 call parse_variable(namespace, 'TDEnergyUpdateIter', default, td%energy_update_iter)
342
343 if (gr%der%boundaries%spiralBC .and. hm%ep%reltype == spin_orbit) then
344 message(1) = "Generalized Bloch theorem cannot be used with spin-orbit coupling."
345 call messages_fatal(1, namespace=namespace)
346 end if
347
348 if (gr%der%boundaries%spiralBC) then
349 if (any(abs(hm%kick%easy_axis(1:2)) > m_epsilon)) then
350 message(1) = "Generalized Bloch theorem cannot be used for an easy axis not along the z direction."
351 call messages_fatal(1, namespace=namespace)
352 end if
353 end if
354
355 !%Variable TDFreezeOrbitals
356 !%Type integer
357 !%Default 0
358 !%Section Time-Dependent
359 !%Description
360 !% (Experimental) You have the possibility of "freezing" a number of orbitals during a time-propagation.
361 !% The Hartree and exchange-correlation potential due to these orbitals (which
362 !% will be the lowest-energy ones) will be added during the propagation, but the orbitals
363 !% will not be propagated.
364 !%Option sae -1
365 !% Single-active-electron approximation. This option is only valid for time-dependent
366 !% calculations (<tt>CalculationMode = td</tt>). Also, the nuclei should not move.
367 !% The idea is that all orbitals except the last one are frozen. The orbitals are to
368 !% be read from a previous ground-state calculation. The active orbital is then treated
369 !% as independent (whether it contains one electron or two) -- although it will
370 !% feel the Hartree and exchange-correlation potentials from the ground-state electronic
371 !% configuration.
372 !%
373 !% It is almost equivalent to setting <tt>TDFreezeOrbitals = N-1</tt>, where <tt>N</tt> is the number
374 !% of orbitals, but not completely.
375 !%End
376 call parse_variable(namespace, 'TDFreezeOrbitals', 0, td%freeze_orbitals)
377
378 if (td%freeze_orbitals /= 0) then
379 call messages_experimental('TDFreezeOrbitals', namespace=namespace)
380
381 if (hm%lda_u_level /= dft_u_none) then
382 call messages_not_implemented('TDFreezeOrbitals with DFT+U', namespace=namespace)
383 end if
384 end if
385
386
387 pop_sub(td_init)
388 nullify(lasers)
389 end subroutine td_init
390
391 ! ---------------------------------------------------------
392 subroutine td_init_run(td, namespace, mc, gr, ions, st, ks, hm, ext_partners, outp, space, from_scratch)
393 type(td_t), intent(inout) :: td
394 type(namespace_t), intent(in) :: namespace
395 type(multicomm_t), intent(inout) :: mc
396 type(grid_t), intent(inout) :: gr
397 type(ions_t), intent(inout) :: ions
398 type(states_elec_t), intent(inout) :: st
399 type(v_ks_t), intent(inout) :: ks
400 type(hamiltonian_elec_t), intent(inout) :: hm
401 type(partner_list_t), intent(in) :: ext_partners
402 type(output_t), intent(inout) :: outp
403 type(electron_space_t), intent(in) :: space
404 logical, intent(inout) :: from_scratch
405 push_sub(td_init_run)
406
407 ! NOTE: please do not change code in this function, but only in functions
408 ! called from here because the logic of this function is replicated in the
409 ! multisystem framework in different places
410
411 call td_allocate_wavefunctions(td, namespace, mc, gr, ions, st, hm, space)
412 call td_init_gaugefield(td, namespace, gr, st, ks, hm, ext_partners, space)
413
414 td%from_scratch = from_scratch
415
416 if (.not. td%from_scratch) then
417 call td_load_restart_from_td(td, namespace, space, mc, gr, ext_partners, st, ks, hm, td%from_scratch)
418 if (td%from_scratch) then
419 message(1) = "Unable to read time-dependent restart information: Starting from scratch"
420 call messages_warning(1, namespace=namespace)
421 end if
422 end if
423
424 if (td%iter >= td%max_iter) then
425 message(1) = "All requested iterations have already been done. Use FromScratch = yes if you want to redo them."
426 call messages_info(1, namespace=namespace)
428 td%iter = td%iter + 1
429 if (ion_dynamics_ions_move(td%ions_dyn) .and. td%recalculate_gs) call restart_end(td%restart_load)
430 pop_sub(td_init_run)
431 return
432 end if
433
434 if (td%from_scratch) then
435 call td_load_restart_from_gs(td, namespace, space, mc, gr, ext_partners, st, ks, hm)
436 end if
437
438 call td_init_with_wavefunctions(td, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, outp, td%from_scratch)
439
440 pop_sub(td_init_run)
441 end subroutine td_init_run
442
443 ! ---------------------------------------------------------
444 subroutine td_allocate_wavefunctions(td, namespace, mc, gr, ions, st, hm, space)
445 type(td_t), intent(inout) :: td
446 type(namespace_t), intent(in) :: namespace
447 type(multicomm_t), intent(inout) :: mc
448 type(grid_t), intent(inout) :: gr
449 type(ions_t), intent(inout) :: ions
450 type(states_elec_t), intent(inout) :: st
451 type(hamiltonian_elec_t), intent(inout) :: hm
452 class(space_t), intent(in) :: space
453
455
456 ! Allocate wavefunctions during time-propagation
457 if (td%dynamics == ehrenfest) then
458 !Note: this is not really clean to do this
459 if (hm%lda_u_level /= dft_u_none .and. states_are_real(st)) then
460 call lda_u_end(hm%lda_u)
461 !complex wfs are required for Ehrenfest
462 call states_elec_allocate_wfns(st, gr, type_cmplx, packed=.true.)
463 call lda_u_init(hm%lda_u, namespace, space, hm%lda_u_level, gr, ions, st, mc, &
464 hm%kpoints, hm%phase%is_allocated())
465 else
466 !complex wfs are required for Ehrenfest
467 call states_elec_allocate_wfns(st, gr, type_cmplx, packed=.true.)
468 end if
469 else
470 call states_elec_allocate_wfns(st, gr, packed=.true.)
471 call scf_init(td%scf, namespace, gr, ions, st, mc, hm, space)
472 end if
473
475 end subroutine td_allocate_wavefunctions
476
477 ! ---------------------------------------------------------
478 subroutine td_init_gaugefield(td, namespace, gr, st, ks, hm, ext_partners, space)
479 type(td_t), intent(inout) :: td
480 type(namespace_t), intent(in) :: namespace
481 type(grid_t), intent(inout) :: gr
482 type(states_elec_t), intent(inout) :: st
483 type(v_ks_t), intent(inout) :: ks
484 type(hamiltonian_elec_t), intent(inout) :: hm
485 type(partner_list_t), intent(in) :: ext_partners
486 class(space_t), intent(in) :: space
487
488 type(gauge_field_t), pointer :: gfield
489
490 push_sub(td_init_gaugefield)
491
492 gfield => list_get_gauge_field(ext_partners)
493 if(associated(gfield)) then
494 if (gauge_field_is_used(gfield)) then
495 !if the gauge field is applied, we need to tell v_ks to calculate the current
496 call v_ks_calculate_current(ks, .true.)
497
498 ! initialize the vector field and update the hamiltonian
499 call gauge_field_init_vec_pot(gfield, st%qtot)
500 call hm%update(gr, namespace, space, ext_partners, time = td%dt*td%iter)
501
502 end if
503 end if
504
505 pop_sub(td_init_gaugefield)
506 end subroutine td_init_gaugefield
507
508 ! ---------------------------------------------------------
509 subroutine td_end(td)
510 type(td_t), intent(inout) :: td
511
512 push_sub(td_end)
513
514 call pes_end(td%pesv)
515 call propagator_elec_end(td%tr) ! clean the evolution method
516 call ion_dynamics_end(td%ions_dyn)
517
518 if (td%dynamics == bo) call scf_end(td%scf)
519
520 pop_sub(td_end)
521 end subroutine td_end
522
523 ! ---------------------------------------------------------
524 subroutine td_end_run(td, st, hm)
525 type(td_t), intent(inout) :: td
526 type(states_elec_t), intent(inout) :: st
527 type(hamiltonian_elec_t), intent(inout) :: hm
528
529 push_sub(td_end_run)
530
531 if (st%pack_states .and. hm%apply_packed()) call st%unpack()
532
533 call restart_end(td%restart_dump)
534 call td_write_end(td%write_handler)
535
536 ! free memory
538 if (ion_dynamics_ions_move(td%ions_dyn) .and. td%recalculate_gs) call restart_end(td%restart_load)
539
540 pop_sub(td_end_run)
541 end subroutine td_end_run
542
543 ! ---------------------------------------------------------
544 subroutine td_run(td, namespace, mc, gr, ions, st, ks, hm, ext_partners, outp, space, from_scratch)
545 type(td_t), intent(inout) :: td
546 type(namespace_t), intent(in) :: namespace
547 type(multicomm_t), intent(inout) :: mc
548 type(grid_t), intent(inout) :: gr
549 type(ions_t), intent(inout) :: ions
550 type(states_elec_t), intent(inout) :: st
551 type(v_ks_t), intent(inout) :: ks
552 type(hamiltonian_elec_t), intent(inout) :: hm
553 type(partner_list_t), intent(in) :: ext_partners
554 type(output_t), intent(inout) :: outp
555 type(electron_space_t), intent(in) :: space
556 logical, intent(inout) :: from_scratch
557
558 logical :: stopping
559 integer :: iter, scsteps
560 real(real64) :: etime
561
562 push_sub(td_run)
563
564 etime = loct_clock()
565 ! This is the time-propagation loop. It starts at t=0 and finishes
566 ! at td%max_iter*dt. The index i runs from 1 to td%max_iter, and
567 ! step "iter" means propagation from (iter-1)*dt to iter*dt.
568 propagation: do iter = td%iter, td%max_iter
569
570 stopping = clean_stop(mc%master_comm) .or. walltimer_alarm(mc%master_comm)
572 call profiling_in("TIME_STEP")
573
574 if (iter > 1) then
575 if (((iter-1)*td%dt <= hm%kick%time) .and. (iter*td%dt > hm%kick%time)) then
576 if (.not. hm%pcm%localf) then
577 call kick_apply(space, gr, st, td%ions_dyn, ions, hm%kick, hm%psolver, hm%kpoints)
578 else
579 call kick_apply(space, gr, st, td%ions_dyn, ions, hm%kick, hm%psolver, hm%kpoints, pcm = hm%pcm)
580 end if
581 call td_write_kick(outp, namespace, space, gr, hm%kick, ions, iter)
582 !We activate the sprial BC only after the kick,
583 !to be sure that the first iteration corresponds to the ground state
584 if (gr%der%boundaries%spiralBC) gr%der%boundaries%spiral = .true.
585 end if
586 end if
587
588 ! time iterate the system, one time step.
589 select case (td%dynamics)
590 case (ehrenfest)
591 call propagator_elec_dt(ks, namespace, space, hm, gr, st, td%tr, iter*td%dt, td%dt, iter, td%ions_dyn, &
592 ions, ext_partners, outp, td%write_handler, scsteps = scsteps, &
593 update_energy = (mod(iter, td%energy_update_iter) == 0) .or. (iter == td%max_iter))
594 case (bo)
595 call propagator_elec_dt_bo(td%scf, namespace, space, gr, ks, st, hm, ions, ext_partners, mc, outp, iter, &
596 td%dt, td%ions_dyn, scsteps)
597 end select
598
599 !Apply mask absorbing boundaries
600 if (hm%abs_boundaries%abtype == mask_absorbing) then
601 if (states_are_real(st)) then
602 call dvmask(gr, hm, st)
603 else
604 call zvmask(gr, hm, st)
605 end if
606 end if
607
608 !Photoelectron stuff
609 if (td%pesv%calc_spm .or. td%pesv%calc_mask .or. td%pesv%calc_flux) then
610 call pes_calc(td%pesv, namespace, space, gr, st, td%dt, iter, gr%der, hm%kpoints, ext_partners, stopping)
611 end if
612
613 call td_write_iter(td%write_handler, namespace, space, outp, gr, st, hm, ions, ext_partners, &
614 hm%kick, ks, td%dt, iter, mc, td%recalculate_gs)
615
616 ! write down data
617 call td_check_point(td, namespace, mc, gr, ions, st, ks, hm, ext_partners, outp, space, &
618 iter, scsteps, etime, stopping, from_scratch)
619
620 ! check if debug mode should be enabled or disabled on the fly
621 call io_debug_on_the_fly(namespace)
622
623 call profiling_out("TIME_STEP")
624 if (stopping) exit
625
626 end do propagation
627
628 pop_sub(td_run)
629 end subroutine td_run
630
631 subroutine td_print_header(namespace)
632 type(namespace_t), intent(in) :: namespace
633
634 push_sub(td_print_header)
635
636 write(message(1), '(a7,1x,a14,a14,a10,a17)') 'Iter ', 'Time ', 'Energy ', 'SC Steps', 'Elapsed Time '
638 call messages_info(1, namespace=namespace)
639 call messages_print_with_emphasis(namespace=namespace)
640
641 pop_sub(td_print_header)
642 end subroutine td_print_header
643
644 ! ---------------------------------------------------------
645 subroutine td_check_point(td, namespace, mc, gr, ions, st, ks, hm, ext_partners, outp, space, &
646 iter, scsteps, etime, stopping, from_scratch)
647 type(td_t), intent(inout) :: td
648 type(namespace_t), intent(in) :: namespace
649 type(multicomm_t), intent(in) :: mc
650 type(grid_t), intent(inout) :: gr
651 type(ions_t), intent(inout) :: ions
652 type(states_elec_t), intent(inout) :: st
653 type(v_ks_t), intent(inout) :: ks
654 type(hamiltonian_elec_t), intent(inout) :: hm
655 type(partner_list_t), intent(in) :: ext_partners
656 type(output_t), intent(in) :: outp
657 type(electron_space_t), intent(in) :: space
658 integer, intent(in) :: iter
659 integer, intent(in) :: scsteps
660 real(real64), intent(inout) :: etime
661 logical, intent(in) :: stopping
662 logical, intent(inout) :: from_scratch
663
664 integer :: ierr
665
666 push_sub(td_check_point)
667
668 call td_print_message(td, namespace, ions, hm, iter, scsteps, etime)
669
670 if (outp%anything_now(iter)) then ! output
671 call td_write_output(namespace, space, gr, st, hm, ks, outp, ions, ext_partners, iter, td%dt)
672 end if
673
674 if (mod(iter, outp%restart_write_interval) == 0 .or. iter == td%max_iter .or. stopping) then ! restart
675 !if (iter == td%max_iter) outp%iter = ii - 1
676 call td_write_data(td%write_handler)
677 call td_dump(td, namespace, space, gr, st, hm, ks, ext_partners, iter, ierr)
678 if (ierr /= 0) then
679 message(1) = "Unable to write time-dependent restart information."
680 call messages_warning(1, namespace=namespace)
681 end if
682
683 call pes_output(td%pesv, namespace, space, gr, st, iter, outp, td%dt, ions)
684
685 if (ion_dynamics_ions_move(td%ions_dyn) .and. td%recalculate_gs) then
686 call messages_print_with_emphasis(msg='Recalculating the ground state.', namespace=namespace)
687 from_scratch = .false.
689 call electrons_ground_state_run(namespace, mc, gr, ions, ext_partners, st, ks, hm, outp, space, from_scratch)
690 call states_elec_allocate_wfns(st, gr, packed=.true.)
691 call td_load(td%restart_load, namespace, space, gr, st, hm, ext_partners, td, ks, ierr)
692 if (ierr /= 0) then
693 message(1) = "Unable to load TD states."
694 call messages_fatal(1, namespace=namespace)
695 end if
696 call density_calc(st, gr, st%rho)
697 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, &
698 calc_eigenval=.true., time = iter*td%dt, calc_energy=.true.)
699 call forces_calculate(gr, namespace, ions, hm, ext_partners, st, ks, t = iter*td%dt, dt = td%dt)
700 call messages_print_with_emphasis(msg="Time-dependent simulation proceeds", namespace=namespace)
701 call td_print_header(namespace)
702 end if
703 end if
704
705 pop_sub(td_check_point)
706 end subroutine td_check_point
707
708 ! ---------------------------------------------------------
709 subroutine td_print_message(td, namespace, ions, hm, iter, scsteps, etime)
710 type(td_t), intent(inout) :: td
711 type(namespace_t), intent(in) :: namespace
712 type(ions_t), intent(inout) :: ions
713 type(hamiltonian_elec_t), intent(inout) :: hm
714 integer, intent(in) :: iter
715 integer, intent(in) :: scsteps
716 real(real64), intent(inout) :: etime
717
718 push_sub(td_print_message)
719
720 write(message(1), '(i7,1x,2f14.6,i10,f14.3)') iter, units_from_atomic(units_out%time, iter*td%dt), &
721 units_from_atomic(units_out%energy, hm%energy%total + ions%kinetic_energy), scsteps, &
722 loct_clock() - etime
723 call messages_info(1, namespace=namespace)
725
726 pop_sub(td_print_message)
727 end subroutine td_print_message
728
729 ! ---------------------------------------------------------
730 subroutine td_update_elapsed_time(etime)
731 real(real64), intent(inout) :: etime
732
733 push_sub(td_update_elapsed_time)
734
735 etime = loct_clock()
736
739
740 ! ---------------------------------------------------------
741 subroutine td_init_with_wavefunctions(td, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, outp, from_scratch)
742 type(td_t), intent(inout) :: td
743 type(namespace_t), intent(in) :: namespace
744 type(electron_space_t), intent(in) :: space
745 type(multicomm_t), intent(in) :: mc
746 type(grid_t), intent(inout) :: gr
747 type(ions_t), intent(inout) :: ions
748 type(partner_list_t), intent(in) :: ext_partners
749 type(states_elec_t), target, intent(inout) :: st
750 type(v_ks_t), intent(inout) :: ks
751 type(hamiltonian_elec_t), intent(inout) :: hm
752 type(output_t), intent(inout) :: outp
753 logical, intent(in) :: from_scratch
754
755 integer :: ierr
756 real(real64) :: x
757 real(real64) :: ndinitial(space%dim)
758 logical :: freeze_hxc, freeze_occ, freeze_u
759 type(restart_t) :: restart, restart_frozen
760 type(gauge_field_t), pointer :: gfield
761 type(lasers_t), pointer :: lasers
763
764 !We activate the sprial BC only after the kick,
765 !to be sure that the first iteration corresponds to the ground state
766 if (gr%der%boundaries%spiralBC) then
767 if ((td%iter-1)*td%dt > hm%kick%time) then
768 gr%der%boundaries%spiral = .true.
769 end if
770 hm%vnl%spin => st%spin
771 hm%phase%spin => st%spin
772 !We fill st%spin. In case of restart, we read it in td_load
773 if (from_scratch) call states_elec_fermi(st, namespace, gr)
774 end if
775
776 if (from_scratch) then
777 ! Initialize the occupation matrices and U for DFT+U
778 ! This must be called before parsing TDFreezeOccupations and TDFreezeU
779 ! in order that the code does properly the initialization.
780 call lda_u_update_occ_matrices(hm%lda_u, namespace, gr, st, hm%hm_base, hm%phase, hm%energy)
781 end if
782
783 if (td%freeze_orbitals > 0) then
784 if (from_scratch) then
785 ! In this case, we first freeze the orbitals, then calculate the Hxc potential.
786 call states_elec_freeze_orbitals(st, namespace, space, gr, mc, hm%kpoints, &
787 td%freeze_orbitals, family_is_mgga(ks%xc_family))
788 else
789 call restart_init(restart, namespace, restart_td, restart_type_load, mc, ierr, mesh=gr)
790 if (ierr == 0) then
791 call td_load_frozen(namespace, restart, space, gr, st, hm, ierr)
792 end if
793 if (ierr /= 0) then
794 td%iter = 0
795 message(1) = "Unable to read frozen restart information."
796 call messages_fatal(1, namespace=namespace)
797 end if
798 call restart_end(restart)
799 end if
800 write(message(1),'(a,i4,a,i4,a)') 'Info: The lowest', td%freeze_orbitals, &
801 ' orbitals have been frozen.', st%nst, ' will be propagated.'
802 call messages_info(1, namespace=namespace)
804 call density_calc(st, gr, st%rho)
805 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_eigenval=.true., time = td%iter*td%dt)
806 else if (td%freeze_orbitals < 0) then
807 ! This means SAE approximation. We calculate the Hxc first, then freeze all
808 ! orbitals minus one.
809 write(message(1),'(a)') 'Info: The single-active-electron approximation will be used.'
810 call messages_info(1, namespace=namespace)
811 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_eigenval=.true., time = td%iter*td%dt)
812 if (from_scratch) then
813 call states_elec_freeze_orbitals(st, namespace, space, gr, mc, hm%kpoints, st%nst-1, family_is_mgga(ks%xc_family))
814 else
815 call messages_not_implemented("TDFreezeOrbials < 0 with FromScratch=no", namespace=namespace)
816 end if
818 call v_ks_freeze_hxc(ks)
819 call density_calc(st, gr, st%rho)
820 else
821 ! Normal run.
822 call density_calc(st, gr, st%rho)
823 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_eigenval=.true., time = td%iter*td%dt)
824 end if
825
826 !%Variable TDFreezeHXC
827 !%Type logical
828 !%Default no
829 !%Section Time-Dependent
830 !%Description
831 !% The electrons are evolved as independent particles feeling the Hartree and
832 !% exchange-correlation potentials from the ground-state electronic configuration.
833 !%End
834 call parse_variable(namespace, 'TDFreezeHXC', .false., freeze_hxc)
835 if (freeze_hxc) then
836 write(message(1),'(a)') 'Info: Freezing Hartree and exchange-correlation potentials.'
837 call messages_info(1, namespace=namespace)
838
839 if (.not. from_scratch) then
840
841 call restart_init(restart_frozen, namespace, restart_gs, restart_type_load, mc, ierr, mesh=gr, exact=.true.)
842 call states_elec_load(restart_frozen, namespace, space, st, gr, hm%kpoints, ierr, label = ": gs")
843 call states_elec_transform(st, namespace, space, restart_frozen, gr, hm%kpoints)
844 call restart_end(restart_frozen)
845
846 call density_calc(st, gr, st%rho)
847 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_eigenval=.true., time = td%iter*td%dt)
848
849 call restart_init(restart_frozen, namespace, restart_td, restart_type_load, mc, ierr, mesh=gr)
850 call states_elec_load(restart_frozen, namespace, space, st, gr, hm%kpoints, ierr, iter=td%iter, label = ": td")
851 call restart_end(restart_frozen)
852 call hm%ks_pot%run_zero_iter(td%tr%vks_old)
853
854 end if
855
856 call v_ks_freeze_hxc(ks)
857
858 end if
859
860 x = minval(st%eigenval(st%st_start, :))
861 if (st%parallel_in_states) then
862 call st%mpi_grp%bcast(x, 1, mpi_double_precision, 0)
863 end if
864 call hm%update_span(gr%spacing(1:space%dim), x, namespace)
865 ! initialize Fermi energy
866 call states_elec_fermi(st, namespace, gr, compute_spin = .not. gr%der%boundaries%spiralBC)
867 call energy_calc_total(namespace, space, hm, gr, st, ext_partners)
868
869 !%Variable TDFreezeDFTUOccupations
870 !%Type logical
871 !%Default no
872 !%Section Time-Dependent
873 !%Description
874 !% The occupation matrices than enters in the DFT+U potential
875 !% are not evolved during the time evolution.
876 !%End
877 call parse_variable(namespace, 'TDFreezeDFTUOccupations', .false., freeze_occ)
878 if (freeze_occ) then
879 write(message(1),'(a)') 'Info: Freezing DFT+U occupation matrices that enters in the DFT+U potential.'
880 call messages_info(1, namespace=namespace)
881 call lda_u_freeze_occ(hm%lda_u)
882
883 !In this case we should reload GS wavefunctions
884 if (hm%lda_u_level /= dft_u_none .and. .not. from_scratch) then
885 call restart_init(restart_frozen, namespace, restart_gs, restart_type_load, mc, ierr, mesh=gr)
886 call lda_u_load(restart_frozen, hm%lda_u, st, hm%energy%dft_u, ierr, occ_only = .true.)
887 call restart_end(restart_frozen)
888 end if
889 end if
890
891 !%Variable TDFreezeU
892 !%Type logical
893 !%Default no
894 !%Section Time-Dependent
895 !%Description
896 !% The effective U of DFT+U is not evolved during the time evolution.
897 !%End
898 call parse_variable(namespace, 'TDFreezeU', .false., freeze_u)
899 if (freeze_u) then
900 write(message(1),'(a)') 'Info: Freezing the effective U of DFT+U.'
901 call messages_info(1, namespace=namespace)
902 call lda_u_freeze_u(hm%lda_u)
903
904 !In this case we should reload GS wavefunctions
905 if (hm%lda_u_level == dft_u_acbn0 .and. .not. from_scratch) then
906 call restart_init(restart_frozen, namespace, restart_gs, restart_type_load, mc, ierr, mesh=gr)
907 call lda_u_load(restart_frozen, hm%lda_u, st, hm%energy%dft_u, ierr, u_only = .true.)
908 call restart_end(restart_frozen)
909 write(message(1),'(a)') 'Loaded GS effective U of DFT+U'
910 call messages_info(1, namespace=namespace)
911 call lda_u_write_u(hm%lda_u, namespace=namespace)
912 call lda_u_write_v(hm%lda_u, namespace=namespace)
913 end if
914 end if
915
916 ! This needs to be called before the calculation of the forces,
917 ! as we need to test of we output the forces or not
918 call td_write_init(td%write_handler, namespace, space, outp, gr, st, hm, ions, ext_partners, &
919 ks, ion_dynamics_ions_move(td%ions_dyn), &
920 list_has_gauge_field(ext_partners), hm%kick, td%iter, td%max_iter, td%dt, mc)
921
922 ! Resets the nondipole integration after laser-file has been written.
923 lasers => list_get_lasers(ext_partners)
924 if(associated(lasers)) then
925 if (lasers_with_nondipole_field(lasers)) then
926 ndinitial(1:space%dim)=m_zero
927 call lasers_set_nondipole_parameters(lasers,ndinitial,m_zero)
928 end if
929 end if
930 nullify(lasers)
931
932 call td_init_ions_and_forces(td, namespace, space, gr, ions, ext_partners, st, ks, hm, outp)
933
934 if (td%scissor > m_epsilon) then
935 call scissor_init(hm%scissor, namespace, space, st, gr, hm%d, hm%kpoints, hm%phase, td%scissor, mc)
936 end if
937
938 if (td%iter == 0) call td_run_zero_iter(td, namespace, space, gr, ions, st, ks, hm, ext_partners, outp, mc)
939
940 gfield => list_get_gauge_field(ext_partners)
941 if(associated(gfield)) then
942 if (gauge_field_is_propagated(gfield)) then
943 if(ks%xc%kernel_lrc_alpha > m_epsilon) then
944 call gauge_field_get_force(gfield, gr, st%d%spin_channels, st%current, ks%xc%kernel_lrc_alpha)
945 call messages_experimental('TD-LRC kernel')
946 else
947 call gauge_field_get_force(gfield, gr, st%d%spin_channels, st%current)
948 endif
949 endif
950 end if
951
952 !call td_check_trotter(td, sys, h)
953 td%iter = td%iter + 1
954
955 call restart_init(td%restart_dump, namespace, restart_td, restart_type_dump, mc, ierr, mesh=gr)
956 if (ion_dynamics_ions_move(td%ions_dyn) .and. td%recalculate_gs) then
957 ! We will also use the TD restart directory as temporary storage during the time propagation
958 call restart_init(td%restart_load, namespace, restart_td, restart_type_load, mc, ierr, mesh=gr)
959 end if
960
961 call messages_print_with_emphasis(msg="Time-Dependent Simulation", namespace=namespace)
962 call td_print_header(namespace)
963
964 if (td%pesv%calc_spm .or. td%pesv%calc_mask .and. from_scratch) then
965 call pes_init_write(td%pesv,gr,st, namespace)
966 end if
967
968 if (st%pack_states .and. hm%apply_packed()) call st%pack()
969
971 end subroutine td_init_with_wavefunctions
972
973 ! ---------------------------------------------------------
974 subroutine td_init_ions_and_forces(td, namespace, space, gr, ions, ext_partners, st, ks, hm, outp)
975 type(td_t), intent(inout) :: td
976 type(namespace_t), intent(in) :: namespace
977 type(electron_space_t), intent(in) :: space
978 type(grid_t), intent(inout) :: gr
979 type(ions_t), intent(inout) :: ions
980 type(partner_list_t), intent(in) :: ext_partners
981 type(states_elec_t), target, intent(inout) :: st
982 type(v_ks_t), intent(inout) :: ks
983 type(hamiltonian_elec_t), intent(inout) :: hm
984 type(output_t), intent(inout) :: outp
985
987
988 ! Calculate initial forces and kinetic energy
989 if (ion_dynamics_ions_move(td%ions_dyn)) then
990 if (td%iter > 0) then
991 call td_read_coordinates(td, namespace, ions)
992 if (ion_dynamics_drive_ions(td%ions_dyn)) then
993 call ion_dynamics_propagate(td%ions_dyn, ions, td%iter*td%dt, td%dt, namespace)
994 end if
995 call hamiltonian_elec_epot_generate(hm, namespace, space, gr, ions, ext_partners, st, time = td%iter*td%dt)
996 ! recompute potential because the ions have moved
997 call v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_eigenval=.true., time = td%iter*td%dt)
998 end if
999
1000 call forces_calculate(gr, namespace, ions, hm, ext_partners, st, ks, t = td%iter*td%dt, dt = td%dt)
1001
1002 call ions%update_kinetic_energy()
1003 else
1004 if (outp%what(option__output__forces) .or. td%write_handler%out(out_separate_forces)%write) then
1005 call forces_calculate(gr, namespace, ions, hm, ext_partners, st, ks, t = td%iter*td%dt, dt = td%dt)
1006 end if
1007 end if
1008
1009 if (outp%what(option__output__stress)) then
1010 call stress_calculate(namespace, gr, hm, st, ions, ks, ext_partners)
1011 end if
1012
1014 end subroutine td_init_ions_and_forces
1015
1016 ! ---------------------------------------------------------
1017 subroutine td_load_restart_from_td(td, namespace, space, mc, gr, ext_partners, st, ks, hm, from_scratch)
1018 type(td_t), intent(inout) :: td
1019 type(namespace_t), intent(in) :: namespace
1020 class(space_t), intent(in) :: space
1021 type(multicomm_t), intent(in) :: mc
1022 type(grid_t), intent(inout) :: gr
1023 type(partner_list_t), intent(in) :: ext_partners
1024 type(states_elec_t), target, intent(inout) :: st
1025 type(v_ks_t), intent(inout) :: ks
1026 type(hamiltonian_elec_t), intent(inout) :: hm
1027 logical, intent(inout) :: from_scratch
1028
1029 integer :: ierr
1030 type(restart_t) :: restart
1031
1032 push_sub(td_load_restart_from_td)
1033
1034 !We redistribute the states before the restarting
1035 if (td%freeze_orbitals > 0) then
1036 call states_elec_freeze_redistribute_states(st, namespace, gr, mc, td%freeze_orbitals)
1037 end if
1038
1039 call restart_init(restart, namespace, restart_td, restart_type_load, mc, ierr, mesh=gr)
1040 if (ierr == 0) then
1041 call td_load(restart, namespace, space, gr, st, hm, ext_partners, td, ks, ierr)
1042 end if
1043 call restart_end(restart)
1044 if (ierr /= 0) then
1045 from_scratch = .true.
1046 td%iter = 0
1047 end if
1048
1050 end subroutine td_load_restart_from_td
1051
1052 ! ---------------------------------------------------------
1053 subroutine td_load_restart_from_gs(td, namespace, space, mc, gr, ext_partners, st, ks, hm)
1054 type(td_t), intent(inout) :: td
1055 type(namespace_t), intent(in) :: namespace
1056 class(space_t), intent(in) :: space
1057 type(multicomm_t), intent(in) :: mc
1058 type(grid_t), intent(inout) :: gr
1059 type(partner_list_t), intent(in) :: ext_partners
1060 type(states_elec_t), target, intent(inout) :: st
1061 type(v_ks_t), intent(inout) :: ks
1062 type(hamiltonian_elec_t), intent(inout) :: hm
1063
1064 integer :: ierr
1065 type(restart_t) :: restart
1066
1068
1069 call restart_init(restart, namespace, restart_gs, restart_type_load, mc, ierr, mesh=gr, exact=.true.)
1070
1071 if (.not. st%only_userdef_istates) then
1072 if (ierr == 0) then
1073 call states_elec_load(restart, namespace, space, st, gr, hm%kpoints, ierr, label = ": gs")
1074 end if
1075 if (ierr /= 0) then
1076 message(1) = 'Unable to read ground-state wavefunctions.'
1077 call messages_fatal(1, namespace=namespace)
1078 end if
1079 end if
1080
1081 ! check if we should deploy user-defined wavefunctions.
1082 ! according to the settings in the input file the routine
1083 ! overwrites orbitals that were read from restart/gs
1084 if (parse_is_defined(namespace, 'UserDefinedStates')) then
1085 call states_elec_read_user_def_orbitals(gr, namespace, space, st)
1086 end if
1087
1088 call states_elec_transform(st, namespace, space, restart, gr, hm%kpoints)
1089 call restart_end(restart)
1090
1092 end subroutine td_load_restart_from_gs
1093
1094 ! ---------------------------------------------------------
1095 subroutine td_run_zero_iter(td, namespace, space, gr, ions, st, ks, hm, ext_partners, outp, mc)
1096 type(td_t), intent(inout) :: td
1097 type(namespace_t), intent(in) :: namespace
1098 type(electron_space_t), intent(in) :: space
1099 type(grid_t), intent(inout) :: gr
1100 type(ions_t), intent(inout) :: ions
1101 type(states_elec_t), intent(inout) :: st
1102 type(v_ks_t), intent(inout) :: ks
1103 type(hamiltonian_elec_t), intent(inout) :: hm
1104 type(partner_list_t), intent(in) :: ext_partners
1105 type(output_t), intent(in) :: outp
1106 type(multicomm_t), intent(in) :: mc
1107
1108 push_sub(td_run_zero_iter)
1109
1110 call td_write_iter(td%write_handler, namespace, space, outp, gr, st, hm, ions, ext_partners, &
1111 hm%kick, ks, td%dt, 0, mc, td%recalculate_gs)
1112
1113 ! I apply the delta electric field *after* td_write_iter, otherwise the
1114 ! dipole matrix elements in write_proj are wrong
1115 if (abs(hm%kick%time) <= m_epsilon) then
1116 if (.not. hm%pcm%localf) then
1117 call kick_apply(space, gr, st, td%ions_dyn, ions, hm%kick, hm%psolver, hm%kpoints)
1118 else
1119 call kick_apply(space, gr, st, td%ions_dyn, ions, hm%kick, hm%psolver, hm%kpoints, pcm = hm%pcm)
1120 end if
1121 call td_write_kick(outp, namespace, space, gr, hm%kick, ions, 0)
1122
1123 !We activate the sprial BC only after the kick
1124 if (gr%der%boundaries%spiralBC) then
1125 gr%der%boundaries%spiral = .true.
1126 end if
1127 end if
1128 call hm%ks_pot%run_zero_iter(td%tr%vks_old)
1129
1130 if (any(outp%output_interval > 0)) then
1131 call td_write_data(td%write_handler)
1132 call td_write_output(namespace, space, gr, st, hm, ks, outp, ions, ext_partners, 0)
1133 end if
1134
1135 pop_sub(td_run_zero_iter)
1136 end subroutine td_run_zero_iter
1137
1138
1139 ! ---------------------------------------------------------
1141 subroutine td_read_coordinates(td, namespace, ions)
1142 type(td_t), intent(in) :: td
1143 type(namespace_t), intent(in) :: namespace
1144 type(ions_t), intent(inout) :: ions
1145
1146 integer :: iatom, iter, iunit
1147
1148 push_sub(td_read_coordinates)
1149
1150 iunit = io_open('td.general/coordinates', namespace, action='read', status='old', die=.false.)
1151 if (iunit == -1) then
1152 message(1) = "Could not open file '"//trim(io_workpath('td.general/coordinates', namespace))//"'."
1153 message(2) = "Starting simulation from initial geometry."
1154 call messages_warning(2, namespace=namespace)
1155 pop_sub(td_read_coordinates)
1156 return
1157 end if
1158
1159 call io_skip_header(iunit)
1160 do iter = 0, td%iter - 1
1161 read(iunit, *) ! skip previous iterations... sorry, but no portable seek in Fortran
1162 end do
1163 read(iunit, '(32x)', advance='no') ! skip the time index.
1164
1165 do iatom = 1, ions%natoms
1166 read(iunit, '(3es24.16)', advance='no') ions%pos(:, iatom)
1167 ions%pos(:, iatom) = units_to_atomic(units_out%length, ions%pos(:, iatom))
1168 end do
1169 do iatom = 1, ions%natoms
1170 read(iunit, '(3es24.16)', advance='no') ions%vel(:, iatom)
1171 ions%vel(:, iatom) = units_to_atomic(units_out%velocity, ions%vel(:, iatom))
1172 end do
1173 do iatom = 1, ions%natoms
1174 read(iunit, '(3es24.16)', advance='no') ions%tot_force(:, iatom)
1175 ions%tot_force(:, iatom) = units_to_atomic(units_out%force, ions%tot_force(:, iatom))
1176 end do
1177
1178 call io_close(iunit)
1179
1180 pop_sub(td_read_coordinates)
1181 end subroutine td_read_coordinates
1182
1183 ! ---------------------------------------------------------
1184 subroutine td_dump(td, namespace, space, gr, st, hm, ks, ext_partners, iter, ierr)
1185 type(td_t), intent(in) :: td
1186 type(namespace_t), intent(in) :: namespace
1187 class(space_t), intent(in) :: space
1188 type(grid_t), intent(in) :: gr
1189 type(states_elec_t), intent(in) :: st
1190 type(hamiltonian_elec_t), intent(in) :: hm
1191 type(v_ks_t), intent(in) :: ks
1192 type(partner_list_t), intent(in) :: ext_partners
1193 integer, intent(in) :: iter
1194 integer, intent(out) :: ierr
1195
1196 type(gauge_field_t), pointer :: gfield
1197 integer :: err, err2
1198
1199 push_sub(td_dump)
1200
1201 ierr = 0
1202
1203 if (restart_skip(td%restart_dump)) then
1204 pop_sub(td_dump)
1205 return
1206 end if
1207
1208 message(1) = "Debug: Writing td restart."
1209 call messages_info(1, namespace=namespace, debug_only=.true.)
1210
1211 ! first write resume file
1212 call states_elec_dump(td%restart_dump, space, st, gr, hm%kpoints, err, iter=iter)
1213 if (err /= 0) ierr = ierr + 1
1214
1215 call states_elec_dump_rho(td%restart_dump, space, st, gr, ierr, iter=iter)
1216 if (err /= 0) ierr = ierr + 1
1217
1218 if (hm%lda_u_level /= dft_u_none) then
1219 call lda_u_dump(td%restart_dump, namespace, hm%lda_u, st, gr, ierr)
1220 if (err /= 0) ierr = ierr + 1
1221 end if
1222
1223 call potential_interpolation_dump(td%tr%vks_old, space, td%restart_dump, gr, st%d%nspin, err2)
1224 if (err2 /= 0) ierr = ierr + 2
1225
1226 call pes_dump(td%pesv, namespace, td%restart_dump, st, gr, err)
1227 if (err /= 0) ierr = ierr + 4
1228
1229 ! Gauge field restart
1230 gfield => list_get_gauge_field(ext_partners)
1231 if(associated(gfield)) then
1232 call gauge_field_dump(td%restart_dump, gfield, ierr)
1233 end if
1235 if (gr%der%boundaries%spiralBC) then
1236 call states_elec_dump_spin(td%restart_dump, st, err)
1237 if (err /= 0) ierr = ierr + 8
1238 end if
1239
1240 if (ks%has_photons) then
1241 call mf_photons_dump(td%restart_dump, ks%pt_mx, gr, td%dt, ks%pt, err)
1242 if (err /= 0) ierr = ierr + 16
1243 end if
1244
1245 if (ks%xc_photon /= 0) then
1246 ! photon-free mean field
1247 call ks%xc_photons%mf_dump(td%restart_dump, err)
1248 if (err /= 0) ierr = ierr + 32
1249 end if
1250
1251 if (allocated(st%frozen_rho)) then
1252 call states_elec_dump_frozen(td%restart_dump, space, st, gr, ierr)
1253 end if
1254 if (err /= 0) ierr = ierr + 64
1255
1256 if (ion_dynamics_ions_move(td%ions_dyn)) then
1257 call ion_dynamics_dump(td%ions_dyn, td%restart_dump, err)
1258 end if
1259 if (err /= 0) ierr = ierr + 128
1260
1261 message(1) = "Debug: Writing td restart done."
1262 call messages_info(1, namespace=namespace, debug_only=.true.)
1263
1264 pop_sub(td_dump)
1265 end subroutine td_dump
1266
1267 ! ---------------------------------------------------------
1268 subroutine td_load(restart, namespace, space, gr, st, hm, ext_partners, td, ks, ierr)
1269 type(restart_t), intent(in) :: restart
1270 type(namespace_t), intent(in) :: namespace
1271 class(space_t), intent(in) :: space
1272 type(grid_t), intent(in) :: gr
1273 type(states_elec_t), intent(inout) :: st
1274 type(hamiltonian_elec_t), intent(inout) :: hm
1275 type(partner_list_t), intent(in) :: ext_partners
1276 type(td_t), intent(inout) :: td
1277 type(v_ks_t), intent(inout) :: ks
1278 integer, intent(out) :: ierr
1279
1280 integer :: err, err2
1281 type(gauge_field_t), pointer :: gfield
1282 push_sub(td_load)
1283
1284 ierr = 0
1285
1286 if (restart_skip(restart)) then
1287 ierr = -1
1288 pop_sub(td_load)
1289 return
1290 end if
1291
1292 message(1) = "Debug: Reading td restart."
1293 call messages_info(1, namespace=namespace, debug_only=.true.)
1294
1295 ! Read states
1296 call states_elec_load(restart, namespace, space, st, gr, hm%kpoints, err, iter=td%iter, label = ": td")
1297 if (err /= 0) then
1298 ierr = ierr + 1
1299 end if
1300
1301 ! read potential from previous interactions
1302 call potential_interpolation_load(td%tr%vks_old, namespace, space, restart, gr, st%d%nspin, err2)
1303 if (err2 /= 0) ierr = ierr + 2
1304
1305 if (hm%lda_u_level /= dft_u_none) then
1306 call lda_u_load(restart, hm%lda_u, st, hm%energy%dft_u, ierr)
1307 if (err /= 0) ierr = ierr + 1
1308 end if
1309
1310
1311 ! read PES restart
1312 if (td%pesv%calc_spm .or. td%pesv%calc_mask .or. td%pesv%calc_flux) then
1313 call pes_load(td%pesv, namespace, restart, st, err)
1314 if (err /= 0) ierr = ierr + 4
1315 end if
1316
1317 ! Gauge field restart
1318 gfield => list_get_gauge_field(ext_partners)
1319 if (associated(gfield)) then
1320 call gauge_field_load(restart, gfield, err)
1321 if (err /= 0) then
1322 ierr = ierr + 8
1323 else
1324 call hm%update(gr, namespace, space, ext_partners, time = td%dt*td%iter)
1325 end if
1326 end if
1327
1328 ! add photon restart
1329 if (ks%has_photons) then
1330 call mf_photons_load(restart, ks%pt_mx, gr, err)
1331 end if
1332 if (err /= 0) ierr = ierr + 16
1333
1334 if (ks%xc_photon /= 0) then
1335 call ks%xc_photons%mf_load(restart, space, err)
1336 end if
1337 if (err /= 0) ierr = ierr + 32
1338
1339 if (gr%der%boundaries%spiralBC) then
1340 call states_elec_load_spin(restart, st, err)
1341 !To ensure back compatibility, if the file is not present, we use the
1342 !current states to get the spins
1343 if (err /= 0) call states_elec_fermi(st, namespace, gr)
1344 end if
1345
1346 if (ion_dynamics_ions_move(td%ions_dyn)) then
1347 call ion_dynamics_load(td%ions_dyn, restart, err)
1348 end if
1349 if (err /= 0) ierr = ierr + 64
1350
1351 message(1) = "Debug: Reading td restart done."
1352 call messages_info(1, namespace=namespace, debug_only=.true.)
1353
1354 pop_sub(td_load)
1355 end subroutine td_load
1356
1357 ! ---------------------------------------------------------
1358 subroutine td_load_frozen(namespace, restart, space, mesh, st, hm, ierr)
1359 type(namespace_t), intent(in) :: namespace
1360 type(restart_t), intent(in) :: restart
1361 class(space_t), intent(in) :: space
1362 class(mesh_t), intent(in) :: mesh
1363 type(states_elec_t), intent(inout) :: st
1364 type(hamiltonian_elec_t), intent(inout) :: hm
1365 integer, intent(out) :: ierr
1366
1367 push_sub(td_load_frozen)
1368
1369 ierr = 0
1370
1371 if (restart_skip(restart)) then
1372 ierr = -1
1373 pop_sub(td_load_frozen)
1374 return
1375 end if
1376
1377 message(1) = "Debug: Reading td frozen restart."
1378 call messages_info(1, namespace=namespace, debug_only=.true.)
1379
1380 safe_allocate(st%frozen_rho(1:mesh%np, 1:st%d%nspin))
1381 if (family_is_mgga(hm%xc%family)) then
1382 safe_allocate(st%frozen_tau(1:mesh%np, 1:st%d%nspin))
1383 safe_allocate(st%frozen_gdens(1:mesh%np, 1:space%dim, 1:st%d%nspin))
1384 safe_allocate(st%frozen_ldens(1:mesh%np, 1:st%d%nspin))
1385 end if
1386
1387 call states_elec_load_frozen(restart, space, st, mesh, ierr)
1388
1389 message(1) = "Debug: Reading td frozen restart done."
1390 call messages_info(1, namespace=namespace, debug_only=.true.)
1391
1392 pop_sub(td_load_frozen)
1393 end subroutine td_load_frozen
1394
1395 ! ---------------------------------------------------------
1396 logical function td_get_from_scratch(td)
1397 type(td_t), intent(in) :: td
1398
1399 push_sub(td_get_from_scratch)
1400
1401 td_get_from_scratch = td%from_scratch
1402
1403 pop_sub(td_get_from_scratch)
1404 end function td_get_from_scratch
1405
1406 ! ---------------------------------------------------------
1407 subroutine td_set_from_scratch(td, from_scratch)
1408 type(td_t), intent(inout) :: td
1409 logical, intent(in) :: from_scratch
1410
1411 push_sub(td_set_from_scratch)
1412
1413 td%from_scratch = from_scratch
1414
1415 pop_sub(td_set_from_scratch)
1416 end subroutine td_set_from_scratch
1417end module td_oct_m
1418
1419!! Local Variables:
1420!! mode: f90
1421!! coding: utf-8
1422!! End:
Prints out to iunit a message in the form: ["InputVariable" = value] where "InputVariable" is given b...
Definition: messages.F90:180
integer, parameter, public mask_absorbing
Module implementing boundary conditions in Octopus.
Definition: boundaries.F90:122
This module handles the calculation mode.
type(calc_mode_par_t), public calc_mode_par
Singleton instance of parallel calculation mode.
integer, parameter, public p_strategy_states
parallelization in states
This module implements a calculator for the density and defines related functions.
Definition: density.F90:120
subroutine, public states_elec_freeze_adjust_qtot(st)
Definition: density.F90:811
subroutine, public states_elec_freeze_redistribute_states(st, namespace, mesh, mc, nn)
Definition: density.F90:783
subroutine, public states_elec_freeze_orbitals(st, namespace, space, gr, mc, kpoints, n, family_is_mgga)
Calculate partial density for frozen orbitals.
Definition: density.F90:652
subroutine, public density_calc(st, gr, density, istin)
Computes the density from the orbitals in st.
Definition: density.F90:609
subroutine, public electrons_ground_state_run(namespace, mc, gr, ions, ext_partners, st, ks, hm, outp, space, fromScratch)
subroutine, public energy_calc_total(namespace, space, hm, gr, st, ext_partners, iunit, full)
This subroutine calculates the total energy of the system. Basically, it adds up the KS eigenvalues,...
integer, parameter, public spin_orbit
Definition: epot.F90:166
logical function, public list_has_gauge_field(partners)
type(gauge_field_t) function, pointer, public list_get_gauge_field(partners)
type(lasers_t) function, pointer, public list_get_lasers(partners)
subroutine, public forces_calculate(gr, namespace, ions, hm, ext_partners, st, ks, vhxc_old, t, dt)
Definition: forces.F90:339
subroutine, public gauge_field_load(restart, gfield, ierr)
subroutine, public gauge_field_get_force(this, gr, spin_channels, current, lrc_alpha)
subroutine, public gauge_field_dump(restart, gfield, ierr)
logical pure function, public gauge_field_is_propagated(this)
logical pure function, public gauge_field_is_used(this)
subroutine, public gauge_field_init_vec_pot(this, qtot)
real(real64), parameter, public m_zero
Definition: global.F90:188
real(real64), parameter, public m_epsilon
Definition: global.F90:204
real(real64), parameter, public m_one
Definition: global.F90:189
This module implements the underlying real-space grid.
Definition: grid.F90:117
subroutine, public zvmask(mesh, hm, st)
subroutine, public hamiltonian_elec_epot_generate(this, namespace, space, gr, ions, ext_partners, st, time)
subroutine, public dvmask(mesh, hm, st)
This module defines classes and functions for interaction partners.
Definition: io.F90:114
subroutine, public io_close(iunit, grp)
Definition: io.F90:418
subroutine, public io_skip_header(iunit)
Definition: io.F90:597
character(len=max_path_len) function, public io_workpath(path, namespace)
Definition: io.F90:270
subroutine, public io_debug_on_the_fly(namespace)
check if debug mode should be enabled or disabled on the fly
Definition: io.F90:486
integer function, public io_open(file, namespace, action, status, form, position, die, recl, grp)
Definition: io.F90:352
subroutine, public ion_dynamics_dump(this, restart, ierr)
subroutine, public ion_dynamics_propagate(this, ions, time, dt, namespace)
logical pure function, public ion_dynamics_ions_move(this)
subroutine, public ion_dynamics_init(this, namespace, ions)
subroutine, public ion_dynamics_load(this, restart, ierr)
subroutine, public ion_dynamics_end(this)
logical pure function, public ion_dynamics_drive_ions(this)
subroutine, public kick_apply(space, mesh, st, ions_dyn, ions, kick, psolver, kpoints, pcm)
Applies the delta-function electric field where k = kick%delta_strength.
Definition: kick.F90:1136
subroutine, public lasers_set_nondipole_parameters(this, ndfield, nd_integration_time)
Set parameters for nondipole SFA calculation.
Definition: lasers.F90:752
logical function, public lasers_with_nondipole_field(lasers)
Check if a nondipole SFA correction should be computed for the given laser.
Definition: lasers.F90:739
subroutine, public laser_write_info(lasers, namespace, dt, max_iter, iunit)
Definition: lasers.F90:930
subroutine, public lda_u_dump(restart, namespace, this, st, mesh, ierr)
Definition: lda_u_io.F90:646
subroutine, public lda_u_write_u(this, iunit, namespace)
Definition: lda_u_io.F90:530
subroutine, public lda_u_load(restart, this, st, dftu_energy, ierr, occ_only, u_only)
Definition: lda_u_io.F90:727
subroutine, public lda_u_write_v(this, iunit, namespace)
Definition: lda_u_io.F90:579
integer, parameter, public dft_u_none
Definition: lda_u.F90:201
subroutine, public lda_u_freeze_occ(this)
Definition: lda_u.F90:887
subroutine, public lda_u_freeze_u(this)
Definition: lda_u.F90:894
subroutine, public lda_u_end(this)
Definition: lda_u.F90:654
subroutine, public lda_u_init(this, namespace, space, level, gr, ions, st, mc, kpoints, has_phase)
Definition: lda_u.F90:281
integer, parameter, public dft_u_acbn0
Definition: lda_u.F90:201
subroutine, public lda_u_update_occ_matrices(this, namespace, mesh, st, hm_base, phase, energy)
Definition: lda_u.F90:797
This module implements fully polymorphic linked lists, and some specializations thereof.
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:118
subroutine, public messages_print_with_emphasis(msg, iunit, namespace)
Definition: messages.F90:920
subroutine, public messages_not_implemented(feature, namespace)
Definition: messages.F90:1113
character(len=512), private msg
Definition: messages.F90:165
subroutine, public messages_warning(no_lines, all_nodes, namespace)
Definition: messages.F90:537
subroutine, public messages_obsolete_variable(namespace, name, rep)
Definition: messages.F90:1045
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
Definition: messages.F90:160
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
Definition: messages.F90:414
subroutine, public messages_input_error(namespace, var, details, row, column)
Definition: messages.F90:713
subroutine, public messages_experimental(name, namespace)
Definition: messages.F90:1085
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
Definition: messages.F90:616
logical function mpi_grp_is_root(grp)
Is the current MPI process of grpcomm, root.
Definition: mpi.F90:430
type(mpi_grp_t), public mpi_world
Definition: mpi.F90:266
This module handles the communicators for the various parallelization strategies.
Definition: multicomm.F90:145
this module contains the low-level part of the output system
Definition: output_low.F90:115
this module contains the output system
Definition: output.F90:115
logical function, public parse_is_defined(namespace, name)
Definition: parser.F90:502
subroutine, public pes_calc(pes, namespace, space, mesh, st, dt, iter, der, kpoints, ext_partners, stopping)
Definition: pes.F90:269
subroutine, public pes_init(pes, namespace, space, mesh, box, st, save_iter, kpoints, abs_boundaries, ext_partners, max_iter, dt)
Definition: pes.F90:175
subroutine, public pes_output(pes, namespace, space, gr, st, iter, outp, dt, ions)
Definition: pes.F90:295
subroutine, public pes_init_write(pes, mesh, st, namespace)
Definition: pes.F90:399
subroutine, public pes_end(pes)
Definition: pes.F90:255
subroutine, public pes_load(pes, namespace, restart, st, ierr)
Definition: pes.F90:359
subroutine, public pes_dump(pes, namespace, restart, st, mesh, ierr)
Definition: pes.F90:319
subroutine, public mf_photons_load(restart, this, gr, ierr)
subroutine, public mf_photons_dump(restart, this, gr, dt, pt_mode, ierr)
subroutine, public potential_interpolation_load(potential_interpolation, namespace, space, restart, mesh, nspin, err2)
subroutine, public potential_interpolation_dump(potential_interpolation, space, restart, mesh, nspin, err2)
subroutine, public profiling_out(label)
Increment out counter and sum up difference between entry and exit time.
Definition: profiling.F90:623
subroutine, public profiling_in(label, exclude)
Increment in counter and save entry time.
Definition: profiling.F90:552
subroutine, public propagator_elec_dt_bo(scf, namespace, space, gr, ks, st, hm, ions, ext_partners, mc, outp, iter, dt, ions_dyn, scsteps)
subroutine, public propagator_elec_init(gr, namespace, st, tr, ks_pot, have_fields, family_is_mgga_with_exc)
subroutine, public propagator_elec_dt(ks, namespace, space, hm, gr, st, tr, time, dt, nt, ions_dyn, ions, ext_partners, outp, write_handler, scsteps, update_energy, qcchi)
Propagates st from time - dt to t. If dt<0, it propagates backwards from t+|dt| to t.
subroutine, public propagator_elec_end(tr)
This module implements the basic propagator framework.
Definition: propagator.F90:117
logical function, public clean_stop(comm)
returns true if a file named stop exists
Definition: restart.F90:276
integer, parameter, public restart_gs
Definition: restart.F90:200
subroutine, public restart_init(restart, namespace, data_type, type, mc, ierr, mesh, dir, exact)
Initializes a restart object.
Definition: restart.F90:516
integer, parameter, public restart_type_dump
Definition: restart.F90:245
logical pure function, public restart_skip(restart)
Returns true if the restart information should neither be read nor written. This might happen because...
Definition: restart.F90:969
integer, parameter, public restart_td
Definition: restart.F90:200
integer, parameter, public restart_type_load
Definition: restart.F90:245
subroutine, public restart_end(restart)
Definition: restart.F90:722
subroutine, public scf_init(scf, namespace, gr, ions, st, mc, hm, space)
Definition: scf.F90:254
subroutine, public scf_end(scf)
Definition: scf.F90:543
subroutine, public scissor_init(this, namespace, space, st, mesh, d, kpoints, phase, gap, mc)
Definition: scissor.F90:162
pure logical function, public states_are_real(st)
subroutine, public states_elec_fermi(st, namespace, mesh, compute_spin)
calculate the Fermi level for the states in this object
subroutine, public states_elec_deallocate_wfns(st)
Deallocates the KS wavefunctions defined within a states_elec_t structure.
subroutine, public states_elec_allocate_wfns(st, mesh, wfs_type, skip, packed)
Allocates the KS wavefunctions defined within a states_elec_t structure.
This module handles reading and writing restart information for the states_elec_t.
subroutine, public states_elec_read_user_def_orbitals(mesh, namespace, space, st)
the routine reads formulas for user-defined wavefunctions from the input file and fills the respectiv...
subroutine, public states_elec_load_frozen(restart, space, st, mesh, ierr)
subroutine, public states_elec_dump(restart, space, st, mesh, kpoints, ierr, iter, lr, st_start_writing, verbose)
subroutine, public states_elec_transform(st, namespace, space, restart, mesh, kpoints, prefix)
subroutine, public states_elec_load(restart, namespace, space, st, mesh, kpoints, ierr, iter, lr, lowest_missing, label, verbose, skip)
returns in ierr: <0 => Fatal error, or nothing read =0 => read all wavefunctions >0 => could only rea...
subroutine, public states_elec_load_spin(restart, st, ierr)
returns in ierr: <0 => Fatal error, or nothing read =0 => read all wavefunctions >0 => could only rea...
subroutine, public states_elec_dump_frozen(restart, space, st, mesh, ierr)
subroutine, public states_elec_dump_rho(restart, space, st, mesh, ierr, iter)
subroutine, public states_elec_dump_spin(restart, st, ierr)
This module implements the calculation of the stress tensor.
Definition: stress.F90:118
subroutine, public stress_calculate(namespace, gr, hm, st, ions, ks, ext_partners)
This computes the total stress on the lattice.
Definition: stress.F90:186
Definition: td.F90:114
subroutine, public td_check_point(td, namespace, mc, gr, ions, st, ks, hm, ext_partners, outp, space, iter, scsteps, etime, stopping, from_scratch)
Definition: td.F90:740
subroutine, public td_end(td)
Definition: td.F90:603
subroutine, public td_load_restart_from_gs(td, namespace, space, mc, gr, ext_partners, st, ks, hm)
Definition: td.F90:1147
subroutine td_print_message(td, namespace, ions, hm, iter, scsteps, etime)
Definition: td.F90:803
subroutine, public td_run_init()
Definition: td.F90:236
subroutine, public td_end_run(td, st, hm)
Definition: td.F90:618
subroutine, public td_init(td, namespace, space, gr, ions, st, ks, hm, ext_partners, outp)
Definition: td.F90:247
subroutine, public td_allocate_wavefunctions(td, namespace, mc, gr, ions, st, hm, space)
Definition: td.F90:538
logical function, public td_get_from_scratch(td)
Definition: td.F90:1490
subroutine td_load_frozen(namespace, restart, space, mesh, st, hm, ierr)
Definition: td.F90:1452
subroutine td_print_header(namespace)
Definition: td.F90:725
integer, parameter, public bo
Definition: td.F90:201
subroutine td_load(restart, namespace, space, gr, st, hm, ext_partners, td, ks, ierr)
Definition: td.F90:1362
subroutine, public td_set_from_scratch(td, from_scratch)
Definition: td.F90:1501
subroutine, public td_dump(td, namespace, space, gr, st, hm, ks, ext_partners, iter, ierr)
Definition: td.F90:1278
subroutine td_read_coordinates(td, namespace, ions)
reads the pos and vel from coordinates file
Definition: td.F90:1235
subroutine, public td_run(td, namespace, mc, gr, ions, st, ks, hm, ext_partners, outp, space, from_scratch)
Definition: td.F90:638
subroutine td_run_zero_iter(td, namespace, space, gr, ions, st, ks, hm, ext_partners, outp, mc)
Definition: td.F90:1189
subroutine, public td_init_gaugefield(td, namespace, gr, st, ks, hm, ext_partners, space)
Definition: td.F90:572
subroutine td_update_elapsed_time(etime)
Definition: td.F90:824
subroutine, public td_init_with_wavefunctions(td, namespace, space, mc, gr, ions, ext_partners, st, ks, hm, outp, from_scratch)
Definition: td.F90:835
subroutine, public td_load_restart_from_td(td, namespace, space, mc, gr, ext_partners, st, ks, hm, from_scratch)
Definition: td.F90:1111
subroutine td_init_ions_and_forces(td, namespace, space, gr, ions, ext_partners, st, ks, hm, outp)
Definition: td.F90:1068
subroutine, public td_init_run(td, namespace, mc, gr, ions, st, ks, hm, ext_partners, outp, space, from_scratch)
Definition: td.F90:486
subroutine, public td_write_output(namespace, space, gr, st, hm, ks, outp, ions, ext_partners, iter, dt)
Definition: td_write.F90:1224
subroutine, public td_write_data(writ)
Definition: td_write.F90:1190
subroutine, public td_write_kick(outp, namespace, space, mesh, kick, ions, iter)
Definition: td_write.F90:333
subroutine, public td_write_end(writ)
Definition: td_write.F90:978
subroutine, public td_write_iter(writ, namespace, space, outp, gr, st, hm, ions, ext_partners, kick, ks, dt, iter, mc, recalculate_gs)
Definition: td_write.F90:1025
subroutine, public td_write_init(writ, namespace, space, outp, gr, st, hm, ions, ext_partners, ks, ions_move, with_gauge_field, kick, iter, max_iter, dt, mc)
Initialize files to write when prograting in time.
Definition: td_write.F90:371
integer, parameter, public out_separate_forces
Definition: td_write.F90:201
type(type_t), public type_cmplx
Definition: types.F90:134
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
Definition: unit.F90:132
This module defines the unit system, used for input and output.
type(unit_system_t), public units_out
type(unit_system_t), public units_inp
the units systems for reading and writing
subroutine, public v_ks_freeze_hxc(ks)
Definition: v_ks.F90:1425
subroutine, public v_ks_calculate_current(this, calc_cur)
Definition: v_ks.F90:1436
subroutine, public v_ks_calc(ks, namespace, space, hm, st, ions, ext_partners, calc_eigenval, time, calc_energy, calc_current, force_semilocal)
Definition: v_ks.F90:730
This module provices a simple timer class which can be used to trigger the writing of a restart file ...
Definition: walltimer.F90:121
logical function, public walltimer_alarm(comm, print)
indicate whether time is up
Definition: walltimer.F90:303
Extension of space that contains the knowledge of the spin dimension.
Description of the grid, containing information on derivatives, stencil, and symmetries.
Definition: grid.F90:168
Describes mesh distribution to nodes.
Definition: mesh.F90:186
Stores all communicators and groups.
Definition: multicomm.F90:206
output handler class
Definition: output_low.F90:164
The states_elec_t class contains all electronic wave functions.
int true(void)