How To Calculate and Write Out Expectation Values¶
Observables in the config file¶
Obtaining time-dependent expectation values for arbitrary observables
while propagating a system is possible by defining observables in the
config file. For example, when propagating a wave packet on a single
potential surface, we may define the following `observables
section <../config_file.html#observables>`__ in the config file.
observables: is_real = T, outfile = exp_values.dat, time_unit = ns
* type = r, column_label='<r>', square='<r^2>', exp_unit = nm
* type = mom, column_label='<p>', square='<p^2>'
After the propagation (i.e. a call to
`prop
<../reference/prop.html#prop>`__), the file exp_values.dat
will have been created and contain 5 columns: the time, and the
expectation values of the position operator, the square of the position
operator, the momentum operator, and the square of the momentum
operator.
The expectation values of the squares may be used to later calculate the
variance of the position or momentum operator. Whether or not these
values are calculated and written out is determined entirely by the
presence of the square
parameter (which gives the column label). If
square
is not given, the expectation value of the square is not
calculated, and the outfile
will only contain a single column for
that observable.
It is possible to write different expectation values to different files
by specifying a different outfile
for each observable. It is
recommended to group all observables whose expectation value go into the
same file into one section in the config file.
A large variety of observables may be defined in the config file,
through the type
parameter. This includes r
and mom
as in
the above example, as well as ham
, and kin
to calculate energy
expectation values (total and kinetic). Moreover, any grid-based
operator may be specified, using the same parameters as in the `ham
section <../config_file.html#ham>`__. For example, the following section
would request the expectation values of the total, kinetic, and
potential energy of a simple harmonic oscillator:
observables: is_real = T, outfile = exp_energies.dat
* type = ham, column_label='<H>'
* type = kin, column_label='<T>'
* type = op, column_label='<V>', &
op_type = pot, op_form = hamos, op_surf = 1, w_0 = 10
The most general definition of an observable is through type=matrix
,
which allows to define an arbitrary operator over the entire Hilbert
space, read as a sparse matrix from file (cf. type=matrix
in the
`ham
section <../config_file.html#ham>`__ of the config file).
Lastly, type=norm
and type=pop
are also available. The former
calculates the norm of the wave function, the latter must be paired with
a specification of exp_surf
and calculates the population on the
given energy surface.
The exp_surf
parameter can be used more generally to calculate
expectation values on a single energy surface, e.g.
observables: is_real = T, outfile = exp_values.dat, time_unit = ns, &
& type = r, exp_unit = nm
* column_label='<r> (tot)'
* column_label='<r> (sf1)', square='<r^2> (sf1)', exp_surf = 1
* column_label='<r> (sf2)', square='<r^2> (sf2)', exp_surf = 2
In this case, before calculating the expectation value, the wave function is projected onto the given surface and then renormalized.
Observables can be non-Hermitian, resulting in complex expectation
values. This is possible for type=ham
if the Hamiltonian is
non-Hermitian, or with type=matrix
, which allows to define
completely arbitrary matrices. In this case, is_real = F
must be
specified. The output file will then contain two columns for the
observable, one for the real part of the expectation value, and one for
the imaginary part. If square
is also defined, the number of columns
for that observable will increase to a total of four.
The information about which expectation values should be calculated is
determined from the config file data by the `init_prop
routine <../reference/prop.html#init_prop>`__. In case any observable
section in the config file uses the label
parameter, only the items
that match gen%ham%label
in init_prop
are taken into
consideration. The information is stored as part of the prop_work
variable.
The calculation and output of the expectation values are handled by the
`prop
routine <../reference/prop.html#prop>`__. Every time prop
is called, the outfile
for any of the observables is overwritten.
While the propagation is running, data is continuously written to each
outfile
. Thus, running tail -f {outfile}
may be a useful way to
monitor the progress of a long-running propagation. On completion of the
prop
routine, the expectation values for the entire time grid
covered by the propagation will be contained in outfile
. For a
backwards propagation the time (first column in outfile
) will also
run backwards.
Accessing expectation values programatically¶
In addition to writing expectation values to file, the `prop
routine <../reference/prop.html#prop>`__ also stores all values in
memory, in the complex matrix prop_work%exp%values
. The number of
rows of the matrix is equal to the number of time grid points. The first
columns of the matrix contains the time grid values. The subsequent
columns contain the expectation values for all observables, in the order
in which they appear in the config file. Any observable that defines a
square
label has an additional column for the expectation value of
the square of the observable. A label for all the columns may be
obtained from the array prop_work%exp%column_labels
.
Note that since prop_work%exp%values
is complex, any observable with
is_real=F
does not result in two columns in the matrix. When
processing columns that are known to contain the expectation values of
Hermitian operators, the values may be converted to real
on demand.
All values are stored in internal units.
Expectation values will be included in prop_work%exp%values
even if
the observable does not define an outfile
. In this case, the
expectation value will not be written to file, but will be available in
memory. It is even possible to not define an outfile
for any of the
observables, in order to not generate any output files, but only
calculate expectation values internally.
The content of prop_work%exp%values
is overwritten whenever the
prop
routine is called. Any values from an earlier call are lost.
Using prop_info_hook
¶
While in most cases, the output of expectation values can entirely be
specified in the config file, it can be useful to calculate some custom
quantities using the prop_info_hook
parameter of the `prop
routine <../reference/prop.html#prop>`__. The parameter specifies a
routine that is called after each propagation step. Inside the routine,
expectation values could be calculated using the routines defined in the
`exp
module <../reference/exp.html>`__ (use qdyn_exp_mod
).
However, the most common usage of the prop_info_hook
is to write out
the wave function after each propagation step. E.g., using an external
routine
subroutine prop_write_psi(psi, grid, ham, ti, work, para)
use qdyn
use qdyn_exp_mod
implicit none
type(state_t), intent(inout) :: psi
type(grid_t), intent(in) :: grid
type(ham_t), intent(inout) :: ham
integer, intent(in) :: ti
type(prop_work_t), intent(inout) :: work
type(para_t), intent(in) :: para
real(idp) :: t, r
type(state_t), save :: state_temp
character(len=datline_l) :: comment
t = para%tgrid%t_start + ti * para%tgrid%dt
call mold_state(state_temp, state)
call exp_r(r, state, state_temp, grid)
comment='# t = '//trim(nr2str(t))//'; <r> = '//trim(nr2str(r))
if (ti == 0) then
call write_psi(psi, grid, 1, 1, 'psi1.dat', &
& append=.false., comment=comment)
else
call write_psi(psi, grid, 1, 1, 'psi1.dat', &
& append=.true., comment=comment)
end if
end subroutine prop_write_psi
the propagation may be called as
call prop(psi, grid, gen, work, para, pulses, prop_info_hook=prop_write_psi)