Skip to content

n_coeffs_deriv sorting #83

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
merged 4 commits into from
May 20, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
67 changes: 34 additions & 33 deletions filter_functions/gradient.py
Original file line number Diff line number Diff line change
Expand Up @@ -262,7 +262,7 @@ def _control_matrix_at_timestep_derivative(
The individual control matrices of all time steps
ctrlmat_g_deriv: ndarray, shape (n_dt, n_nops, d**2, n_ctrl, n_omega)
The corresponding derivative with respect to the control
strength :math:`\frac{\partial\mathcal{B}_{\alpha j}^{(g)}(\omega)}
strength :math:`\frac{\partial\mathcal{B}_{\alpha j}^{(g)}(\omega)}`

Notes
-----
Expand Down Expand Up @@ -384,8 +384,6 @@ def calculate_derivative_of_control_matrix_from_scratch(
n_opers: Sequence[Operator],
n_coeffs: Sequence[Coefficients],
c_opers: Sequence[Operator],
all_identifiers: Sequence[str],
control_identifiers: Optional[Sequence[str]] = None,
n_coeffs_deriv: Optional[Sequence[Coefficients]] = None,
intermediates: Optional[Dict[str, ndarray]] = None
) -> ndarray:
Expand Down Expand Up @@ -421,28 +419,16 @@ def calculate_derivative_of_control_matrix_from_scratch(
n_coeffs: array_like, shape (n_nops, n_dt)
The sensitivities of the system to the noise operators given by
*n_opers* at the given time step.
c_opers: array_like, shape (n_cops, d, d)
Control operators :math:`H_k`.
all_identifiers: array_like, shape (n_cops)
Identifiers of all control operators.
control_identifiers: Sequence[str], shape (n_ctrl), Optional
Sequence of strings with the control identifiers to distinguish
between accessible control and drift Hamiltonian. The default is
None.
c_opers: array_like, shape (n_ctrl, d, d)
Control operators :math:`H_k` with respect to which the
derivative is computed.
n_coeffs_deriv: array_like, shape (n_nops, n_ctrl, n_dt)
The derivatives of the noise susceptibilities by the control
amplitudes. Defaults to None.
intermediates: Dict[str, ndarray], optional
Optional dictionary containing intermediate results of the
calculation of the control matrix.

Raises
------
ValueError
If the given identifiers *control_identifier* are not subset of
the total identifiers *all_identifiers* of all control
operators.

Returns
-------
ctrlmat_deriv: ndarray, shape (n_ctrl, n_omega, n_dt, n_nops, d**2)
Expand Down Expand Up @@ -470,24 +456,16 @@ def calculate_derivative_of_control_matrix_from_scratch(
_liouville_derivative
_control_matrix_at_timestep_derivative
"""
# Distinction between control and drift operators and only
# calculate the derivatives in control direction
try:
idx = util.get_indices_from_identifiers(all_identifiers, control_identifiers)
except ValueError as err:
raise ValueError('Given control identifiers have to be a subset of (drift+control) ' +
'Hamiltonian!') from err

d = eigvecs.shape[-1]
n_dt = len(dt)
n_ctrl = len(idx)
n_ctrl = len(c_opers)
n_nops = len(n_opers)
n_omega = len(omega)

# Precompute some transformations or grab from cache if possible
basis_transformed = numeric._transform_by_unitary(eigvecs[:, None], basis[None],
out=np.empty((n_dt, d**2, d, d), complex))
c_opers_transformed = numeric._transform_hamiltonian(eigvecs, c_opers[idx]).swapaxes(0, 1)
c_opers_transformed = numeric._transform_hamiltonian(eigvecs, c_opers).swapaxes(0, 1)
if not intermediates:
# None or empty
n_opers_transformed = numeric._transform_hamiltonian(eigvecs, n_opers,
Expand Down Expand Up @@ -575,6 +553,7 @@ def infidelity_derivative(
spectrum: Coefficients,
omega: Coefficients,
control_identifiers: Optional[Sequence[str]] = None,
n_oper_identifiers: Optional[Sequence[str]] = None,
n_coeffs_deriv: Optional[Sequence[Coefficients]] = None
) -> ndarray:
r"""Calculate the entanglement infidelity derivative of the
Expand All @@ -599,11 +578,29 @@ def infidelity_derivative(
omega: array_like, shape (n_omega,)
The frequencies at which the integration is to be carried out.
control_identifiers: Sequence[str], shape (n_ctrl,)
Sequence of strings with the control identifiern to distinguish
between accessible control and drift Hamiltonian.
Sequence of strings with the control identifiers to
distinguish between control and drift Hamiltonian. The
default is None, in which case the derivative is computed
for all known non-noise operators.
n_oper_identifiers: Sequence[str], shape (n_nops,)
Sequence of strings with the noise identifiers for which to
compute the derivative contribution. The default is None, in
which case it is computed for all known noise operators.
n_coeffs_deriv: array_like, shape (n_nops, n_ctrl, n_dt)
The derivatives of the noise susceptibilities by the control
amplitudes. Defaults to None.
amplitudes. The rows and columns should be in the same order
as the corresponding identifiers above. Defaults to None, in
which case the coefficients are assumed to be constant and
hence their derivative vanishing.

.. warning::

Internally, control and noise terms of the Hamiltonian
are stored alphanumerically sorted by their identifiers.
If the noise and/or control identifiers above are not
explicitly given, the rows and/or columns of this
parameter need to be sorted in the same fashion.


Raises
------
Expand All @@ -616,7 +613,9 @@ def infidelity_derivative(
infid_deriv: ndarray, shape (n_nops, n_dt, n_ctrl)
Array with the derivative of the infidelity for each noise
source taken for each control direction at each time step
:math:`\frac{\partial I_e}{\partial u_h(t_{g'})}`.
:math:`\frac{\partial I_e}{\partial u_h(t_{g'})}`. Sorted in
the same fashion as `n_coeffs_deriv` or, if not given,
alphanumerically by the identifiers.

Notes
-----
Expand Down Expand Up @@ -658,7 +657,9 @@ def infidelity_derivative(
https://doi.org/10.1016/S0375-9601(02)01272-0
"""
spectrum = numeric._parse_spectrum(spectrum, omega, range(len(pulse.n_opers)))
filter_function_deriv = pulse.get_filter_function_derivative(omega, control_identifiers,
filter_function_deriv = pulse.get_filter_function_derivative(omega,
control_identifiers,
n_oper_identifiers,
n_coeffs_deriv)

integrand = np.einsum('ao,atho->atho', spectrum, filter_function_deriv)
Expand Down
84 changes: 68 additions & 16 deletions filter_functions/pulse_sequence.py
Original file line number Diff line number Diff line change
Expand Up @@ -160,17 +160,23 @@ class PulseSequence:
Attributes
----------
c_opers: ndarray, shape (n_cops, d, d)
Control operators
Control operators. Note that they are stored sorted by their
corresponding identifiers.
n_opers: ndarray, shape (n_nops, d, d)
Noise operators
Noise operators. Note that they are stored sorted by their
corresponding identifiers.
c_oper_identifers: sequence of str
Identifiers for the control operators of the system
Identifiers for the control operators of the system. Stored
sorted.
n_oper_identifers: sequence of str
Identifiers for the noise operators of the system
Identifiers for the noise operators of the system. Stored
sorted.
c_coeffs: ndarray, shape (n_cops, n_dt)
Control parameters in units of :math:`\hbar`
Control parameters in units of :math:`\hbar`. Note that they
are stored sorted by their corresponding identifiers.
n_coeffs: ndarray, shape (n_nops, n_dt)
Noise sensitivities in units of :math:`\hbar`
Noise sensitivities in units of :math:`\hbar`. Note that they
are stored sorted by their corresponding identifiers.
dt: ndarray, shape (n_dt,)
Time steps
t: ndarray, shape (n_dt + 1,)
Expand Down Expand Up @@ -247,8 +253,10 @@ def __init__(self, *args, **kwargs) -> None:
self.basis = None

# Parse the input arguments and set attributes
attributes = ('c_opers', 'c_oper_identifiers', 'c_coeffs', 'n_opers',
'n_oper_identifiers', 'n_coeffs', 'dt', 'd', 'basis')
# TODO: Jesus, this is in need of refactoring.
attributes = ('c_opers', 'c_oper_identifiers', 'c_coeffs',
'n_opers', 'n_oper_identifiers', 'n_coeffs',
'dt', 'd', 'basis')
if not args:
# Bypass args parsing and directly set necessary attributes
values = (kwargs[attr] for attr in attributes)
Expand Down Expand Up @@ -855,6 +863,7 @@ def get_filter_function_derivative(
self,
omega: Coefficients,
control_identifiers: Optional[Sequence[str]] = None,
n_oper_identifiers: Optional[Sequence[str]] = None,
n_coeffs_deriv: Optional[Sequence[Coefficients]] = None
) -> ndarray:
r"""Calculate the pulse sequence's filter function derivative.
Expand All @@ -864,27 +873,70 @@ def get_filter_function_derivative(
omega: array_like, shape (n_omega,)
Frequencies at which the pulse control matrix is to be
evaluated.
control_identifiers: Sequence[str]
Sequence of strings with the control identifiern to
control_identifiers: Sequence[str], shape (n_ctrl,)
Sequence of strings with the control identifiers to
distinguish between control and drift Hamiltonian. The
default is None.
default is None, in which case the derivative is computed
for all known non-noise operators.
n_oper_identifiers: Sequence[str], shape (n_nops,)
Sequence of strings with the noise identifiers for which to
compute the derivative contribution. The default is None, in
which case it is computed for all known noise operators.
n_coeffs_deriv: array_like, shape (n_nops, n_ctrl, n_dt)
The derivatives of the noise susceptibilities by the control
amplitudes. Defaults to None.
amplitudes. The rows and columns should be in the same order
as the corresponding identifiers above. Defaults to None, in
which case the coefficients are assumed to be constant and
hence their derivative vanishing.

.. warning::

Internally, control and noise terms of the Hamiltonian
are stored alphanumerically sorted by their identifiers.
If the noise and/or control identifiers above are not
explicitly given, the rows and/or columns of this
parameter need to be sorted in the same fashion.

Returns
-------
filter_function_deriv: ndarray, shape (n_nops, n_t, n_ctrl, n_omega)
The regular filter functions' derivatives for variation in
each control contribution.
each control contribution. Sorted in the same fashion as
`n_coeffs_deriv` or, if not given, alphanumerically by the
identifiers.

"""
control_matrix = self.get_control_matrix(omega, cache_intermediates=True)
c_idx = util.get_indices_from_identifiers(self.c_oper_identifiers, control_identifiers)
n_idx = util.get_indices_from_identifiers(self.n_oper_identifiers, n_oper_identifiers)

if n_coeffs_deriv is not None:
# TODO 05/22: walrus once support for 3.7 is dropped.
actual_shape = np.shape(n_coeffs_deriv)
required_shape = (len(n_idx), len(c_idx), len(self))
if actual_shape != required_shape:
raise ValueError(f'Expected n_coeffs_deriv to be of shape {required_shape}, '
f'not {actual_shape}. Did you forget to specify identifiers?')
else:
# Do nothing; n_coeffs_deriv specifies the sorting order. If identifiers
# are given, we sort everything else by them. If not, n_coeffs_deriv is
# expected to be sorted in accordance with the opers and coeffs.
pass

# Check if we can pass on intermediates.
intermediates = dict()
# TODO 05/22: walrus once support for 3.7 is dropped.
n_opers_transformed = self._intermediates.get('n_opers_transformed')
first_order_integral = self._intermediates.get('first_order_integral')
if n_opers_transformed is not None:
intermediates['n_opers_transformed'] = n_opers_transformed[n_idx]
if first_order_integral is not None:
intermediates['first_order_integral'] = first_order_integral

control_matrix = self.get_control_matrix(omega, cache_intermediates=True)[n_idx]
control_matrix_deriv = gradient.calculate_derivative_of_control_matrix_from_scratch(
omega, self.propagators, self.eigvals, self.eigvecs, self.basis, self.t, self.dt,
self.n_opers, self.n_coeffs, self.c_opers, self.c_oper_identifiers,
control_identifiers, n_coeffs_deriv, self._intermediates
self.n_opers[n_idx], self.n_coeffs[n_idx], self.c_opers[c_idx], n_coeffs_deriv,
intermediates
)
return gradient.calculate_filter_function_derivative(control_matrix, control_matrix_deriv)

Expand Down
13 changes: 3 additions & 10 deletions tests/gradient_testutil.py
Original file line number Diff line number Diff line change
Expand Up @@ -21,10 +21,6 @@ def deriv_exchange_interaction(eps):
return j_0 / eps_0 * np.exp(eps / eps_0)


def deriv_2_exchange_interaction(eps):
return j_0 / eps_0 ** 2 * np.exp(eps / eps_0)


def one_over_f_noise(f):
spectrum = np.divide(S_0, f, where=(f != 0))
spectrum[f == 0] = spectrum[np.abs(f).argmin()]
Expand All @@ -40,8 +36,6 @@ def create_sing_trip_pulse_seq(eps, dbz, *args):

H_n = [
[sigma_z, deriv_exchange_interaction(eps[0])]
# [sigma_z, np.ones_like(eps[0])]

]

dt = time_step * np.ones(n_time_steps)
Expand Down Expand Up @@ -84,18 +78,17 @@ def finite_diff_infid(u_ctrl_central, u_drift, d, pulse_sequence_builder,

"""
pulse = pulse_sequence_builder(u_ctrl_central, u_drift, d)
all_id = pulse.c_oper_identifiers
if c_id is None:
c_id = all_id[:len(u_ctrl_central)]
c_id = pulse.c_oper_identifiers[:len(u_ctrl_central)]

# Make sure we test for zero frequency case (possible convergence issues)
omega = ff.util.get_sample_frequencies(pulse=pulse, n_samples=n_freq_samples, spacing='log',
include_quasistatic=True)
spectrum = spectral_noise_density(omega)

gradient = np.empty((len(pulse.n_coeffs), len(pulse.dt), len(c_id)))
gradient = np.empty(pulse.n_coeffs.shape + (len(c_id),))

for g in range(len(pulse.dt)):
for g in range(len(pulse)):
for k in range(len(c_id)):
u_plus = u_ctrl_central.copy()
u_plus[k, g] += delta_u
Expand Down
Loading