Octopus
target_tdlocal.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
20#include "global.h"
21
24 use debug_oct_m
26 use global_oct_m
27 use grid_oct_m
28 use io_oct_m
29 use ions_oct_m
32 use, intrinsic :: iso_fortran_env
33 use mesh_oct_m
38 use parser_oct_m
40 use space_oct_m
43 use td_oct_m
44 use unit_oct_m
46
47
48 implicit none
49
50 private
51 public :: &
59
60
61contains
62
63 ! ----------------------------------------------------------------------
65 subroutine target_init_tdlocal(gr, namespace, tg, td)
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
70
71 type(block_t) :: blk
72 push_sub(target_init_tdlocal)
73
74 tg%move_ions = td%ions_dyn%ions_move()
75 tg%dt = td%dt
76
77 !%Variable OCTTdTarget
78 !%Type block
79 !%Section Calculation Modes::Optimal Control
80 !%Description
81 !% (Experimental) If <tt>OCTTargetOperator = oct_tg_td_local</tt>, then you must supply
82 !% a OCTTdTarget block. The block should only contain one element, a string cotaining the
83 !% definition of the time-dependent local target, <i>i.e.</i> a function of x,y,z and t that
84 !% is to be maximized along the evolution.
85 !%End
86 if (parse_block(namespace, 'OCTTdTarget', blk) == 0) then
87 call parse_block_string(blk, 0, 0, tg%td_local_target, convert_to_c=.true.)
88 safe_allocate(tg%rho(1:gr%np))
89 call parse_block_end(blk)
90 else
91 message(1) = 'If OCTTargetOperator = oct_tg_td_local, you must supply a OCTTdTarget block.'
92 call messages_fatal(1, namespace=namespace)
93 end if
94 safe_allocate(tg%td_fitness(0:td%max_iter))
95 tg%td_fitness = m_zero
96 call target_build_tdlocal(tg, gr, m_zero)
97
98 pop_sub(target_init_tdlocal)
99 end subroutine target_init_tdlocal
100
101
102 ! ----------------------------------------------------------------------
104 subroutine target_end_tdlocal(tg)
105 type(target_t), intent(inout) :: tg
106
107 push_sub(target_end_tdlocal)
108
109 safe_deallocate_a(tg%rho)
110 safe_deallocate_a(tg%td_fitness)
111
112 pop_sub(target_end_tdlocal)
113 end subroutine target_end_tdlocal
114
115
116 ! ----------------------------------------------------------------------
117 subroutine target_output_tdlocal(tg, namespace, space, gr, dir, ions, outp)
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
125
126 integer :: ierr
127 push_sub(target_output_tdlocal)
128
129 call io_mkdir(trim(dir), namespace)
130 call target_build_tdlocal(tg, gr, m_zero)
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)
133
134 pop_sub(target_output_tdlocal)
135 end subroutine target_output_tdlocal
136 ! ----------------------------------------------------------------------
137
138
139 ! ----------------------------------------------------------------------
141 real(real64) function target_j1_tdlocal(tg) result(j1)
142 type(target_t), intent(in) :: tg
143
144 integer :: maxiter
145 push_sub(target_j1_tdlocal)
146
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))
151
152
153 pop_sub(target_j1_tdlocal)
154 end function target_j1_tdlocal
155
156
157 ! ----------------------------------------------------------------------
159 subroutine target_chi_tdlocal(chi_out)
160 type(states_elec_t), intent(inout) :: chi_out
161
162 integer :: ik, ib
163 push_sub(target_chi_tdlocal)
164
165 !We assume that there is no time-independent operator.
166
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))
170 end do
171 end do
172
173 pop_sub(target_chi_tdlocal)
174 end subroutine target_chi_tdlocal
175
176
177 ! ---------------------------------------------------------
180 subroutine target_tdcalc_tdlocal(tg, gr, psi, time)
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
185
186 complex(real64), allocatable :: opsi(:, :), zpsi(:, :)
187 integer :: ist, ip
188 push_sub(target_tdcalc_tdlocal)
189
190 tg%td_fitness(time) = m_zero
191
192 safe_allocate(zpsi(1:gr%np, 1:psi%d%dim))
193
194 !!!! WARNING Here one should build the time-dependent target.
195 select case (psi%d%ispin)
196 case (unpolarized)
197 assert(psi%nik == 1)
198 safe_allocate(opsi(1:gr%np_part, 1))
199 opsi = m_z0
200 do ist = psi%st_start, psi%st_end
201
202 call states_elec_get_state(psi, gr, ist, 1, zpsi)
203
204 do ip = 1, gr%np
205 opsi(ip, 1) = tg%rho(ip)*zpsi(ip, 1)
206 end do
207
208 tg%td_fitness(time) = tg%td_fitness(time) + psi%occ(ist, 1)*real(zmf_dotp(gr, psi%d%dim, zpsi, opsi), real64)
209
210 end do
211 safe_deallocate_a(opsi)
212 case (spin_polarized)
213 message(1) = 'Error in target.target_tdcalc: spin_polarized.'
214 call messages_fatal(1)
215 case (spinors)
216 message(1) = 'Error in target.target_tdcalc: spinors.'
217 call messages_fatal(1)
218 end select
219
220 safe_deallocate_a(zpsi)
221
222 pop_sub(target_tdcalc_tdlocal)
223 end subroutine target_tdcalc_tdlocal
224 ! ----------------------------------------------------------------------
225
226
227 !----------------------------------------------------------
228 subroutine target_build_tdlocal(tg, gr, time)
229 type(target_t), intent(inout) :: tg
230 type(grid_t), intent(in) :: gr
231 real(real64), intent(in) :: time
232
233 integer :: ip
234 real(real64) :: xx(gr%box%dim), rr, re, im
235
237
238 do ip = 1, gr%np
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)
241 tg%rho(ip) = re
242 end do
243
244 pop_sub(target_build_tdlocal)
245 end subroutine target_build_tdlocal
246 !----------------------------------------------------------
247
248end module target_tdlocal_oct_m
249
250!! Local Variables:
251!! mode: f90
252!! coding: utf-8
253!! End:
This module implements common operations on batches of mesh functions.
Definition: batch_ops.F90:118
real(real64), parameter, public m_zero
Definition: global.F90:191
real(real64), parameter, public m_half
Definition: global.F90:197
This module implements the underlying real-space grid.
Definition: grid.F90:119
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.
Definition: io.F90:116
subroutine, public io_mkdir(fname, namespace, parents)
Definition: io.F90:360
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:120
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
Definition: messages.F90:162
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
Definition: messages.F90:410
this module contains the low-level part of the output system
Definition: output_low.F90:117
subroutine, public parse_block_string(blk, l, c, res, convert_to_c)
Definition: parser.F90:810
integer function, public parse_block(namespace, name, blk, check_varinfo_)
Definition: parser.F90:615
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)
Definition: td.F90:116
brief This module defines the class unit_t which is used by the unit_systems_oct_m module.
Definition: unit.F90:134
This module defines the unit system, used for input and output.
type(unit_system_t), public units_out
int true(void)