Octopus
target.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 target_oct_m
23 use debug_oct_m
24 use epot_oct_m
25 use global_oct_m
26 use grid_oct_m
28 use ions_oct_m
30 use, intrinsic :: iso_fortran_env
32 use lasers_oct_m
33 use mesh_oct_m
40 use parser_oct_m
43 use space_oct_m
60 use td_oct_m
61 use types_oct_m
62 use unit_oct_m
65
66 implicit none
67
68 private
69 public :: &
70 target_t, &
73 target_end, &
76 target_inh, &
79 target_j1, &
80 target_chi, &
84
85contains
86
87 ! ---------------------------------------------------------
92 subroutine target_init_propagation(tg)
93 type(target_t), intent(inout) :: tg
95
96 select case (tg%type)
97 case (oct_tg_hhgnew)
98 tg%vel = m_z0
99 tg%gvec = m_z0
100 tg%acc = m_z0
101 end select
102
104 end subroutine target_init_propagation
105
106
107 ! ----------------------------------------------------------------------
109 subroutine target_get_state(tg, st)
110 type(target_t), intent(in) :: tg
111 type(states_elec_t), intent(inout) :: st
112
113 push_sub(target_get_state)
114 call states_elec_copy(st, tg%st)
115
116 pop_sub(target_get_state)
117 end subroutine target_get_state
118 ! ----------------------------------------------------------------------
119
120
121 ! ----------------------------------------------------------------------
123 subroutine target_init(gr, kpoints, namespace, space, ions, qcs, td, w0, tg, oct, ep, mc)
124 type(grid_t), intent(in) :: gr
125 type(kpoints_t), intent(in) :: kpoints
126 type(namespace_t), intent(in) :: namespace
127 class(space_t), intent(in) :: space
128 type(ions_t), intent(in) :: ions
129 type(opt_control_state_t), intent(inout) :: qcs
130 type(td_t), intent(in) :: td
131 real(real64), intent(in) :: w0
132 type(target_t), intent(inout) :: tg
133 type(oct_t), intent(in) :: oct
134 type(epot_t), intent(inout) :: ep
135 type(multicomm_t), intent(in) :: mc
136
137 integer :: ierr
138 type(states_elec_t), pointer :: stin
139 type(restart_t) :: restart
140
141 push_sub(target_init)
142
143 stin => opt_control_point_qs(qcs)
144
145 !%Variable OCTTargetOperator
146 !%Type integer
147 !%Section Calculation Modes::Optimal Control
148 !%Default oct_tg_gstransformation
149 !%Description
150 !% The variable <tt>OCTTargetOperator</tt> prescribes which kind of target functional is
151 !% to be used.
152 !%Option oct_tg_groundstate 1
153 !% The target operator is a projection operator on the ground state, <i>i.e.</i> the
154 !% objective is to populate the ground state as much as possible.
155 !%Option oct_tg_excited 2
156 !% (Experimental) The target operator is an "excited state". This means that the target operator
157 !% is a linear combination of Slater determinants, each one formed by replacing
158 !% in the ground-state Slater determinant one occupied state with one excited
159 !% state (<i>i.e.</i> "single excitations"). The description of which excitations are
160 !% used, and with which weights, should be given in a file called
161 !% <tt>oct-excited-state-target</tt>.
162 !% See the documentation of subroutine <tt>excited_states_elec_init</tt> in the source
163 !% code in order to use this feature.
164 !%Option oct_tg_gstransformation 3
165 !% The target operator is a projection operator on a transformation of the ground-state
166 !% orbitals defined by the block <tt>OCTTargetTransformStates</tt>.
167 !%Option oct_tg_userdefined 4
168 !% (Experimental) Allows to define target state by using <tt>OCTTargetUserdefined</tt>.
169 !%Option oct_tg_jdensity 5
170 !% (Experimental)
171 !%Option oct_tg_local 6
172 !% (Experimental) The target operator is a local operator.
173 !%Option oct_tg_td_local 7
174 !% (Experimental) The target operator is a time-dependent local operator.
175 !%Option oct_tg_exclude_state 8
176 !% (Experimental) Target operator is the projection onto the complement of a given state, given by the
177 !% block <tt>OCTTargetTransformStates</tt>. This means that the target operator is the unity
178 !% operator minus the projector onto that state.
179 !%Option oct_tg_hhg 9
180 !% (Experimental) The target is the optimization of the HHG yield. You must supply the <tt>OCTOptimizeHarmonicSpectrum</tt>
181 !% block, and it attempts to optimize the maximum of the spectrum around each harmonic peak. You may
182 !% use only one of the gradient-less optimization schemes.
183 !%Option oct_tg_velocity 10
184 !% (Experimental) The target is a function of the velocities of the nuclei at the end of the influence of
185 !% the external field, defined by <tt>OCTVelocityTarget</tt>
186 !%Option oct_tg_hhgnew 12
187 !% (Experimental) The target is the optimization of the HHG yield. You must supply the
188 !% <tt>OCTHarmonicWeight</tt> string. It attempts to optimize the integral of the harmonic spectrum multiplied
189 !% by some user-defined weight function.
190 !%Option oct_tg_classical 13
191 !% (Experimental)
192 !%Option oct_tg_spin 14
193 !% (Experimental)
194 !%End
195 call parse_variable(namespace, 'OCTTargetOperator', oct_tg_gstransformation, tg%type)
196 if (tg%type == oct_tg_excited) call messages_experimental('OCTTargetOperator = oct_tg_excited', namespace=namespace)
197 if (tg%type == oct_tg_userdefined) call messages_experimental('OCTTargetOperator = oct_tg_userdefined')
198 if (tg%type == oct_tg_jdensity) call messages_experimental('OCTTargetOperator = oct_tg_jdensity', namespace=namespace)
199 if (tg%type == oct_tg_local) call messages_experimental('OCTTargetOperator = oct_tg_local', namespace=namespace)
200 if (tg%type == oct_tg_td_local) call messages_experimental('OCTTargetOperator = oct_tg_td_local', namespace=namespace)
201 if (tg%type == oct_tg_exclude_state) call messages_experimental('OCTTargetOperator = oct_tg_exclude_state', namespace=namespace)
202 if (tg%type == oct_tg_hhg) call messages_experimental('OCTTargetOperator = oct_tg_hhg', namespace=namespace)
203 if (tg%type == oct_tg_velocity) call messages_experimental('OCTTargetOperator = oct_tg_velocity', namespace=namespace)
204 if (tg%type == oct_tg_hhgnew) call messages_experimental('OCTTargetOperator = oct_tg_hhgnew', namespace=namespace)
205 if (tg%type == oct_tg_classical) call messages_experimental('OCTTargetOperator = oct_tg_classical', namespace=namespace)
206 if (tg%type == oct_tg_spin) call messages_experimental('OCTTargetOperator = oct_tg_spin', namespace=namespace)
207
208
209 if (.not. varinfo_valid_option('OCTTargetOperator', tg%type)) then
210 call messages_input_error(namespace, 'OCTTargetOperator')
211 end if
212
213 call states_elec_copy(tg%st, stin)
216 call restart_init(restart, namespace, restart_gs, restart_type_load, mc, ierr, mesh=gr, exact=.true.)
217 if (ierr /= 0) then
218 message(1) = "Could not read gs for OCTTargetOperator."
219 call messages_fatal(1, namespace=namespace)
220 end if
221
222 tg%curr_functional = oct_no_curr
223
224 select case (tg%type)
225 case (oct_tg_groundstate)
226 call target_init_groundstate(gr, namespace, space, tg, td, restart, kpoints)
227 case (oct_tg_excited)
228 call messages_experimental('OCTTargetOperator = oct_tg_excited', namespace=namespace)
229 call target_init_excited(gr, namespace, space, tg, td, restart, kpoints)
231 call target_init_exclude(gr, namespace, space, tg, td, restart, kpoints)
233 call target_init_gstransformation(gr, namespace, space, tg, td, restart, kpoints)
234 case (oct_tg_userdefined)
235 call target_init_userdefined(gr, namespace, tg, td)
236 case (oct_tg_jdensity)
237 call target_init_density(gr, kpoints, namespace, space, tg, stin, td, restart)
238 case (oct_tg_local)
239 call target_init_local(gr, namespace, tg, td)
240 case (oct_tg_td_local)
241 call target_init_tdlocal(gr, namespace, tg, td)
242 case (oct_tg_hhg)
243 call target_init_hhg(tg, namespace, td, w0)
244 case (oct_tg_hhgnew)
245 call messages_experimental('OCTTargetOperator = oct_tg_hhgnew', namespace=namespace)
246 call target_init_hhgnew(gr, namespace, tg, td, ions, ep)
247 case (oct_tg_velocity)
248 call target_init_velocity(gr, namespace, ions, tg, oct, td, ep)
249 case (oct_tg_classical)
250 call messages_experimental('OCTTargetOperator = oct_tg_classical', namespace=namespace)
251 call target_init_classical(ions, namespace, tg, td, oct)
252 case (oct_tg_spin)
253 call messages_experimental('OCTTargetOperator = oct_tg_spin', namespace=namespace)
254 call target_init_spin(tg, namespace)
255 case default
256 write(message(1),'(a)') "Target Operator not properly defined."
257 call messages_fatal(1, namespace=namespace)
258 end select
259
260 call restart_end(restart)
261
262 nullify(stin)
263 pop_sub(target_init)
264 end subroutine target_init
265 ! ----------------------------------------------------------------------
266
267
268 ! ----------------------------------------------------------------------
269 subroutine target_end(tg, oct)
270 type(target_t), intent(inout) :: tg
271 type(oct_t), intent(in) :: oct
272
273 push_sub(target_end)
274
275 call states_elec_end(tg%st)
276
277 select case (tg%type)
279 call target_end_exclude()
280 case (oct_tg_jdensity)
281 call target_end_density(tg)
282 case (oct_tg_local)
283 call target_end_local(tg)
284 case (oct_tg_td_local)
285 call target_end_tdlocal(tg)
286 case (oct_tg_hhg)
287 call target_end_hhg(tg)
288 case (oct_tg_hhgnew)
289 call target_end_hhgnew(tg, oct)
290 case (oct_tg_velocity)
291 call target_end_velocity(tg, oct)
292 case (oct_tg_classical)
293 call target_end_classical(tg)
294 end select
295
296 pop_sub(target_end)
297 end subroutine target_end
298 ! ----------------------------------------------------------------------
299
300
301 ! ----------------------------------------------------------------------
302 subroutine target_output(tg, namespace, space, gr, dir, ions, hm, outp)
303 type(target_t), intent(inout) :: tg
304 type(namespace_t), intent(in) :: namespace
305 class(space_t), intent(in) :: space
306 type(grid_t), intent(in) :: gr
307 character(len=*), intent(in) :: dir
308 type(ions_t), intent(in) :: ions
309 type(hamiltonian_elec_t), intent(in) :: hm
310 type(output_t), intent(in) :: outp
311
312 push_sub(target_output)
313
314 select case (tg%type)
315 case (oct_tg_groundstate)
316 call target_output_groundstate(tg, namespace, space, gr, dir, ions, hm, outp)
317 case (oct_tg_excited)
318 call target_output_excited(tg, namespace, space, gr, dir, ions, hm, outp)
320 call target_output_exclude(tg, namespace, space, gr, dir, ions, hm, outp)
322 call target_output_gstransformation(tg, namespace, space, gr, dir, ions, hm, outp)
323 case (oct_tg_userdefined)
324 call target_output_userdefined(tg, namespace, space, gr, dir, ions, hm, outp)
325 case (oct_tg_jdensity)
326 call target_output_density(tg, namespace, space, gr, dir, ions, outp)
327 case (oct_tg_local)
328 call target_output_local(tg, namespace, space, gr, dir, ions, outp)
329 case (oct_tg_td_local)
330 call target_output_tdlocal(tg, namespace, space, gr, dir, ions, outp)
331 case (oct_tg_hhg)
332 call target_output_hhg(tg, namespace, space, gr, dir, ions, hm, outp)
333 case (oct_tg_hhgnew)
334 call target_output_hhg(tg, namespace, space, gr, dir, ions, hm, outp)
335 case (oct_tg_velocity)
336 call target_output_velocity(tg, namespace, space, gr, dir, ions, hm, outp)
337 case (oct_tg_classical)
339 end select
340
341 pop_sub(target_output)
342 end subroutine target_output
343 ! ----------------------------------------------------------------------
344
345
346 ! ---------------------------------------------------------
350 subroutine target_tdcalc(tg, namespace, space, hm, gr, ions, ext_partners, psi, time, max_time)
351 type(target_t), intent(inout) :: tg
352 type(namespace_t), intent(in) :: namespace
353 class(space_t), intent(in) :: space
354 type(hamiltonian_elec_t), intent(inout) :: hm
355 type(grid_t), intent(in) :: gr
356 type(ions_t), intent(inout) :: ions
357 type(partner_list_t), intent(in) :: ext_partners
358 type(states_elec_t), intent(inout) :: psi
359 integer, intent(in) :: time
360 integer, intent(in) :: max_time
361
363
364 push_sub(target_tdcalc)
365
366 tg%td_fitness(time) = m_zero
367
368 select case (tg%type)
369 case (oct_tg_hhgnew)
370 call target_tdcalc_hhgnew(tg, gr, psi, time, max_time)
371 case (oct_tg_velocity)
372 call target_tdcalc_velocity(tg, hm, gr, ions, psi, time, max_time)
373 case (oct_tg_td_local)
374 call target_tdcalc_tdlocal(tg, gr, psi, time)
375 case (oct_tg_hhg)
376 call target_tdcalc_hhg(tg, namespace, space, hm, gr, ions, ext_partners, psi, time)
377 case (oct_tg_jdensity)
378 call target_tdcalc_density(tg, gr, hm%kpoints, psi, time)
379 case default
380 message(1) = 'Error in target.target_tdcalc: default.'
381 call messages_fatal(1, namespace=namespace)
382 end select
383
384 pop_sub(target_tdcalc)
385 end subroutine target_tdcalc
386 ! ----------------------------------------------------------------------
387
388
389
390 ! ---------------------------------------------------------------
393 subroutine target_inh(psi, gr, kpoints, tg, time, inh, iter)
394 type(states_elec_t), intent(inout) :: psi
395 type(grid_t), intent(in) :: gr
396 type(kpoints_t), intent(in) :: kpoints
397 type(target_t), intent(inout) :: tg
398 real(real64), intent(in) :: time
399 type(states_elec_t), intent(inout) :: inh
400 integer, intent(in) :: iter
401
402 integer :: ik, ist, ip, idim, ib
403 complex(real64), allocatable :: zpsi(:)
404 complex(real64) :: gvec(gr%box%dim)
405
406 push_sub(target_inh)
407
408 safe_allocate(zpsi(1:gr%np))
409
410 select case (tg%type)
411 case (oct_tg_td_local)
412
413 call target_build_tdlocal(tg, gr, time)
414
415 do ik = inh%d%kpt%start, inh%d%kpt%end
416 do ist = inh%st_start, inh%st_end
417 do idim = 1, inh%d%dim
418 call states_elec_get_state(psi, gr, idim, ist, ik, zpsi)
419 zpsi(1:gr%np) = -psi%occ(ist, ik)*tg%rho(1:gr%np)*zpsi(1:gr%np)
420 call states_elec_set_state(inh, gr, idim, ist, ik, zpsi)
421 end do
422 end do
423 end do
424
425 case (oct_tg_hhgnew)
426 gvec(:) = real(tg%gvec(iter + 1, :), real64)
427
428 do ik = inh%d%kpt%start, inh%d%kpt%end
429 do ist = inh%st_start, inh%st_end
430 do idim = 1, inh%d%dim
431 call states_elec_get_state(psi, gr, idim, ist, ik, zpsi)
432 do ip = 1, gr%np
433 zpsi(ip) = -psi%occ(ist, ik)*m_two*sum(tg%grad_local_pot(1, ip, 1:gr%box%dim)*gvec(:))*zpsi(ip)
434 end do
435 call states_elec_set_state(inh, gr, idim, ist, ik, zpsi)
436 end do
437 end do
438 end do
439
440 case (oct_tg_velocity)
441
442 do ik = inh%d%kpt%start, inh%d%kpt%end
443 do ist = inh%st_start, inh%st_end
444 do idim = 1, inh%d%dim
445 call states_elec_get_state(psi, gr, idim, ist, ik, zpsi)
446 do ip = 1, gr%np
447 zpsi(ip) = -psi%occ(ist, ik)*tg%rho(ip)*zpsi(ip)
448 end do
449 call states_elec_set_state(inh, gr, idim, ist, ik, zpsi)
450 end do
451 end do
452 end do
453
454 case (oct_tg_jdensity)
455
456 do ik = inh%d%kpt%start, inh%d%kpt%end
457 do ib = inh%group%block_start, inh%group%block_end
458 call batch_set_zero(inh%group%psib(ib, ik))
459 end do
460 end do
461
462 if (abs(nint(time/tg%dt)) >= tg%strt_iter_curr_tg) then
463 call chi_current(tg, gr, kpoints, -1.0_real64, psi, inh)
464 end if
465
466 case default
467 write(message(1),'(a)') 'Internal error in target_inh'
468 call messages_fatal(1)
469
470 end select
471
472 safe_deallocate_a(zpsi)
473
474 pop_sub(target_inh)
475 end subroutine target_inh
476 !----------------------------------------------------------
477
478
479 ! ---------------------------------------------------------
484 real(real64) function target_j1(tg, namespace, gr, kpoints, qcpsi, ions) result(j1)
485 type(target_t), intent(inout) :: tg
486 type(namespace_t), intent(in) :: namespace
487 type(grid_t), intent(in) :: gr
488 type(kpoints_t), intent(in) :: kpoints
489 type(opt_control_state_t), intent(inout) :: qcpsi
490 type(ions_t), optional, intent(in) :: ions
491
492 type(states_elec_t), pointer :: psi
493
494 psi => opt_control_point_qs(qcpsi)
495
496 push_sub(target_j1)
497
498 j1 = m_zero
499 select case (tg%type)
500 case (oct_tg_groundstate)
501 j1 = target_j1_groundstate(tg, gr, psi)
502 case (oct_tg_excited)
503 j1 = target_j1_excited(tg, namespace, gr, psi)
505 j1 = target_j1_gstransformation(tg, gr, psi)
506 case (oct_tg_userdefined)
507 j1 = target_j1_userdefined(tg, gr, psi)
508 case (oct_tg_jdensity)
509 j1 = target_j1_density(gr, kpoints, tg, psi)
510 case (oct_tg_local)
511 j1 = target_j1_local(gr, tg, psi)
512 case (oct_tg_td_local)
513 j1 = target_j1_tdlocal(tg)
515 j1 = target_j1_exclude(gr, tg, psi)
516 case (oct_tg_hhg)
517 j1 = target_j1_hhg(tg, namespace)
518 case (oct_tg_hhgnew)
519 j1 = target_j1_hhgnew(gr, tg)
520 case (oct_tg_velocity)
521 j1 = target_j1_velocity(tg, ions)
522 case (oct_tg_classical)
523 j1 = target_j1_classical(tg, qcpsi)
524 case (oct_tg_spin)
525 j1 = target_j1_spin(tg, gr, psi)
526 end select
527
528 nullify(psi)
529 pop_sub(target_j1)
530 end function target_j1
531 ! ---------------------------------------------------------
532
533
534 ! ---------------------------------------------------------
536 subroutine target_chi(tg, namespace, gr, kpoints, qcpsi_in, qcchi_out, ions)
537 type(target_t), intent(inout) :: tg
538 type(namespace_t), intent(in) :: namespace
539 type(grid_t), intent(in) :: gr
540 type(kpoints_t), intent(in) :: kpoints
541 type(opt_control_state_t), target, intent(inout) :: qcpsi_in
542 type(opt_control_state_t), target, intent(inout) :: qcchi_out
543 type(ions_t), intent(in) :: ions
544
545 real(real64), pointer :: q(:, :), p(:, :)
546 type(states_elec_t), pointer :: psi_in, chi_out
547 push_sub(target_chi)
548
549 psi_in => opt_control_point_qs(qcpsi_in)
550 chi_out => opt_control_point_qs(qcchi_out)
551
552 select case (tg%type)
553 case (oct_tg_groundstate)
554
555 call target_chi_groundstate(tg, gr, psi_in, chi_out)
556 case (oct_tg_excited)
557 call target_chi_excited(tg, namespace, gr, psi_in, chi_out)
558 case (oct_tg_gstransformation)
559 call target_chi_gstransformation(tg, gr, psi_in, chi_out)
560 case (oct_tg_userdefined)
561 call target_chi_userdefined(tg, gr, psi_in, chi_out)
562 case (oct_tg_jdensity)
563 call target_chi_density(tg, gr, kpoints, psi_in, chi_out)
564 case (oct_tg_local)
565 call target_chi_local(tg, gr, psi_in, chi_out)
566 case (oct_tg_td_local)
567 call target_chi_tdlocal(chi_out)
568 case (oct_tg_exclude_state)
569 call target_chi_exclude(tg, gr, psi_in, chi_out)
570 case (oct_tg_hhg)
571 call target_chi_hhg(chi_out)
572 case (oct_tg_hhgnew)
573 call target_chi_hhg(chi_out)
574 case (oct_tg_velocity)
575 call target_chi_velocity(gr, tg, chi_out, ions)
576 case (oct_tg_classical)
577 call target_chi_classical(tg, qcpsi_in, qcchi_out, ions)
578 case (oct_tg_spin)
579 call target_chi_spin(tg, gr, psi_in, chi_out)
580 end select
581
582 ! Unless the target is "classical", the co-state classical variables are zero at time t=T.
583 if (tg%type .ne. oct_tg_classical) then
584 q => opt_control_point_q(qcchi_out)
585 p => opt_control_point_p(qcchi_out)
586 q = m_zero
587 p = m_zero
588 nullify(q)
589 nullify(p)
590 end if
591
592 nullify(psi_in)
593 nullify(chi_out)
594 pop_sub(target_chi)
595 end subroutine target_chi
596
597
598end module target_oct_m
599
600!! Local Variables:
601!! mode: f90
602!! coding: utf-8
603!! End:
This module implements common operations on batches of mesh functions.
Definition: batch_ops.F90:116
subroutine, public batch_set_zero(this, np, async)
fill all mesh functions of the batch with zero
Definition: batch_ops.F90:242
real(real64), parameter, public m_two
Definition: global.F90:190
real(real64), parameter, public m_zero
Definition: global.F90:188
complex(real64), parameter, public m_z0
Definition: global.F90:198
This module implements the underlying real-space grid.
Definition: grid.F90:117
This module defines classes and functions for interaction partners.
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:118
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
This module handles the communicators for the various parallelization strategies.
Definition: multicomm.F90:145
This module contains the definition of the oct_t data type, which contains some of the basic informat...
This module holds the "opt_control_state_t" datatype, which contains a quantum-classical state.
type(states_elec_t) function, pointer, public opt_control_point_qs(ocs)
this module contains the low-level part of the output system
Definition: output_low.F90:115
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_load
Definition: restart.F90:245
subroutine, public restart_end(restart)
Definition: restart.F90:722
subroutine, public states_elec_end(st)
finalize the states_elec_t 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.
subroutine, public states_elec_copy(stout, stin, exclude_wfns, exclude_eigenval, special)
make a (selective) copy of a states_elec_t object
real(real64) function, public target_j1_classical(tg, qcpsi)
subroutine, public target_init_classical(ions, namespace, tg, td, oct)
subroutine, public target_end_classical(tg)
subroutine, public target_output_classical
subroutine, public chi_current(tg, gr, kpoints, factor, psi_in, chi)
real(real64) function, public target_j1_density(gr, kpoints, tg, psi)
subroutine, public target_end_density(tg)
subroutine, public target_init_density(gr, kpoints, namespace, space, tg, stin, td, restart)
subroutine, public target_output_density(tg, namespace, space, mesh, dir, ions, outp)
subroutine, public target_tdcalc_density(tg, gr, kpoints, psi, time)
subroutine, public target_init_excited(mesh, namespace, space, tg, td, restart, kpoints)
subroutine, public target_output_excited(tg, namespace, space, gr, dir, ions, hm, outp)
real(real64) function, public target_j1_excited(tg, namespace, gr, psi)
subroutine, public target_init_exclude(mesh, namespace, space, tg, td, restart, kpoints)
subroutine, public target_end_exclude()
real(real64) function, public target_j1_exclude(gr, tg, psi)
subroutine, public target_output_exclude(tg, namespace, space, gr, dir, ions, hm, outp)
subroutine, public target_output_groundstate(tg, namespace, space, gr, dir, ions, hm, outp)
subroutine, public target_init_groundstate(mesh, namespace, space, tg, td, restart, kpoints)
real(real64) function, public target_j1_groundstate(tg, gr, psi)
real(real64) function, public target_j1_gstransformation(tg, gr, psi)
subroutine, public target_output_gstransformation(tg, namespace, space, gr, dir, ions, hm, outp)
subroutine, public target_init_gstransformation(gr, namespace, space, tg, td, restart, kpoints)
subroutine, public target_end_hhg(tg)
Definition: target_hhg.F90:342
subroutine, public target_init_hhg(tg, namespace, td, w0)
Definition: target_hhg.F90:167
subroutine, public target_tdcalc_hhgnew(tg, gr, psi, time, max_time)
Definition: target_hhg.F90:476
subroutine, public target_tdcalc_hhg(tg, namespace, space, hm, gr, ions, ext_partners, psi, time)
Definition: target_hhg.F90:544
real(real64) function, public target_j1_hhg(tg, namespace)
Definition: target_hhg.F90:402
subroutine, public target_init_hhgnew(gr, namespace, tg, td, ions, ep)
Definition: target_hhg.F90:246
real(real64) function, public target_j1_hhgnew(gr, tg)
Definition: target_hhg.F90:434
subroutine, public target_end_hhgnew(tg, oct)
Definition: target_hhg.F90:379
subroutine, public target_output_hhg(tg, namespace, space, gr, dir, ions, hm, outp)
Definition: target_hhg.F90:357
real(real64) function, public target_j1_local(mesh, tg, psi)
subroutine, public target_init_local(gr, namespace, tg, td)
subroutine, public target_output_local(tg, namespace, space, mesh, dir, ions, outp)
subroutine, public target_end_local(tg)
integer, parameter, public oct_tg_velocity
Definition: target_low.F90:212
integer, parameter, public oct_tg_hhgnew
Definition: target_low.F90:212
integer, parameter, public oct_tg_excited
Definition: target_low.F90:212
integer, parameter, public oct_tg_hhg
Definition: target_low.F90:212
integer pure function, public target_curr_functional(tg)
Definition: target_low.F90:273
integer, parameter, public oct_tg_spin
Definition: target_low.F90:212
integer pure function, public target_type(tg)
Definition: target_low.F90:265
integer, parameter, public oct_targetmode_td
Definition: target_low.F90:227
integer, parameter, public oct_tg_groundstate
Definition: target_low.F90:212
integer, parameter, public oct_tg_classical
Definition: target_low.F90:212
integer, parameter, public oct_no_curr
Definition: target_low.F90:231
integer, parameter, public oct_tg_exclude_state
Definition: target_low.F90:212
integer pure function, public target_mode(tg)
Definition: target_low.F90:243
integer, parameter, public oct_tg_jdensity
Definition: target_low.F90:212
integer, parameter, public oct_tg_gstransformation
Definition: target_low.F90:212
logical pure function, public target_move_ions(tg)
Definition: target_low.F90:281
integer, parameter, public oct_tg_local
Definition: target_low.F90:212
integer, parameter, public oct_tg_td_local
Definition: target_low.F90:212
integer, parameter, public oct_tg_userdefined
Definition: target_low.F90:212
subroutine, public target_init_propagation(tg)
This routine performs all the things that must be initialized prior to a forward evolution,...
Definition: target.F90:186
real(real64) function, public target_j1(tg, namespace, gr, kpoints, qcpsi, ions)
Calculates the J1 functional, i.e.: in the time-independent case, or else in the time-dependent cas...
Definition: target.F90:578
subroutine, public target_init(gr, kpoints, namespace, space, ions, qcs, td, w0, tg, oct, ep, mc)
The target is initialized, mainly by reading from the inp file.
Definition: target.F90:217
subroutine, public target_get_state(tg, st)
This just copies the states_elec_t variable present in target, into st.
Definition: target.F90:203
subroutine, public target_inh(psi, gr, kpoints, tg, time, inh, iter)
Calculates the inhomogeneous term that appears in the equation for chi, and places it into inh.
Definition: target.F90:487
subroutine, public target_output(tg, namespace, space, gr, dir, ions, hm, outp)
Definition: target.F90:396
subroutine, public target_tdcalc(tg, namespace, space, hm, gr, ions, ext_partners, psi, time, max_time)
Calculates, at a given point in time marked by the integer index, the integrand of the target functio...
Definition: target.F90:444
subroutine, public target_chi(tg, namespace, gr, kpoints, qcpsi_in, qcchi_out, ions)
Calculate .
Definition: target.F90:630
subroutine, public target_end(tg, oct)
Definition: target.F90:363
subroutine, public target_init_spin(tg, namespace)
real(real64) function, public target_j1_spin(tg, gr, psi)
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)
real(real64) function, public target_j1_tdlocal(tg)
subroutine, public target_build_tdlocal(tg, gr, time)
real(real64) function, public target_j1_userdefined(tg, gr, psi)
subroutine, public target_init_userdefined(gr, namespace, tg, td)
subroutine, public target_output_userdefined(tg, namespace, space, gr, dir, ions, hm, outp)
real(real64) function, public target_j1_velocity(tg, ions)
subroutine, public target_output_velocity(tg, namespace, space, gr, dir, ions, hm, outp)
subroutine, public target_tdcalc_velocity(tg, hm, gr, ions, psi, time, max_time)
subroutine, public target_init_velocity(gr, namespace, ions, tg, oct, td, ep)
subroutine, public target_end_velocity(tg, oct)
Definition: td.F90:114
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.
Description of the grid, containing information on derivatives, stencil, and symmetries.
Definition: grid.F90:168
!brief The oct_t datatype stores the basic information about how the OCT run is done.
This is the datatype that contains the objects that are propagated: in principle this could be both t...
output handler class
Definition: output_low.F90:164
The states_elec_t class contains all electronic wave functions.
int true(void)