32 use,
intrinsic :: iso_fortran_env
66 type(grid_t),
intent(in) :: gr
67 type(namespace_t),
intent(in) :: namespace
68 type(target_t),
intent(inout) :: tg
69 type(td_t),
intent(in) :: td
74 tg%move_ions = td%ions_dyn%ions_move()
86 if (
parse_block(namespace,
'OCTTdTarget', blk) == 0)
then
88 safe_allocate(tg%rho(1:gr%np))
91 message(1) =
'If OCTTargetOperator = oct_tg_td_local, you must supply a OCTTdTarget block.'
94 safe_allocate(tg%td_fitness(0:td%max_iter))
105 type(target_t),
intent(inout) :: tg
109 safe_deallocate_a(tg%rho)
110 safe_deallocate_a(tg%td_fitness)
118 type(target_t),
intent(inout) :: tg
119 type(namespace_t),
intent(in) :: namespace
120 class(space_t),
intent(in) :: space
121 type(grid_t),
intent(in) :: gr
122 character(len=*),
intent(in) :: dir
123 type(ions_t),
intent(in) :: ions
124 type(output_t),
intent(in) :: outp
131 call dio_function_output(outp%how(0), trim(dir),
'td_local_target', namespace, space, gr, &
132 tg%rho,
units_out%length**(-space%dim), ierr, pos=ions%pos, atoms=ions%atom)
141 real(real64) function target_j1_tdlocal(tg) result(j1)
142 type(target_t),
intent(in) :: tg
145 push_sub(target_j1_tdlocal)
147 maxiter =
size(tg%td_fitness) - 1
148 j1 =
m_half * tg%dt * tg%td_fitness(0) + &
149 m_half * tg%dt * tg%td_fitness(maxiter) + &
150 tg%dt * sum(tg%td_fitness(1:maxiter-1))
153 pop_sub(target_j1_tdlocal)
160 type(states_elec_t),
intent(inout) :: chi_out
167 do ik = chi_out%d%kpt%start, chi_out%d%kpt%end
168 do ib = chi_out%group%block_start, chi_out%group%block_end
169 call batch_set_zero(chi_out%group%psib(ib, ik))
181 type(target_t),
intent(inout) :: tg
182 type(grid_t),
intent(in) :: gr
183 type(states_elec_t),
intent(in) :: psi
184 integer,
intent(in) :: time
186 complex(real64),
allocatable :: opsi(:, :), zpsi(:, :)
190 tg%td_fitness(time) = m_zero
192 safe_allocate(zpsi(1:gr%np, 1:psi%d%dim))
195 select case (psi%d%ispin)
198 safe_allocate(opsi(1:gr%np_part, 1))
200 do ist = psi%st_start, psi%st_end
202 call states_elec_get_state(psi, gr, ist, 1, zpsi)
205 opsi(ip, 1) = tg%rho(ip)*zpsi(ip, 1)
208 tg%td_fitness(time) = tg%td_fitness(time) + psi%occ(ist, 1)*real(zmf_dotp(gr, psi%d%dim, zpsi, opsi), real64)
211 safe_deallocate_a(opsi)
212 case (spin_polarized)
213 message(1) =
'Error in target.target_tdcalc: spin_polarized.'
214 call messages_fatal(1)
216 message(1) =
'Error in target.target_tdcalc: spinors.'
217 call messages_fatal(1)
220 safe_deallocate_a(zpsi)
229 type(target_t),
intent(inout) :: tg
230 type(grid_t),
intent(in) :: gr
231 real(real64),
intent(in) :: time
234 real(real64) :: xx(gr%box%dim), rr, re, im
239 call mesh_r(gr, ip, rr, coords = xx)
240 call parse_expression(re, im, gr%box%dim, xx, rr, time, tg%td_local_target)
This module implements common operations on batches of mesh functions.
real(real64), parameter, public m_zero
real(real64), parameter, public m_half
This module implements the underlying real-space grid.
subroutine, public dio_function_output(how, dir, fname, namespace, space, mesh, ff, unit, ierr, pos, atoms, grp, root)
Top-level IO routine for functions defined on the mesh.
subroutine, public io_mkdir(fname, namespace, parents)
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
this module contains the low-level part of the output system
subroutine, public parse_block_string(blk, l, c, res, convert_to_c)
integer function, public parse_block(namespace, name, blk, check_varinfo_)
subroutine, public target_output_tdlocal(tg, namespace, space, gr, dir, ions, outp)
subroutine, public target_init_tdlocal(gr, namespace, tg, td)
subroutine, public target_end_tdlocal(tg)
subroutine, public target_tdcalc_tdlocal(tg, gr, psi, time)
subroutine, public target_chi_tdlocal(chi_out)
real(real64) function, public target_j1_tdlocal(tg)
subroutine, public target_build_tdlocal(tg, gr, time)
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
This module defines the unit system, used for input and output.
type(unit_system_t), public units_out