Compare commits
15 Commits
Author | SHA1 | Date | |
---|---|---|---|
885508e104 | |||
6193ecb9c9 | |||
5ad442750e | |||
9b1538b3c6 | |||
7b277fdc85 | |||
e5fc1207a8 | |||
387a607e09 | |||
9e6d1df559 | |||
45031857f2 | |||
64181eeef2 | |||
d682d50554 | |||
e9e34162a3 | |||
32a7812a43 | |||
7c39475742 | |||
1dd7569fde |
35
CHANGELOG.md
35
CHANGELOG.md
@ -2,6 +2,41 @@
|
||||
|
||||
All notable changes to this project will be documented in this file. See [standard-version](https://github.com/conventional-changelog/standard-version) for commit guidelines.
|
||||
|
||||
## [1.5.0](https://gitea.deepak.science:2222/physics/pdme/compare/1.4.0...1.5.0) (2024-05-17)
|
||||
|
||||
|
||||
### Features
|
||||
|
||||
* adds mcmc chain that returns number of repeats ([6193ecb](https://gitea.deepak.science:2222/physics/pdme/commit/6193ecb9c9f7a21d24e860987a7107549a4b2fa7))
|
||||
|
||||
## [1.4.0](https://gitea.deepak.science:2222/physics/pdme/compare/1.3.0...1.4.0) (2024-05-17)
|
||||
|
||||
|
||||
### Features
|
||||
|
||||
* adds relative squared diff calc utility method ([9b1538b](https://gitea.deepak.science:2222/physics/pdme/commit/9b1538b3c63bfaf2a779bb109cd160a8d7887195))
|
||||
|
||||
## [1.3.0](https://gitea.deepak.science:2222/physics/pdme/compare/1.2.0...1.3.0) (2024-05-17)
|
||||
|
||||
|
||||
### Features
|
||||
|
||||
* adds utility function for sorting samples by frequency for subspace simulation ([e5fc120](https://gitea.deepak.science:2222/physics/pdme/commit/e5fc1207a8b7d5b67208ad825907baa442eec648))
|
||||
|
||||
## [1.2.0](https://gitea.deepak.science:2222/physics/pdme/compare/1.1.0...1.2.0) (2024-05-03)
|
||||
|
||||
|
||||
### Features
|
||||
|
||||
* adds pdme fast calc for e field xs ([9e6d1df](https://gitea.deepak.science:2222/physics/pdme/commit/9e6d1df559e58998851a1c2bf24fcc46d8c1b148))
|
||||
|
||||
## [1.1.0](https://gitea.deepak.science:2222/physics/pdme/compare/1.0.0...1.1.0) (2024-05-02)
|
||||
|
||||
|
||||
### Features
|
||||
|
||||
* adds both electric potential and electric field x sources, makes some fast util tests use the slower explicit versions as double check ([e9e3416](https://gitea.deepak.science:2222/physics/pdme/commit/e9e34162a3b84faad5c18ddeda327c2f7f5ac5aa))
|
||||
|
||||
## [1.0.0](https://gitea.deepak.science:2222/physics/pdme/compare/0.9.3...1.0.0) (2024-04-29)
|
||||
|
||||
|
||||
|
@ -5,7 +5,7 @@
|
||||
[](https://jenkins.deepak.science/job/gitea-physics/job/pdme/job/master/)
|
||||

|
||||

|
||||

|
||||

|
||||
|
||||
This repo has library code for evaluating dipole models.
|
||||
|
||||
|
5029
diagnosis1.nb
5029
diagnosis1.nb
File diff suppressed because it is too large
Load Diff
10
justfile
10
justfile
@ -34,6 +34,16 @@ test: fmt
|
||||
poetry run pytest
|
||||
fi
|
||||
|
||||
# update all test snapshots, use if snapshots are out of date
|
||||
update-snapshots:
|
||||
#!/usr/bin/env bash
|
||||
set -euxo pipefail
|
||||
if [[ "${DO_NIX_CUSTOM:=0}" -eq 1 ]]; then
|
||||
pytest --snapshot-update
|
||||
else
|
||||
poetry run pytest --snapshot-update
|
||||
fi
|
||||
|
||||
# format code
|
||||
fmt:
|
||||
#!/usr/bin/env bash
|
||||
|
46
pdme/calculations/__init__.py
Normal file
46
pdme/calculations/__init__.py
Normal file
@ -0,0 +1,46 @@
|
||||
"""
|
||||
This module is a canonical source of the accurate expressions we want to use for calculating our noise.
|
||||
No reference to class or anything, just a straight set of math functions.
|
||||
"""
|
||||
import numpy
|
||||
|
||||
|
||||
def telegraph_beta(f: float, w: float) -> float:
|
||||
"""
|
||||
This function represents the frequency component of analytic telegraph noise.
|
||||
|
||||
We're assuming we care about the one-sided PSD where we are ignoring negative frequencies.
|
||||
This matches with experimental data from say Connors et al., and I think is better than keeping with one-sided.
|
||||
Note that this means that it will only be comparable then with time series data assuming one-sided!
|
||||
|
||||
Don't bikeshed yet, if we care about two-sided things for any reason down the line divide this by two or just change it then.
|
||||
"""
|
||||
return 2 * w / ((numpy.pi * f) ** 2 + w**2)
|
||||
|
||||
|
||||
def electric_potential(p: numpy.ndarray, s: numpy.ndarray, r: numpy.ndarray) -> float:
|
||||
"""
|
||||
Gives the electric potential of a defect with dipole moment p, located at position s,
|
||||
as measured from position r.
|
||||
|
||||
p, s, r, are numpy arrays of length 3
|
||||
"""
|
||||
diff = r - s
|
||||
return (p.dot(diff) / (numpy.linalg.norm(diff) ** 3)).item()
|
||||
|
||||
|
||||
def electric_field(
|
||||
p: numpy.ndarray, s: numpy.ndarray, r: numpy.ndarray
|
||||
) -> numpy.ndarray:
|
||||
"""
|
||||
Gives the electric field of a defect with dipole moment p, located at position s,
|
||||
as measured from position r.
|
||||
|
||||
p, s, r, are numpy arrays of length 3
|
||||
|
||||
Returns an array of length 3, ideally.
|
||||
"""
|
||||
diff = r - s
|
||||
norm_diff = numpy.linalg.norm(diff)
|
||||
|
||||
return ((3 * (p.dot(diff) * diff) / (norm_diff**2)) - p) / (norm_diff**3)
|
@ -7,6 +7,7 @@ from pdme.measurement.dot_pair_measure import (
|
||||
DotPairMeasurement,
|
||||
DotPairRangeMeasurement,
|
||||
)
|
||||
import pdme.calculations
|
||||
from pdme.measurement.input_types import DotInput, DotPairInput
|
||||
|
||||
|
||||
@ -38,10 +39,27 @@ class OscillatingDipole:
|
||||
self.p = numpy.array(self.p)
|
||||
self.s = numpy.array(self.s)
|
||||
|
||||
def s_at_position(self, r: numpy.ndarray, f: float) -> float:
|
||||
def _alpha_electric_potential(self, r: numpy.ndarray) -> float:
|
||||
"""
|
||||
Returns the electric potential of this dipole at position r.
|
||||
"""
|
||||
return pdme.calculations.electric_potential(self.p, self.s, r)
|
||||
|
||||
def _alpha_electric_field(self, r: numpy.ndarray) -> numpy.ndarray:
|
||||
"""
|
||||
Returns the electric field of this dipole at position r.
|
||||
"""
|
||||
return pdme.calculations.electric_field(self.p, self.s, r)
|
||||
|
||||
def _b(self, f: float) -> float:
|
||||
return pdme.calculations.telegraph_beta(f, self.w)
|
||||
|
||||
def s_electric_potential_at_position(self, r: numpy.ndarray, f: float) -> float:
|
||||
"""
|
||||
Returns the noise potential at a point r, at some frequency f.
|
||||
|
||||
Specifically for electric potential!
|
||||
|
||||
Parameters
|
||||
----------
|
||||
r : numpy.ndarray
|
||||
@ -50,17 +68,49 @@ class OscillatingDipole:
|
||||
f : float
|
||||
The dot frequency to sample.
|
||||
"""
|
||||
return (self._alpha(r)) ** 2 * self._b(f)
|
||||
return (self._alpha_electric_potential(r)) ** 2 * self._b(f)
|
||||
|
||||
def _alpha(self, r: numpy.ndarray) -> float:
|
||||
diff = r - self.s
|
||||
return self.p.dot(diff) / (numpy.linalg.norm(diff) ** 3)
|
||||
def s_electric_potential_for_dot_pair(
|
||||
self, r1: numpy.ndarray, r2: numpy.ndarray, f: float
|
||||
) -> float:
|
||||
"""
|
||||
This is specifically the analytic cpsd for electric potential noise.
|
||||
This should be deprecated
|
||||
"""
|
||||
return (
|
||||
self._alpha_electric_potential(r1)
|
||||
* self._alpha_electric_potential(r2)
|
||||
* self._b(f)
|
||||
)
|
||||
|
||||
def _b(self, f: float) -> float:
|
||||
return (1 / numpy.pi) * (self.w / (f**2 + self.w**2))
|
||||
def s_electric_fieldx_at_position(self, r: numpy.ndarray, f: float) -> float:
|
||||
"""
|
||||
Returns the noise potential at a point r, at some frequency f.
|
||||
|
||||
def s_for_dot_pair(self, r1: numpy.ndarray, r2: numpy.ndarray, f: float) -> float:
|
||||
return self._alpha(r1) * self._alpha(r2) * self._b(f)
|
||||
Specifically for electric potential!
|
||||
|
||||
Parameters
|
||||
----------
|
||||
r : numpy.ndarray
|
||||
The position of the dot.
|
||||
|
||||
f : float
|
||||
The dot frequency to sample.
|
||||
"""
|
||||
return (self._alpha_electric_field(r)[0]) ** 2 * self._b(f)
|
||||
|
||||
def s_electric_fieldx_for_dot_pair(
|
||||
self, r1: numpy.ndarray, r2: numpy.ndarray, f: float
|
||||
) -> float:
|
||||
"""
|
||||
This is specifically the analytic cpsd for electric potential noise.
|
||||
This should be deprecated
|
||||
"""
|
||||
return (
|
||||
self._alpha_electric_field(r1)[0]
|
||||
* self._alpha_electric_field(r2)[0]
|
||||
* self._b(f)
|
||||
)
|
||||
|
||||
def to_flat_array(self) -> numpy.ndarray:
|
||||
return numpy.concatenate([self.p, self.s, numpy.array([self.w])])
|
||||
@ -78,69 +128,97 @@ class OscillatingDipoleArrangement:
|
||||
def __init__(self, dipoles: Sequence[OscillatingDipole]):
|
||||
self.dipoles = dipoles
|
||||
|
||||
def get_dot_measurement(self, dot_input: DotInput) -> DotMeasurement:
|
||||
def get_potential_dot_measurement(self, dot_input: DotInput) -> DotMeasurement:
|
||||
r = numpy.array(dot_input[0])
|
||||
f = dot_input[1]
|
||||
return DotMeasurement(
|
||||
sum([dipole.s_at_position(r, f) for dipole in self.dipoles]), r, f
|
||||
sum(
|
||||
[
|
||||
dipole.s_electric_potential_at_position(r, f)
|
||||
for dipole in self.dipoles
|
||||
]
|
||||
),
|
||||
r,
|
||||
f,
|
||||
)
|
||||
|
||||
def get_dot_pair_measurement(
|
||||
def get_potential_dot_pair_measurement(
|
||||
self, dot_pair_input: DotPairInput
|
||||
) -> DotPairMeasurement:
|
||||
r1 = numpy.array(dot_pair_input[0])
|
||||
r2 = numpy.array(dot_pair_input[1])
|
||||
f = dot_pair_input[2]
|
||||
return DotPairMeasurement(
|
||||
sum([dipole.s_for_dot_pair(r1, r2, f) for dipole in self.dipoles]),
|
||||
sum(
|
||||
[
|
||||
dipole.s_electric_potential_for_dot_pair(r1, r2, f)
|
||||
for dipole in self.dipoles
|
||||
]
|
||||
),
|
||||
r1,
|
||||
r2,
|
||||
f,
|
||||
)
|
||||
|
||||
def get_dot_measurements(
|
||||
def get_potential_dot_measurements(
|
||||
self, dot_inputs: Sequence[DotInput]
|
||||
) -> List[DotMeasurement]:
|
||||
"""
|
||||
For a series of points, each with three coordinates and a frequency, return a list of the corresponding DotMeasurements.
|
||||
"""
|
||||
return [self.get_dot_measurement(dot_input) for dot_input in dot_inputs]
|
||||
return [
|
||||
self.get_potential_dot_measurement(dot_input) for dot_input in dot_inputs
|
||||
]
|
||||
|
||||
def get_dot_pair_measurements(
|
||||
def get_potential_dot_pair_measurements(
|
||||
self, dot_pair_inputs: Sequence[DotPairInput]
|
||||
) -> List[DotPairMeasurement]:
|
||||
"""
|
||||
For a series of pairs of points, each with three coordinates and a frequency, return a list of the corresponding DotPairMeasurements.
|
||||
"""
|
||||
return [
|
||||
self.get_dot_pair_measurement(dot_pair_input)
|
||||
self.get_potential_dot_pair_measurement(dot_pair_input)
|
||||
for dot_pair_input in dot_pair_inputs
|
||||
]
|
||||
|
||||
def get_percent_range_dot_measurement(
|
||||
def get_percent_range_potential_dot_measurement(
|
||||
self, dot_input: DotInput, low_percent: float, high_percent: float
|
||||
) -> DotRangeMeasurement:
|
||||
r = numpy.array(dot_input[0])
|
||||
f = dot_input[1]
|
||||
return DotRangeMeasurement(
|
||||
low_percent * sum([dipole.s_at_position(r, f) for dipole in self.dipoles]),
|
||||
high_percent * sum([dipole.s_at_position(r, f) for dipole in self.dipoles]),
|
||||
low_percent
|
||||
* sum(
|
||||
[
|
||||
dipole.s_electric_potential_at_position(r, f)
|
||||
for dipole in self.dipoles
|
||||
]
|
||||
),
|
||||
high_percent
|
||||
* sum(
|
||||
[
|
||||
dipole.s_electric_potential_at_position(r, f)
|
||||
for dipole in self.dipoles
|
||||
]
|
||||
),
|
||||
r,
|
||||
f,
|
||||
)
|
||||
|
||||
def get_percent_range_dot_measurements(
|
||||
def get_percent_range_potential_dot_measurements(
|
||||
self, dot_inputs: Sequence[DotInput], low_percent: float, high_percent: float
|
||||
) -> List[DotRangeMeasurement]:
|
||||
"""
|
||||
For a series of pairs of points, each with three coordinates and a frequency, and also a lower error range and upper error range, return a list of the corresponding DotPairRangeMeasurements.
|
||||
"""
|
||||
return [
|
||||
self.get_percent_range_dot_measurement(dot_input, low_percent, high_percent)
|
||||
self.get_percent_range_potential_dot_measurement(
|
||||
dot_input, low_percent, high_percent
|
||||
)
|
||||
for dot_input in dot_inputs
|
||||
]
|
||||
|
||||
def get_percent_range_dot_pair_measurement(
|
||||
def get_percent_range_potential_dot_pair_measurement(
|
||||
self, pair_input: DotPairInput, low_percent: float, high_percent: float
|
||||
) -> DotPairRangeMeasurement:
|
||||
r1 = numpy.array(pair_input[0])
|
||||
@ -148,15 +226,25 @@ class OscillatingDipoleArrangement:
|
||||
f = pair_input[2]
|
||||
return DotPairRangeMeasurement(
|
||||
low_percent
|
||||
* sum([dipole.s_for_dot_pair(r1, r2, f) for dipole in self.dipoles]),
|
||||
* sum(
|
||||
[
|
||||
dipole.s_electric_potential_for_dot_pair(r1, r2, f)
|
||||
for dipole in self.dipoles
|
||||
]
|
||||
),
|
||||
high_percent
|
||||
* sum([dipole.s_for_dot_pair(r1, r2, f) for dipole in self.dipoles]),
|
||||
* sum(
|
||||
[
|
||||
dipole.s_electric_potential_for_dot_pair(r1, r2, f)
|
||||
for dipole in self.dipoles
|
||||
]
|
||||
),
|
||||
r1,
|
||||
r2,
|
||||
f,
|
||||
)
|
||||
|
||||
def get_percent_range_dot_pair_measurements(
|
||||
def get_percent_range_potential_dot_pair_measurements(
|
||||
self,
|
||||
pair_inputs: Sequence[DotPairInput],
|
||||
low_percent: float,
|
||||
@ -166,7 +254,7 @@ class OscillatingDipoleArrangement:
|
||||
For a series of pairs of points, each with three coordinates and a frequency, and also a lower error range and upper error range, return a list of the corresponding DotPairRangeMeasurements.
|
||||
"""
|
||||
return [
|
||||
self.get_percent_range_dot_pair_measurement(
|
||||
self.get_percent_range_potential_dot_pair_measurement(
|
||||
pair_input, low_percent, high_percent
|
||||
)
|
||||
for pair_input in pair_inputs
|
||||
|
@ -99,3 +99,65 @@ class DipoleModel:
|
||||
else:
|
||||
chain.append((numpy.squeeze(current_cost).item(), current))
|
||||
return chain
|
||||
|
||||
def get_repeat_counting_mcmc_chain(
|
||||
self,
|
||||
seed,
|
||||
cost_function,
|
||||
chain_length,
|
||||
threshold_cost: float,
|
||||
stdevs: pdme.subspace_simulation.MCMCStandardDeviation,
|
||||
initial_cost: Optional[float] = None,
|
||||
rng_arg: Optional[numpy.random.Generator] = None,
|
||||
) -> Tuple[int, List[Tuple[float, numpy.ndarray]]]:
|
||||
"""
|
||||
performs constrained markov chain monte carlo starting on seed parameter.
|
||||
The cost function given is used as a constrained to condition the chain;
|
||||
a new state is only accepted if cost_function(state) < cost_function(previous_state).
|
||||
The stdevs passed in are the stdevs we're expected to use.
|
||||
|
||||
Because we're using this for subspace simulation where our proposal function is not too important, we're in good shape.
|
||||
Note that for our adaptive stdevs to work, there's an unwritten contract that we sort each dipole in the state by frequency (increasing).
|
||||
|
||||
The seed is a list of dipoles, and each chain state is a list of dipoles as well.
|
||||
|
||||
initial_cost is a performance guy that lets you pre-populate the initial cost used to define the condition.
|
||||
Probably premature optimisation.
|
||||
|
||||
Chain has type of [ (cost: float, state: dipole_ndarray ) ] format,
|
||||
returning (repeat_count, chain) to keep track of number of repeats
|
||||
"""
|
||||
_logger.debug(
|
||||
f"Starting Markov Chain Monte Carlo with seed: {seed} for chain length {chain_length} and provided stdevs {stdevs}"
|
||||
)
|
||||
chain: List[Tuple[float, numpy.ndarray]] = []
|
||||
if initial_cost is None:
|
||||
current_cost = cost_function(numpy.array([seed]))
|
||||
else:
|
||||
current_cost = initial_cost
|
||||
current = seed
|
||||
repeat_event_count = 0
|
||||
for _ in range(chain_length):
|
||||
dips = []
|
||||
for dipole_index, dipole in enumerate(current):
|
||||
_logger.debug(dipole_index)
|
||||
_logger.debug(dipole)
|
||||
stdev = stdevs[dipole_index]
|
||||
tentative_dip = self.markov_chain_monte_carlo_proposal(
|
||||
dipole, stdev, rng_arg
|
||||
)
|
||||
|
||||
dips.append(tentative_dip)
|
||||
dips_array = pdme.subspace_simulation.sort_array_of_dipoles_by_frequency(
|
||||
dips
|
||||
)
|
||||
tentative_cost = cost_function(numpy.array([dips_array]))[0]
|
||||
if tentative_cost < threshold_cost:
|
||||
chain.append((numpy.squeeze(tentative_cost).item(), dips_array))
|
||||
current = dips_array
|
||||
current_cost = tentative_cost
|
||||
else:
|
||||
# repeating a sample, increase count
|
||||
repeat_event_count += 1
|
||||
chain.append((numpy.squeeze(current_cost).item(), current))
|
||||
return (repeat_event_count, chain)
|
||||
|
@ -41,11 +41,30 @@ def sort_array_of_dipoles_by_frequency(configuration) -> numpy.ndarray:
|
||||
Say we have a situation of 2 dipoles, and we've created 8 samples. Then we'll have an (8, 2, 7) numpy array.
|
||||
For each of the 8 samples, we want the 2 dipoles to be in order of frequency.
|
||||
|
||||
This just sorts each sample, the 2x7 array.
|
||||
|
||||
Utility function.
|
||||
"""
|
||||
return numpy.array(sorted(configuration, key=lambda l: l[6]))
|
||||
|
||||
|
||||
def sort_array_of_dipoleses_by_frequency(configurations) -> numpy.ndarray:
|
||||
"""
|
||||
Say we have a situation of 2 dipoles, and we've created 8 samples. Then we'll have an (8, 2, 7) numpy array.
|
||||
For each of the 8 samples, we want the 2 dipoles to be in order of frequency.
|
||||
|
||||
This is the wrapper that sorts everything.
|
||||
|
||||
Utility function.
|
||||
"""
|
||||
return numpy.array(
|
||||
[
|
||||
sort_array_of_dipoles_by_frequency(configuration)
|
||||
for configuration in configurations
|
||||
]
|
||||
)
|
||||
|
||||
|
||||
__all__ = [
|
||||
"DipoleStandardDeviation",
|
||||
"MCMCStandardDeviation",
|
||||
|
@ -18,3 +18,13 @@ def proportional_costs_vs_actual_measurement(
|
||||
dot_inputs_array, dipoles_to_test
|
||||
)
|
||||
return proportional_cost(actual_measurement_array, vals)
|
||||
|
||||
|
||||
def relative_square_diffs(
|
||||
approx: numpy.ndarray, target: numpy.ndarray
|
||||
) -> numpy.ndarray:
|
||||
# Assume that both approx and target are arrays of length m
|
||||
# Approx can broadcast if additional indexes to the left
|
||||
# diffs.shape = [ m ]
|
||||
diffs = (approx - target) ** 2 / (target**2)
|
||||
return diffs.sum(axis=-1)
|
||||
|
@ -45,7 +45,7 @@ def fast_s_nonlocal(
|
||||
_logger.debug(f"alphses1: {alphses1}")
|
||||
_logger.debug(f"alphses2: {alphses2}")
|
||||
|
||||
bses = (1 / numpy.pi) * (ws[:, None] / (f1s**2 + ws[:, None] ** 2))
|
||||
bses = 2 * ws[:, None] / ((numpy.pi * f1s) ** 2 + ws[:, None] ** 2)
|
||||
if _logger.isEnabledFor(logging.DEBUG):
|
||||
_logger.debug(f"bses: {bses}")
|
||||
|
||||
@ -58,44 +58,78 @@ def fast_s_nonlocal_dipoleses(
|
||||
"""
|
||||
No error correction here baby.
|
||||
"""
|
||||
|
||||
# We're going to annotate the indices on this class.
|
||||
# Let's define some indices:
|
||||
# A -> index of dipoleses configurations
|
||||
# measurement_index -> if we have 100 frequencies for example, indexes which one of them it is
|
||||
# j -> within a particular configuration, indexes dipole j
|
||||
# If we need to use numbers, let's use A -> 2, j -> 10, measurement_index -> 9 for consistency with
|
||||
# my other notes
|
||||
# cart -> {x, y, z} is a cartesian axis
|
||||
|
||||
# ps, ss have shape [A, j, cart]
|
||||
ps = dipoleses[:, :, 0:3]
|
||||
ss = dipoleses[:, :, 3:6]
|
||||
# ws shape [A, j]
|
||||
ws = dipoleses[:, :, 6]
|
||||
|
||||
_logger.debug(f"ps: {ps}")
|
||||
_logger.debug(f"ss: {ss}")
|
||||
_logger.debug(f"ws: {ws}")
|
||||
|
||||
# rs have shape [meas_idx, {}, cart], where the inner index goes away leaving
|
||||
# [meas, cart]
|
||||
r1s = dot_pair_inputs[:, 0, 0:3]
|
||||
r2s = dot_pair_inputs[:, 1, 0:3]
|
||||
# fs have index [meas_idx], this makes sense
|
||||
f1s = dot_pair_inputs[:, 0, 3]
|
||||
f2s = dot_pair_inputs[:, 1, 3]
|
||||
|
||||
if (f1s != f2s).all():
|
||||
raise ValueError(f"Dot pair frequencies are inconsistent: {dot_pair_inputs}")
|
||||
|
||||
# r1s have shape [meas, cart], adding the none makes it
|
||||
# r1s[:, None].shape = [meas, 1, cart]
|
||||
# ss[:, None, :].shape = [A, 1, j, cart]
|
||||
# subtracting broadcasts by matching from the right to the left, giving a final shape of
|
||||
# diffses.shape [A, meas, j, cart]
|
||||
diffses1 = r1s[:, None] - ss[:, None, :]
|
||||
diffses2 = r2s[:, None] - ss[:, None, :]
|
||||
if _logger.isEnabledFor(logging.DEBUG):
|
||||
_logger.debug(f"diffses1: {diffses1}")
|
||||
_logger.debug(f"diffses2: {diffses2}")
|
||||
|
||||
# norming on the cartesian axis, which is axis 3 as seen above
|
||||
# norms.shape [A, meas, j]
|
||||
norms1 = numpy.linalg.norm(diffses1, axis=3) ** 3
|
||||
norms2 = numpy.linalg.norm(diffses2, axis=3) ** 3
|
||||
if _logger.isEnabledFor(logging.DEBUG):
|
||||
_logger.debug(f"norms1: {norms1}")
|
||||
_logger.debug(f"norms2: {norms2}")
|
||||
|
||||
# diffses shape [A, meas, j, cart]
|
||||
# ps shape [A, j, cart]
|
||||
# so we're summing over the d axis, the cartesian one.
|
||||
# final shape of numerator is [A, meas, j]
|
||||
# denom shape is [A, meas, j]
|
||||
# final shape stayes [A, meas, j]
|
||||
alphses1 = numpy.einsum("abcd,acd->abc", diffses1, ps) / norms1
|
||||
alphses2 = numpy.einsum("abcd,acd->abc", diffses2, ps) / norms2
|
||||
if _logger.isEnabledFor(logging.DEBUG):
|
||||
_logger.debug(f"alphses1: {alphses1}")
|
||||
_logger.debug(f"alphses2: {alphses2}")
|
||||
|
||||
bses = (1 / numpy.pi) * (ws[:, None, :] / (f1s[:, None] ** 2 + ws[:, None, :] ** 2))
|
||||
# ws shape [A, j], so numerator has shape [A, 1, j]
|
||||
# f1s shape is [meas], so first term of denom is [meas, 1]
|
||||
# ws[:, None, :].shape [A, 1, j] so breadcasting the sum in denom gives
|
||||
# denom.shape [A meas, j]
|
||||
# so now overall shape is [A, meas, j]
|
||||
bses = 2 * ws[:, None, :] / ((numpy.pi * f1s[:, None]) ** 2 + ws[:, None, :] ** 2)
|
||||
if _logger.isEnabledFor(logging.DEBUG):
|
||||
_logger.debug(f"bses: {bses}")
|
||||
|
||||
# so our output shape is [A, meas, j]
|
||||
_logger.debug(f"Raw pair calc: [{alphses1 * alphses2 * bses}]")
|
||||
return numpy.einsum("...j->...", alphses1 * alphses2 * bses)
|
||||
|
||||
@ -152,15 +186,15 @@ def fast_s_spin_qubit_tarucha_nonlocal_dipoleses(
|
||||
diffses1 = r1s[:, None] - ss[:, None, :]
|
||||
diffses2 = r2s[:, None] - ss[:, None, :]
|
||||
if _logger.isEnabledFor(logging.DEBUG):
|
||||
_logger.warning(f"diffses1: {diffses1}")
|
||||
_logger.warning(f"diffses2: {diffses2}")
|
||||
_logger.debug(f"diffses1: {diffses1}")
|
||||
_logger.debug(f"diffses2: {diffses2}")
|
||||
|
||||
# norms takes out axis 3, the last one, giving [A, measurement_idx, j]
|
||||
norms1 = numpy.linalg.norm(diffses1, axis=3)
|
||||
norms2 = numpy.linalg.norm(diffses2, axis=3)
|
||||
if _logger.isEnabledFor(logging.DEBUG):
|
||||
_logger.warning(f"norms1: {norms1}")
|
||||
_logger.warning(f"norms2: {norms2}")
|
||||
_logger.debug(f"norms1: {norms1}")
|
||||
_logger.debug(f"norms2: {norms2}")
|
||||
|
||||
# _logger.info(f"norms1: {norms1}")
|
||||
# _logger.info(f"norms1 shape: {norms1.shape}")
|
||||
@ -215,14 +249,19 @@ def fast_s_spin_qubit_tarucha_nonlocal_dipoleses(
|
||||
- ps[:, numpy.newaxis, :, 0]
|
||||
) / (norms2**3)
|
||||
if _logger.isEnabledFor(logging.DEBUG):
|
||||
_logger.warning(f"alphses1: {alphses1}")
|
||||
_logger.warning(f"alphses2: {alphses2}")
|
||||
_logger.debug(f"alphses1: {alphses1}")
|
||||
_logger.debug(f"alphses2: {alphses2}")
|
||||
|
||||
bses = (1 / numpy.pi) * (ws[:, None, :] / (f1s[:, None] ** 2 + ws[:, None, :] ** 2))
|
||||
# ws has shape (A, j), so it becomes (A, 1, j) in numerator with the new axis
|
||||
# f1s has shape (m), so we get in the denominator (m, 1) + (A, 1, j)
|
||||
# This becomes (A, m, j)
|
||||
bses = 2 * ws[:, None, :] / ((numpy.pi * f1s[:, None]) ** 2 + ws[:, None, :] ** 2)
|
||||
if _logger.isEnabledFor(logging.DEBUG):
|
||||
_logger.warning(f"bses: {bses}")
|
||||
_logger.debug(f"bses: {bses}")
|
||||
|
||||
_logger.warning(f"Raw pair calc: [{alphses1 * alphses2 * bses}]")
|
||||
# alphas have (A, 1, j), betas have (A, m, j)
|
||||
# Final result is (A, m, j)
|
||||
_logger.debug(f"Raw pair calc: [{alphses1 * alphses2 * bses}]")
|
||||
return numpy.einsum("...j->...", alphses1 * alphses2 * bses)
|
||||
|
||||
|
||||
|
@ -12,27 +12,44 @@ def fast_vs_for_dipoles(
|
||||
No error correction here baby.
|
||||
Expects dot_inputs to be numpy array of [rx, ry, rz, f] entries, so a n by 4 where n is number of measurement points.
|
||||
"""
|
||||
# name indexes:
|
||||
# A: dipole config index
|
||||
# cart: cartesian index
|
||||
# m: measurement index
|
||||
|
||||
# [A, cart]
|
||||
ps = dipoles[:, 0:3]
|
||||
ss = dipoles[:, 3:6]
|
||||
# [A]
|
||||
ws = dipoles[:, 6]
|
||||
|
||||
_logger.debug(f"ps: {ps}")
|
||||
_logger.debug(f"ss: {ss}")
|
||||
_logger.debug(f"ws: {ws}")
|
||||
|
||||
# [m, cart]
|
||||
rs = dot_inputs[:, 0:3]
|
||||
# [m]
|
||||
fs = dot_inputs[:, 3]
|
||||
|
||||
# [m, cart] - [A, 1, cart]
|
||||
# [A, m, cart]
|
||||
diffses = rs - ss[:, None]
|
||||
|
||||
_logger.debug(f"diffses: {diffses}")
|
||||
# [A, m]
|
||||
norms = numpy.linalg.norm(diffses, axis=2) ** 3
|
||||
_logger.debug(f"norms: {norms}")
|
||||
# [A, m, cart] [A, cart] -> [A, m]
|
||||
ases = (numpy.einsum("...ji, ...i", diffses, ps) / norms) ** 2
|
||||
_logger.debug(f"ases: {ases}")
|
||||
|
||||
bses = (1 / numpy.pi) * (ws[:, None] / (fs**2 + ws[:, None] ** 2))
|
||||
# [A, 1], denom [m] + [A, 1] -> [A, m]
|
||||
# [A, m]
|
||||
bses = 2 * ws[:, None] / ((numpy.pi * fs) ** 2 + ws[:, None] ** 2)
|
||||
|
||||
_logger.debug(f"bses: {bses}")
|
||||
# returns shape [A, m]
|
||||
return ases * bses
|
||||
|
||||
|
||||
@ -64,11 +81,126 @@ def fast_vs_for_dipoleses(
|
||||
ases = (numpy.einsum("abcd,acd->abc", diffses, ps) / norms) ** 2
|
||||
_logger.debug(f"ases: {ases}")
|
||||
|
||||
bses = (1 / numpy.pi) * (ws[:, None, :] / (fs[:, None] ** 2 + ws[:, None, :] ** 2))
|
||||
bses = 2 * ws[:, None, :] / ((numpy.pi * fs[:, None]) ** 2 + ws[:, None, :] ** 2)
|
||||
_logger.debug(f"bses: {bses}")
|
||||
return numpy.einsum("...j->...", ases * bses)
|
||||
|
||||
|
||||
def fast_efieldxs_for_dipoles(
|
||||
dot_inputs: numpy.ndarray, dipoles: numpy.ndarray
|
||||
) -> numpy.ndarray:
|
||||
"""
|
||||
No error correction here baby.
|
||||
Expects dot_inputs to be numpy array of [rx, ry, rz, f] entries, so a n by 4 where n is number of measurement points.
|
||||
"""
|
||||
# name indexes:
|
||||
# A: dipole config index
|
||||
# j: dipole index within a config
|
||||
# cart: cartesian index
|
||||
# m: measurement index
|
||||
|
||||
# Indexes [A, cart]
|
||||
ps = dipoles[:, 0:3]
|
||||
ss = dipoles[:, 3:6]
|
||||
# Indexes [A]
|
||||
ws = dipoles[:, 6]
|
||||
|
||||
# Indexes [m, cart]
|
||||
rs = dot_inputs[:, 0:3]
|
||||
# Indexes [m]
|
||||
fs = dot_inputs[:, 3]
|
||||
|
||||
# Indexes [m, cart] - [A, 1, cart]
|
||||
# Broadcasting from right
|
||||
# diffses.indexes [A, m, cart]
|
||||
diffses = rs - ss[:, None, :]
|
||||
|
||||
# [A, m, cart][2] = cart
|
||||
# norms takes out axis 2, the last one, giving [A, m]
|
||||
norms = numpy.linalg.norm(diffses, axis=2)
|
||||
|
||||
# long story but this ends up becoming (A, 1, j)
|
||||
# Some evidence is looking at ps term, which has shape (A, 1, j, cart=0) becoming (A, 1, j)
|
||||
|
||||
# [A, m, cart] einsum [A, cart] explicitly gives [A, m]
|
||||
dot_products = numpy.einsum("amc,ac->am", diffses, ps) / (norms**2)
|
||||
|
||||
# [m, A] * [cart, m, A] -> [cart, m, A], transpose that and you get [A, m, cart]
|
||||
projections = numpy.transpose(
|
||||
numpy.transpose(dot_products) * numpy.transpose(diffses)
|
||||
)
|
||||
|
||||
# numerator [A, m, cart] - [A, 1, cart] -> [A, m, cart][:, :, 0] -> [A, m]
|
||||
alphas = (3 * projections - ps[:, numpy.newaxis])[:, :, 0] / norms**3
|
||||
|
||||
# [A, m]
|
||||
ases = alphas**2
|
||||
|
||||
# [A, 1], denom [m] + [A, 1] -> [A, m]
|
||||
# [A, m]
|
||||
bses = 2 * ws[:, None] / ((numpy.pi * fs) ** 2 + ws[:, None] ** 2)
|
||||
|
||||
# return shape [A, m, j]
|
||||
return ases * bses
|
||||
|
||||
|
||||
def fast_efieldxs_for_dipoleses(
|
||||
dot_inputs: numpy.ndarray, dipoleses: numpy.ndarray
|
||||
) -> numpy.ndarray:
|
||||
"""
|
||||
No error correction here baby.
|
||||
Expects dot_inputs to be numpy array of [rx, ry, rz, f] entries, so a n by 4 where n is number of measurement points.
|
||||
|
||||
Dipoleses are expected to be array of arrays of arrays:
|
||||
list of sets of dipoles which are part of a single arrangement to be added together.
|
||||
"""
|
||||
# name indexes:
|
||||
# A: dipole config index
|
||||
# j: dipole index within a config
|
||||
# cart: cartesian index
|
||||
# m: measurement index
|
||||
|
||||
# Indexes [A, j, cart]
|
||||
ps = dipoleses[:, :, 0:3]
|
||||
ss = dipoleses[:, :, 3:6]
|
||||
# Indexes [A, j]
|
||||
ws = dipoleses[:, :, 6]
|
||||
|
||||
# Indexes [m, cart]
|
||||
rs = dot_inputs[:, 0:3]
|
||||
# Indexes [m]
|
||||
fs = dot_inputs[:, 3]
|
||||
|
||||
# Indexes [m, 1, cart] - [A, 1, j, cart]
|
||||
# Broadcasting from right
|
||||
# diffses.indexes [A, m, j, cart]
|
||||
diffses = rs[:, None] - ss[:, None, :]
|
||||
|
||||
# norms takes out axis 3, the last one, giving [A, m, j]
|
||||
norms = numpy.linalg.norm(diffses, axis=3)
|
||||
|
||||
# long story but this ends up becoming (A, 1, j)
|
||||
# Some evidence is looking at ps term, which has shape (A, 1, j, cart=0) becoming (A, 1, j)
|
||||
alphas = (
|
||||
(
|
||||
3
|
||||
* numpy.transpose(
|
||||
numpy.transpose(
|
||||
numpy.einsum("abcd,acd->abc", diffses, ps) / (norms**2)
|
||||
)
|
||||
* numpy.transpose(diffses)
|
||||
)[:, :, :, 0]
|
||||
)
|
||||
- ps[:, numpy.newaxis, :, 0]
|
||||
) / (norms**3)
|
||||
ases = alphas**2
|
||||
# bses.shape [A, m, j)]
|
||||
bses = 2 * ws[:, None, :] / ((numpy.pi * fs[:, None]) ** 2 + ws[:, None, :] ** 2)
|
||||
|
||||
# return shape [A, m, j]
|
||||
return numpy.einsum("...j->...", ases * bses)
|
||||
|
||||
|
||||
def fast_vs_for_asymmetric_dipoleses(
|
||||
dot_inputs: numpy.ndarray, dipoleses: numpy.ndarray, temp: numpy.ndarray
|
||||
) -> numpy.ndarray:
|
||||
@ -92,6 +224,9 @@ def fast_vs_for_asymmetric_dipoleses(
|
||||
|
||||
diffses = rs[:, None] - ss[:, None, :]
|
||||
|
||||
_logger.warning(
|
||||
"This method is very likely to be broken, and should not be used without more thought"
|
||||
)
|
||||
w1s = numpy.exp(-e1s / temp) * raw_ws
|
||||
w2s = numpy.exp(-e2s / temp) * raw_ws
|
||||
|
||||
@ -105,7 +240,7 @@ def fast_vs_for_asymmetric_dipoleses(
|
||||
|
||||
ases = (numpy.einsum("abcd,acd->abc", diffses, ps) / norms) ** 2
|
||||
|
||||
bses = (1 / numpy.pi) * (ws[:, None, :] / (fs[:, None] ** 2 + ws[:, None, :] ** 2))
|
||||
bses = ws[:, None, :] / ((numpy.pi * fs[:, None]) ** 2 + ws[:, None, :] ** 2)
|
||||
|
||||
return numpy.einsum("...j->...", ases * bses)
|
||||
|
||||
|
@ -1,6 +1,6 @@
|
||||
[tool.poetry]
|
||||
name = "pdme"
|
||||
version = "1.0.0"
|
||||
version = "1.5.0"
|
||||
description = "Python dipole model evaluator"
|
||||
authors = ["Deepak <dmallubhotla+github@gmail.com>"]
|
||||
license = "GPL-3.0-only"
|
||||
|
81
tests/calculations/__snapshots__/test_calculations.ambr
Normal file
81
tests/calculations/__snapshots__/test_calculations.ambr
Normal file
@ -0,0 +1,81 @@
|
||||
# serializer version: 1
|
||||
# name: test_multiple_electric_field_alphas
|
||||
list([
|
||||
list([
|
||||
0.009579434215333742,
|
||||
0.007714411624737791,
|
||||
0.0035604976729559034,
|
||||
]),
|
||||
list([
|
||||
0.005036717012002481,
|
||||
0.006259221141129298,
|
||||
-0.009682232702684382,
|
||||
]),
|
||||
list([
|
||||
0.01503516326014651,
|
||||
0.012028130608117207,
|
||||
0.009021097956087907,
|
||||
]),
|
||||
list([
|
||||
0.0033871215792458027,
|
||||
0.0038182097802407235,
|
||||
-0.005604146612933966,
|
||||
]),
|
||||
list([
|
||||
0.007409666906433089,
|
||||
0.006778734778841714,
|
||||
0.003467602973242188,
|
||||
]),
|
||||
list([
|
||||
0.004730939083250055,
|
||||
0.0046546336141653774,
|
||||
-0.011766303332857395,
|
||||
]),
|
||||
list([
|
||||
-0.010766266772985707,
|
||||
-0.012631289363581657,
|
||||
-0.010851040527103707,
|
||||
]),
|
||||
list([
|
||||
-0.008410828408392494,
|
||||
-0.020391368873835285,
|
||||
0.009877833363344673,
|
||||
]),
|
||||
list([
|
||||
-0.015035163260146511,
|
||||
-0.018042195912175815,
|
||||
-0.021049228564205113,
|
||||
]),
|
||||
list([
|
||||
-0.005850482727788205,
|
||||
-0.012193637685284892,
|
||||
0.0054809785555068454,
|
||||
]),
|
||||
list([
|
||||
-0.007682229585552562,
|
||||
-0.010584517372472879,
|
||||
-0.009352937859414517,
|
||||
]),
|
||||
list([
|
||||
-0.009767100042838822,
|
||||
-0.020831393060117175,
|
||||
0.014024945217763873,
|
||||
]),
|
||||
])
|
||||
# ---
|
||||
# name: test_multiple_electric_potential_alphas
|
||||
list([
|
||||
0.05035560994609065,
|
||||
-0.03369221379873504,
|
||||
0.07216878364870323,
|
||||
-0.024633611485424027,
|
||||
0.04496779459769236,
|
||||
-0.03108684810509794,
|
||||
-0.08019597139562586,
|
||||
0.08293468011996319,
|
||||
-0.10825317547305485,
|
||||
0.060044427995721066,
|
||||
-0.06935710692186449,
|
||||
0.08733923991432278,
|
||||
])
|
||||
# ---
|
102
tests/calculations/test_calculations.py
Normal file
102
tests/calculations/test_calculations.py
Normal file
@ -0,0 +1,102 @@
|
||||
import pdme.calculations
|
||||
import pytest
|
||||
import numpy
|
||||
import numpy.testing
|
||||
|
||||
|
||||
# generated in mathematica to compare here
|
||||
beta_test_data = [
|
||||
[-2, -2, 0.8072976151],
|
||||
[-1, -2, 0.008105366193],
|
||||
[0, -2, 0.00008105691406],
|
||||
[1, -2, 8.105694659e-7],
|
||||
[2, -2, 8.105694691e-9],
|
||||
[-2, -1, 5.768008783],
|
||||
[-1, -1, 0.08072976151],
|
||||
[0, -1, 0.0008105366193],
|
||||
[1, -1, 8.105691406e-6],
|
||||
[2, -1, 8.105694659e-8],
|
||||
[-2, 0, 1.951840272],
|
||||
[-1, 0, 0.5768008783],
|
||||
[0, 0, 0.008072976151],
|
||||
[1, 0, 0.00008105366193],
|
||||
[2, 0, 8.105691406e-7],
|
||||
[-2, 1, 0.1999506642],
|
||||
[-1, 1, 0.1951840272],
|
||||
[0, 1, 0.05768008783],
|
||||
[1, 1, 0.0008072976151],
|
||||
[2, 1, 8.105366193e-6],
|
||||
[-2, 2, 0.01999995065],
|
||||
[-1, 2, 0.01999506642],
|
||||
[0, 2, 0.01951840272],
|
||||
[1, 2, 0.005768008783],
|
||||
[2, 2, 0.00008072976151],
|
||||
]
|
||||
|
||||
|
||||
@pytest.mark.parametrize("f, w, expected", beta_test_data)
|
||||
def test_calculations_beta_func(f, w, expected):
|
||||
# there's nothing special about the 5 * 10^ f and 10^w passing in logs
|
||||
# this was just to get a variety of orders of magnitude in results.
|
||||
actual = pdme.calculations.telegraph_beta(5 * 10**f, 10**w)
|
||||
numpy.testing.assert_allclose(actual, expected, atol=0, rtol=1e-8)
|
||||
|
||||
|
||||
def test_multiple_electric_potential_alphas(snapshot):
|
||||
"""
|
||||
Manually compare these with mathematica stuff because manually including a list is a bit annoying.
|
||||
Basically just visually compare the snapshot values to make sure they're actually correct.
|
||||
"""
|
||||
dipole_ps = [
|
||||
numpy.array([1, 2, 3]),
|
||||
numpy.array([-4, -3, -2]),
|
||||
]
|
||||
dipole_ss = [
|
||||
numpy.array([0, 0, 1]),
|
||||
numpy.array([1, 1, 1]),
|
||||
numpy.array([0, -0.5, 0.5]),
|
||||
]
|
||||
|
||||
test_rs = [
|
||||
numpy.array([5, 5, 5]),
|
||||
numpy.array([-4, -6, 2]),
|
||||
]
|
||||
|
||||
actuals = [
|
||||
pdme.calculations.electric_potential(p, s, r)
|
||||
for p in dipole_ps
|
||||
for s in dipole_ss
|
||||
for r in test_rs
|
||||
]
|
||||
|
||||
assert actuals == snapshot
|
||||
|
||||
|
||||
def test_multiple_electric_field_alphas(snapshot):
|
||||
"""
|
||||
Manually compare these with mathematica stuff because manually including a list is a bit annoying.
|
||||
Basically just visually compare the snapshot values to make sure they're actually correct.
|
||||
"""
|
||||
dipole_ps = [
|
||||
numpy.array([1, 2, 3]),
|
||||
numpy.array([-4, -3, -2]),
|
||||
]
|
||||
dipole_ss = [
|
||||
numpy.array([0, 0, 1]),
|
||||
numpy.array([1, 1, 1]),
|
||||
numpy.array([0, -0.5, 0.5]),
|
||||
]
|
||||
|
||||
test_rs = [
|
||||
numpy.array([5, 5, 5]),
|
||||
numpy.array([-4, -6, 2]),
|
||||
]
|
||||
|
||||
actuals = [
|
||||
pdme.calculations.electric_field(p, s, r).tolist()
|
||||
for p in dipole_ps
|
||||
for s in dipole_ss
|
||||
for r in test_rs
|
||||
]
|
||||
|
||||
assert actuals == snapshot
|
@ -2,23 +2,30 @@ import numpy
|
||||
import pdme.measurement
|
||||
|
||||
|
||||
def test_static_dipole():
|
||||
d1 = pdme.measurement.OscillatingDipole((1, 2, 3), (4, 5, 6), 7)
|
||||
d2 = pdme.measurement.OscillatingDipole((2, 5, 3), (4, -5, -6), 2)
|
||||
def test_static_dipole_electric_potential():
|
||||
d1 = pdme.measurement.OscillatingDipole(
|
||||
numpy.array((1, 2, 3)), numpy.array([4, 5, 6]), 7
|
||||
)
|
||||
d2 = pdme.measurement.OscillatingDipole(
|
||||
numpy.array([2, 5, 3]), numpy.array([4, -5, -6]), 2
|
||||
)
|
||||
dipoles = pdme.measurement.OscillatingDipoleArrangement([d1, d2])
|
||||
|
||||
dot_position1 = (-1, -1, -1)
|
||||
dot_position1 = numpy.array([-1, -1, -1])
|
||||
dot_frequency1 = 11
|
||||
expected_v1 = 0.00001421963287022476
|
||||
expected_v2 = 0.00001107180225755457
|
||||
|
||||
expected_v1 = 0.00001221710876727626
|
||||
expected_v2 = 7.257229625870065e-6
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
d1.s_at_position(dot_position1, dot_frequency1),
|
||||
d1.s_electric_potential_at_position(dot_position1, dot_frequency1),
|
||||
expected_v1,
|
||||
err_msg="Voltage at dot isn't as expected.",
|
||||
)
|
||||
|
||||
dot_measurements = dipoles.get_dot_measurements([(dot_position1, dot_frequency1)])
|
||||
dot_measurements = dipoles.get_potential_dot_measurements(
|
||||
[(dot_position1, dot_frequency1)]
|
||||
)
|
||||
assert len(dot_measurements) == 1, "Should have only had one dot measurement."
|
||||
measurement = dot_measurements[0]
|
||||
numpy.testing.assert_allclose(
|
||||
@ -38,24 +45,57 @@ def test_static_dipole():
|
||||
)
|
||||
|
||||
|
||||
def test_dipole_dot_pair():
|
||||
d1 = pdme.measurement.OscillatingDipole((1, 2, 3), (4, 5, 6), 7)
|
||||
dipoles = pdme.measurement.OscillatingDipoleArrangement([d1])
|
||||
def test_static_dipole_electric_field_x():
|
||||
d1 = pdme.measurement.OscillatingDipole(
|
||||
numpy.array((1, 2, 3)), numpy.array([4, 5, 6]), 7
|
||||
)
|
||||
d2 = pdme.measurement.OscillatingDipole(
|
||||
numpy.array([2, 5, 3]), numpy.array([4, -5, -6]), 2
|
||||
)
|
||||
# dipoles = pdme.measurement.OscillatingDipoleArrangement([d1, d2])
|
||||
|
||||
dot_position1 = (-1, -1, -1)
|
||||
dot_position2 = (1, 2, 3)
|
||||
dot_frequency = 8
|
||||
expected_sij = 0.000083328037100902801698
|
||||
dot_position1 = numpy.array([-1, -1, -1])
|
||||
dot_frequency1 = 11
|
||||
|
||||
# these should be true for electric field
|
||||
expected_v1 = 1.479556451978925e-7
|
||||
expected_v2 = 6.852024308908262e-7
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
d1.s_for_dot_pair(dot_position1, dot_position2, dot_frequency),
|
||||
d1.s_electric_fieldx_at_position(dot_position1, dot_frequency1),
|
||||
expected_v1,
|
||||
err_msg="Fieldx noise at dot isn't as expected.",
|
||||
)
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
d2.s_electric_fieldx_at_position(dot_position1, dot_frequency1),
|
||||
expected_v2,
|
||||
err_msg="Fieldx at dot isn't as expected.",
|
||||
)
|
||||
|
||||
|
||||
def test_dipole_dot_pair():
|
||||
d1 = pdme.measurement.OscillatingDipole(
|
||||
numpy.array([1, 2, 3]), numpy.array([4, 5, 6]), 7
|
||||
)
|
||||
dipoles = pdme.measurement.OscillatingDipoleArrangement([d1])
|
||||
|
||||
dot_position1 = numpy.array([-1, -1, -1])
|
||||
dot_position2 = numpy.array([1, 2, 3])
|
||||
dot_frequency = 8
|
||||
expected_sij = 0.00008692058236162933
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
d1.s_electric_potential_for_dot_pair(
|
||||
dot_position1, dot_position2, dot_frequency
|
||||
),
|
||||
expected_sij,
|
||||
err_msg="Sij for dot pair isn't as expected.",
|
||||
)
|
||||
numpy.testing.assert_allclose(
|
||||
[
|
||||
m.v
|
||||
for m in dipoles.get_dot_pair_measurements(
|
||||
for m in dipoles.get_potential_dot_pair_measurements(
|
||||
[(dot_position1, dot_position2, dot_frequency)]
|
||||
)
|
||||
],
|
||||
@ -65,15 +105,17 @@ def test_dipole_dot_pair():
|
||||
|
||||
|
||||
def test_range_pairs():
|
||||
d1 = pdme.measurement.OscillatingDipole((1, 2, 3), (4, 5, 6), 7)
|
||||
d1 = pdme.measurement.OscillatingDipole(
|
||||
numpy.array([1, 2, 3]), numpy.array([4, 5, 6]), 7
|
||||
)
|
||||
dipoles = pdme.measurement.OscillatingDipoleArrangement([d1])
|
||||
|
||||
dot_position1 = (-1, -1, -1)
|
||||
dot_position2 = (1, 2, 3)
|
||||
dot_position1 = numpy.array([-1, -1, -1])
|
||||
dot_position2 = numpy.array([1, 2, 3])
|
||||
dot_frequency = 8
|
||||
expected_sij = 0.000083328037100902801698
|
||||
expected_sij = 0.00008692058236162933
|
||||
|
||||
actuals = dipoles.get_percent_range_dot_pair_measurements(
|
||||
actuals = dipoles.get_percent_range_potential_dot_pair_measurements(
|
||||
[(dot_position1, dot_position2, dot_frequency)], 0.5, 1.5
|
||||
)
|
||||
assert len(actuals) == 1, "should have only been one pair"
|
||||
@ -86,22 +128,26 @@ def test_range_pairs():
|
||||
|
||||
|
||||
def test_range_dipole_measurements():
|
||||
d1 = pdme.measurement.OscillatingDipole((1, 2, 3), (4, 5, 6), 7)
|
||||
d2 = pdme.measurement.OscillatingDipole((2, 5, 3), (4, -5, -6), 2)
|
||||
d1 = pdme.measurement.OscillatingDipole(
|
||||
numpy.array([1, 2, 3]), numpy.array([4, 5, 6]), 7
|
||||
)
|
||||
d2 = pdme.measurement.OscillatingDipole(
|
||||
numpy.array([2, 5, 3]), numpy.array([4, -5, -6]), 2
|
||||
)
|
||||
dipoles = pdme.measurement.OscillatingDipoleArrangement([d1, d2])
|
||||
|
||||
dot_position1 = (-1, -1, -1)
|
||||
dot_position1 = numpy.array([-1, -1, -1])
|
||||
dot_frequency1 = 11
|
||||
expected_v1 = 0.00001421963287022476
|
||||
expected_v2 = 0.00001107180225755457
|
||||
expected_v1 = 0.00001221710876727626
|
||||
expected_v2 = 7.257229625870065e-6
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
d1.s_at_position(dot_position1, dot_frequency1),
|
||||
d1.s_electric_potential_at_position(dot_position1, dot_frequency1),
|
||||
expected_v1,
|
||||
err_msg="Voltage at dot isn't as expected.",
|
||||
)
|
||||
|
||||
range_dot_measurements = dipoles.get_percent_range_dot_measurements(
|
||||
range_dot_measurements = dipoles.get_percent_range_potential_dot_measurements(
|
||||
[(dot_position1, dot_frequency1)], 0.5, 1.5
|
||||
)
|
||||
assert len(range_dot_measurements) == 1, "Should have only had one dot measurement."
|
||||
|
@ -2,52 +2,52 @@
|
||||
# name: test_log_spaced_fixedxy_orientation_mcmc_basic
|
||||
list([
|
||||
tuple(
|
||||
3984.461796565,
|
||||
3984.4199552664,
|
||||
array([[ 9.55610128, 2.94634152, 0. , 9.21529051, -2.46576127,
|
||||
2.42481096, 9.19034554]]),
|
||||
),
|
||||
tuple(
|
||||
8583.9908787152,
|
||||
8583.8032728084,
|
||||
array([[ 9.99991539, 0.04113671, 0. , 8.71258954, -2.26599865,
|
||||
2.60452102, 6.37042214]]),
|
||||
),
|
||||
tuple(
|
||||
6215.6376616016,
|
||||
6215.5671525452,
|
||||
array([[ 9.81950685, -1.89137124, 0. , 8.90637055, -2.48043039,
|
||||
2.28444435, 8.84239221]]),
|
||||
),
|
||||
tuple(
|
||||
424.7332846598,
|
||||
424.7282747232,
|
||||
array([[ 1.00028483, 9.94984574, 0. , 8.53064898, -2.59230757,
|
||||
2.33774773, 8.6714416 ]]),
|
||||
),
|
||||
tuple(
|
||||
300.9220380849,
|
||||
300.920463669,
|
||||
array([[ 1.4003442 , 9.90146636, 0. , 8.05557992, -2.6753126 ,
|
||||
2.65915755, 13.02021385]]),
|
||||
),
|
||||
tuple(
|
||||
2400.0107277085,
|
||||
2400.0053792245,
|
||||
array([[ 9.97761813, 0.66868263, 0. , 8.69171028, -2.73145011,
|
||||
2.90140456, 19.94999593]]),
|
||||
),
|
||||
tuple(
|
||||
5001.4620511303,
|
||||
5001.4154377966,
|
||||
array([[ 9.93976109, -1.09596962, 0. , 8.95245025, -2.59409162,
|
||||
2.90140456, 9.75535945]]),
|
||||
),
|
||||
tuple(
|
||||
195.2198074488,
|
||||
195.2191433934,
|
||||
array([[ 0.20690762, 9.99785923, 0. , 9.59636585, -2.83240984,
|
||||
2.90140456, 16.14771567]]),
|
||||
),
|
||||
tuple(
|
||||
2698.258844498,
|
||||
2698.2459654601,
|
||||
array([[-9.68130127, -2.50447712, 0. , 8.94823619, -2.92889659,
|
||||
2.77065328, 13.63173263]]),
|
||||
),
|
||||
tuple(
|
||||
1193.6985473944,
|
||||
1193.6653292625,
|
||||
array([[-6.16597091, -7.87278875, 0. , 9.62210721, -2.75993924,
|
||||
2.77065328, 5.64553534]]),
|
||||
),
|
||||
|
@ -2,52 +2,52 @@
|
||||
# name: test_log_spaced_free_orientation_mcmc_basic
|
||||
list([
|
||||
tuple(
|
||||
3167.6711268743,
|
||||
3167.6603875494,
|
||||
array([[ 9.60483896, -1.41627817, -2.3960853 , -4.76615152, -1.80902942,
|
||||
2.11809123, 16.17452242]]),
|
||||
),
|
||||
tuple(
|
||||
3167.6711268743,
|
||||
3167.6603875494,
|
||||
array([[ 9.60483896, -1.41627817, -2.3960853 , -4.76615152, -1.80902942,
|
||||
2.11809123, 16.17452242]]),
|
||||
),
|
||||
tuple(
|
||||
3167.6711268743,
|
||||
3167.6603875494,
|
||||
array([[ 9.60483896, -1.41627817, -2.3960853 , -4.76615152, -1.80902942,
|
||||
2.11809123, 16.17452242]]),
|
||||
),
|
||||
tuple(
|
||||
736.0306527136,
|
||||
2320.6090682509,
|
||||
array([[ 4.1660069 , -8.11557337, 4.0965663 , -4.35968351, -1.97945216,
|
||||
2.43615641, 12.92143144]]),
|
||||
),
|
||||
tuple(
|
||||
736.0306527136,
|
||||
2320.6090682509,
|
||||
array([[ 4.1660069 , -8.11557337, 4.0965663 , -4.35968351, -1.97945216,
|
||||
2.43615641, 12.92143144]]),
|
||||
),
|
||||
tuple(
|
||||
736.0306527136,
|
||||
2320.6090682509,
|
||||
array([[ 4.1660069 , -8.11557337, 4.0965663 , -4.35968351, -1.97945216,
|
||||
2.43615641, 12.92143144]]),
|
||||
),
|
||||
tuple(
|
||||
2248.0779986277,
|
||||
2248.0760851141,
|
||||
array([[-1.71755535, -5.59925137, 8.10545419, -4.03306318, -1.81098441,
|
||||
2.77407111, 32.28020575]]),
|
||||
),
|
||||
tuple(
|
||||
1663.310672736,
|
||||
1663.3101472167,
|
||||
array([[-5.16785855, 2.7558756 , 8.10545419, -3.34620897, -1.74763642,
|
||||
2.42770463, 52.98214008]]),
|
||||
),
|
||||
tuple(
|
||||
1329.2704143918,
|
||||
1329.2703365134,
|
||||
array([[ -1.39600464, 9.69718343, -2.00394725, -2.59147366,
|
||||
-1.91246681, 2.07361175, 123.01833742]]),
|
||||
),
|
||||
tuple(
|
||||
355.7695591897,
|
||||
355.7695549121,
|
||||
array([[ 9.76047401, 0.84696075, -2.00394725, -3.04310053,
|
||||
-1.99338573, 2.1185589 , 271.35743739]]),
|
||||
),
|
||||
|
@ -2,52 +2,52 @@
|
||||
# name: test_log_spaced_fixed_orientation_mcmc_basic
|
||||
list([
|
||||
tuple(
|
||||
50.5683119299,
|
||||
50.5632126772,
|
||||
array([[ 0. , 0. , 10. , -2.3960853 , 4.23246234,
|
||||
2.26169242, 39.39900844]]),
|
||||
),
|
||||
tuple(
|
||||
50.5683119299,
|
||||
50.5632126772,
|
||||
array([[ 0. , 0. , 10. , -2.3960853 , 4.23246234,
|
||||
2.26169242, 39.39900844]]),
|
||||
),
|
||||
tuple(
|
||||
47.408654554,
|
||||
47.4038652151,
|
||||
array([[ 0. , 0. , 10. , -2.03666518, 4.14084039,
|
||||
2.21309317, 47.82371559]]),
|
||||
),
|
||||
tuple(
|
||||
47.408654554,
|
||||
47.4038652151,
|
||||
array([[ 0. , 0. , 10. , -2.03666518, 4.14084039,
|
||||
2.21309317, 47.82371559]]),
|
||||
),
|
||||
tuple(
|
||||
47.408654554,
|
||||
47.4038652151,
|
||||
array([[ 0. , 0. , 10. , -2.03666518, 4.14084039,
|
||||
2.21309317, 47.82371559]]),
|
||||
),
|
||||
tuple(
|
||||
47.408654554,
|
||||
47.4038652151,
|
||||
array([[ 0. , 0. , 10. , -2.03666518, 4.14084039,
|
||||
2.21309317, 47.82371559]]),
|
||||
),
|
||||
tuple(
|
||||
22.9327902847,
|
||||
22.9304785991,
|
||||
array([[ 0. , 0. , 10. , -1.63019717, 3.97041764,
|
||||
2.53115835, 38.2051999 ]]),
|
||||
),
|
||||
tuple(
|
||||
28.8119773322,
|
||||
28.8090658953,
|
||||
array([[ 0. , 0. , 10. , -1.14570315, 4.07709911,
|
||||
2.48697441, 49.58615195]]),
|
||||
),
|
||||
tuple(
|
||||
28.8119773322,
|
||||
28.8090658953,
|
||||
array([[ 0. , 0. , 10. , -1.14570315, 4.07709911,
|
||||
2.48697441, 49.58615195]]),
|
||||
),
|
||||
tuple(
|
||||
40.9740600543,
|
||||
40.9699102596,
|
||||
array([[ 0. , 0. , 10. , -0.50178755, 3.83878089,
|
||||
2.93560796, 82.07827571]]),
|
||||
),
|
||||
|
@ -43,7 +43,7 @@ def get_cost_function():
|
||||
dot_input_array = pdme.measurement.input_types.dot_inputs_to_array(dot_inputs)
|
||||
|
||||
actual_dipoles = model.get_dipoles(0, numpy.random.default_rng(SEED_TO_USE))
|
||||
actual_measurements = actual_dipoles.get_dot_measurements(dot_inputs)
|
||||
actual_measurements = actual_dipoles.get_potential_dot_measurements(dot_inputs)
|
||||
actual_measurements_array = numpy.array([m.v for m in actual_measurements])
|
||||
|
||||
def cost_to_use(sample_dipoleses: numpy.ndarray) -> numpy.ndarray:
|
||||
|
@ -43,7 +43,7 @@ def get_cost_function():
|
||||
dot_input_array = pdme.measurement.input_types.dot_inputs_to_array(dot_inputs)
|
||||
|
||||
actual_dipoles = model.get_dipoles(0, numpy.random.default_rng(SEED_TO_USE))
|
||||
actual_measurements = actual_dipoles.get_dot_measurements(dot_inputs)
|
||||
actual_measurements = actual_dipoles.get_potential_dot_measurements(dot_inputs)
|
||||
actual_measurements_array = numpy.array([m.v for m in actual_measurements])
|
||||
|
||||
def cost_to_use(sample_dipoleses: numpy.ndarray) -> numpy.ndarray:
|
||||
|
@ -47,7 +47,7 @@ def get_cost_function():
|
||||
dot_input_array = pdme.measurement.input_types.dot_inputs_to_array(dot_inputs)
|
||||
|
||||
actual_dipoles = model.get_dipoles(0, numpy.random.default_rng(SEED_TO_USE))
|
||||
actual_measurements = actual_dipoles.get_dot_measurements(dot_inputs)
|
||||
actual_measurements = actual_dipoles.get_potential_dot_measurements(dot_inputs)
|
||||
actual_measurements_array = numpy.array([m.v for m in actual_measurements])
|
||||
|
||||
def cost_to_use(sample_dipoleses: numpy.ndarray) -> numpy.ndarray:
|
||||
|
@ -30,3 +30,65 @@
|
||||
]),
|
||||
])
|
||||
# ---
|
||||
# name: test_sort_dipoleses_by_freq
|
||||
list([
|
||||
list([
|
||||
list([
|
||||
100.0,
|
||||
200.0,
|
||||
300.0,
|
||||
400.0,
|
||||
500.0,
|
||||
600.0,
|
||||
0.07,
|
||||
]),
|
||||
list([
|
||||
1.0,
|
||||
2.0,
|
||||
3.0,
|
||||
4.0,
|
||||
5.0,
|
||||
6.0,
|
||||
7.0,
|
||||
]),
|
||||
list([
|
||||
10.0,
|
||||
200.0,
|
||||
30.0,
|
||||
41.0,
|
||||
315.0,
|
||||
0.31,
|
||||
100.0,
|
||||
]),
|
||||
]),
|
||||
list([
|
||||
list([
|
||||
1.0,
|
||||
1.0,
|
||||
1.0,
|
||||
1.0,
|
||||
1.0,
|
||||
1.0,
|
||||
100.0,
|
||||
]),
|
||||
list([
|
||||
22.0,
|
||||
22.2,
|
||||
2.2,
|
||||
222.0,
|
||||
22.0,
|
||||
2.0,
|
||||
200.0,
|
||||
]),
|
||||
list([
|
||||
33.0,
|
||||
33.3,
|
||||
33.0,
|
||||
3.3,
|
||||
0.33,
|
||||
0.3,
|
||||
300.0,
|
||||
]),
|
||||
]),
|
||||
])
|
||||
# ---
|
||||
|
@ -1,4 +1,5 @@
|
||||
import pdme.subspace_simulation
|
||||
import pdme.subspace_simulation.mcmc_costs
|
||||
import numpy
|
||||
|
||||
|
||||
@ -8,3 +9,23 @@ def test_proportional_costs(snapshot):
|
||||
|
||||
actual_result = pdme.subspace_simulation.proportional_cost(a, b).tolist()
|
||||
assert actual_result == snapshot
|
||||
|
||||
|
||||
def test_squared_costs_manual():
|
||||
target = numpy.array([100, 400, 900])
|
||||
approx1 = numpy.array([0, 400, 800])
|
||||
approx2 = numpy.array([200, 400, 600])
|
||||
|
||||
expected1 = 1.0123456790123457
|
||||
expected2 = 1.1111111111111111
|
||||
|
||||
actual1 = pdme.subspace_simulation.mcmc_costs.relative_square_diffs(approx1, target)
|
||||
assert actual1 == expected1
|
||||
|
||||
actual2 = pdme.subspace_simulation.mcmc_costs.relative_square_diffs(approx2, target)
|
||||
assert actual2 == expected2
|
||||
|
||||
combined_actual = pdme.subspace_simulation.mcmc_costs.relative_square_diffs(
|
||||
numpy.array([approx1, approx2]), target
|
||||
)
|
||||
numpy.testing.assert_allclose(combined_actual, [expected1, expected2], rtol=1e-14)
|
||||
|
@ -13,3 +13,28 @@ def test_sort_dipoles_by_freq(snapshot):
|
||||
|
||||
actual_sorted = pdme.subspace_simulation.sort_array_of_dipoles_by_frequency(orig)
|
||||
assert actual_sorted.tolist() == snapshot
|
||||
|
||||
|
||||
def test_sort_dipoleses_by_freq(snapshot):
|
||||
sample_1 = numpy.array(
|
||||
[
|
||||
[1, 2, 3, 4, 5, 6, 7],
|
||||
[100, 200, 300, 400, 500, 600, 0.07],
|
||||
[10, 200, 30, 41, 315, 0.31, 100],
|
||||
]
|
||||
)
|
||||
|
||||
sample_2 = numpy.array(
|
||||
[
|
||||
[1, 1, 1, 1, 1, 1, 100],
|
||||
[33, 33.3, 33, 3.3, 0.33, 0.3, 300],
|
||||
[22, 22.2, 2.2, 222, 22, 2, 200],
|
||||
]
|
||||
)
|
||||
|
||||
original_samples = numpy.array([sample_1, sample_2])
|
||||
|
||||
actual_sorted = pdme.subspace_simulation.sort_array_of_dipoleses_by_frequency(
|
||||
original_samples
|
||||
)
|
||||
assert actual_sorted.tolist() == snapshot
|
||||
|
@ -14,12 +14,12 @@
|
||||
# name: test_fast_spin_qubit_frequency_tarucha_calc
|
||||
list([
|
||||
list([
|
||||
0.24485375860291592,
|
||||
0.009068657726033923,
|
||||
1.1471308805198364,
|
||||
0.042486328908142086,
|
||||
]),
|
||||
list([
|
||||
0.0001499877991005419,
|
||||
0.00011239982757520363,
|
||||
0.0008292459978410822,
|
||||
0.0006214312613006971,
|
||||
]),
|
||||
])
|
||||
# ---
|
||||
|
@ -1,9 +1,24 @@
|
||||
import numpy
|
||||
import pdme.util.fast_nonlocal_spectrum
|
||||
import pdme.measurement
|
||||
import logging
|
||||
import pytest
|
||||
|
||||
|
||||
def dipole_from_array(arr: numpy.ndarray) -> pdme.measurement.OscillatingDipole:
|
||||
return pdme.measurement.OscillatingDipole(arr[0:3], arr[3:6], arr[6])
|
||||
|
||||
|
||||
def s_potential_from_arrays(
|
||||
dipole_array: numpy.ndarray, dotf_pair_array: numpy.ndarray
|
||||
) -> float:
|
||||
dipole = dipole_from_array(dipole_array)
|
||||
r1 = dotf_pair_array[0][0:3]
|
||||
f1 = dotf_pair_array[0][3]
|
||||
r2 = dotf_pair_array[1][0:3]
|
||||
return dipole.s_electric_potential_for_dot_pair(r1, r2, f1)
|
||||
|
||||
|
||||
def test_fast_nonlocal_calc():
|
||||
d1 = [1, 2, 3, 4, 5, 6, 7]
|
||||
d2 = [1, 2, 3, 4, 5, 6, 8]
|
||||
@ -16,21 +31,11 @@ def test_fast_nonlocal_calc():
|
||||
[[[-1, -2, -3, 11], [-1, 2, 5, 11]], [[-1, -2, -3, 6], [2, 4, 6, 6]]]
|
||||
)
|
||||
# expected_ij is for pair i, dipole j
|
||||
expected_11 = 0.000021124454334947546213
|
||||
expected_12 = 0.000022184755131682365135
|
||||
expected_13 = 0.0000053860643617855849275
|
||||
expected_14 = -0.0000023069501696755220764
|
||||
expected_21 = 0.00022356021100884617796
|
||||
expected_22 = 0.00021717277640859343002
|
||||
expected_23 = 0.000017558321044891869169
|
||||
expected_24 = -0.000034714318479634499683
|
||||
|
||||
expected = numpy.array(
|
||||
[
|
||||
[expected_11, expected_21],
|
||||
[expected_12, expected_22],
|
||||
[expected_13, expected_23],
|
||||
[expected_14, expected_24],
|
||||
[s_potential_from_arrays(dipole_array, dot_pair) for dot_pair in dot_pairs]
|
||||
for dipole_array in dipoles
|
||||
]
|
||||
)
|
||||
|
||||
|
@ -1,11 +1,26 @@
|
||||
import numpy
|
||||
import pdme.util.fast_nonlocal_spectrum
|
||||
import pdme.measurement
|
||||
import logging
|
||||
import pytest
|
||||
|
||||
_logger = logging.getLogger(__name__)
|
||||
|
||||
|
||||
def dipole_from_array(arr: numpy.ndarray) -> pdme.measurement.OscillatingDipole:
|
||||
return pdme.measurement.OscillatingDipole(arr[0:3], arr[3:6], arr[6])
|
||||
|
||||
|
||||
def s_potential_from_arrays(
|
||||
dipole_array: numpy.ndarray, dotf_pair_array: numpy.ndarray
|
||||
) -> float:
|
||||
dipole = dipole_from_array(dipole_array)
|
||||
r1 = dotf_pair_array[0][0:3]
|
||||
f1 = dotf_pair_array[0][3]
|
||||
r2 = dotf_pair_array[1][0:3]
|
||||
return dipole.s_electric_potential_for_dot_pair(r1, r2, f1)
|
||||
|
||||
|
||||
def test_fast_nonlocal_calc_multidipole():
|
||||
d1 = [1, 2, 3, 4, 5, 6, 7]
|
||||
d2 = [1, 2, 3, 4, 5, 6, 8]
|
||||
@ -18,19 +33,19 @@ def test_fast_nonlocal_calc_multidipole():
|
||||
[[[-1, -2, -3, 11], [-1, 2, 5, 11]], [[-1, -2, -3, 6], [2, 4, 6, 6]]]
|
||||
)
|
||||
# expected_ij is for pair i, dipole j
|
||||
expected_11 = 0.000021124454334947546213
|
||||
expected_12 = 0.000022184755131682365135
|
||||
expected_13 = 0.0000053860643617855849275
|
||||
expected_14 = -0.0000023069501696755220764
|
||||
expected_21 = 0.00022356021100884617796
|
||||
expected_22 = 0.00021717277640859343002
|
||||
expected_23 = 0.000017558321044891869169
|
||||
expected_24 = -0.000034714318479634499683
|
||||
|
||||
expected = numpy.array(
|
||||
[
|
||||
[expected_11 + expected_12, expected_21 + expected_22],
|
||||
[expected_13 + expected_14, expected_23 + expected_24],
|
||||
[
|
||||
sum(
|
||||
[
|
||||
s_potential_from_arrays(dipole_array, dot_pair)
|
||||
for dipole_array in dipoles
|
||||
]
|
||||
)
|
||||
for dot_pair in dot_pairs
|
||||
]
|
||||
for dipoles in dipoleses
|
||||
]
|
||||
)
|
||||
|
||||
|
@ -1,6 +1,30 @@
|
||||
import numpy
|
||||
import pdme.util.fast_v_calc
|
||||
|
||||
import pdme.measurement
|
||||
|
||||
|
||||
def dipole_from_array(arr: numpy.ndarray) -> pdme.measurement.OscillatingDipole:
|
||||
return pdme.measurement.OscillatingDipole(arr[0:3], arr[3:6], arr[6])
|
||||
|
||||
|
||||
def s_potential_from_arrays(
|
||||
dipole_array: numpy.ndarray, dotf_array: numpy.ndarray
|
||||
) -> float:
|
||||
dipole = dipole_from_array(dipole_array)
|
||||
r = dotf_array[0:3]
|
||||
f = dotf_array[3]
|
||||
return dipole.s_electric_potential_at_position(r, f)
|
||||
|
||||
|
||||
def s_electric_field_x_from_arrays(
|
||||
dipole_array: numpy.ndarray, dotf_array: numpy.ndarray
|
||||
) -> float:
|
||||
dipole = dipole_from_array(dipole_array)
|
||||
r = dotf_array[0:3]
|
||||
f = dotf_array[3]
|
||||
return dipole.s_electric_fieldx_at_position(r, f)
|
||||
|
||||
|
||||
def test_fast_v_calc():
|
||||
d1 = [1, 2, 3, 4, 5, 6, 7]
|
||||
@ -9,11 +33,11 @@ def test_fast_v_calc():
|
||||
dipoles = numpy.array([d1, d2])
|
||||
|
||||
dot_inputs = numpy.array([[-1, -1, -1, 11], [2, 3, 1, 5.5]])
|
||||
# expected_ij is for dot i, dipole j
|
||||
expected_11 = 0.00001421963287022476
|
||||
expected_12 = 0.00001107180225755457
|
||||
expected_21 = 0.000345021108583681380388722
|
||||
expected_22 = 0.0000377061050587914705139781
|
||||
|
||||
expected_11 = s_potential_from_arrays(dipoles[0], dot_inputs[0])
|
||||
expected_12 = s_potential_from_arrays(dipoles[1], dot_inputs[0])
|
||||
expected_21 = s_potential_from_arrays(dipoles[0], dot_inputs[1])
|
||||
expected_22 = s_potential_from_arrays(dipoles[1], dot_inputs[1])
|
||||
|
||||
expected = numpy.array([[expected_11, expected_21], [expected_12, expected_22]])
|
||||
|
||||
@ -31,11 +55,11 @@ def test_fast_v_calc_multidipoles():
|
||||
dipoles = numpy.array([[d1, d2]])
|
||||
|
||||
dot_inputs = numpy.array([[-1, -1, -1, 11], [2, 3, 1, 5.5]])
|
||||
# expected_ij is for dot i, dipole j
|
||||
expected_11 = 0.00001421963287022476
|
||||
expected_12 = 0.00001107180225755457
|
||||
expected_21 = 0.000345021108583681380388722
|
||||
expected_22 = 0.0000377061050587914705139781
|
||||
|
||||
expected_11 = s_potential_from_arrays(dipoles[0][0], dot_inputs[0])
|
||||
expected_12 = s_potential_from_arrays(dipoles[0][1], dot_inputs[0])
|
||||
expected_21 = s_potential_from_arrays(dipoles[0][0], dot_inputs[1])
|
||||
expected_22 = s_potential_from_arrays(dipoles[0][1], dot_inputs[1])
|
||||
|
||||
expected = numpy.array([[expected_11 + expected_12, expected_21 + expected_22]])
|
||||
|
||||
@ -48,7 +72,7 @@ def test_fast_v_calc_multidipoles():
|
||||
|
||||
def test_fast_v_calc_big_multidipole():
|
||||
|
||||
dipoles = numpy.array(
|
||||
dipoleses = numpy.array(
|
||||
[
|
||||
[
|
||||
[1, 1, 5, 6, 3, 1, 1],
|
||||
@ -72,27 +96,113 @@ def test_fast_v_calc_big_multidipole():
|
||||
]
|
||||
)
|
||||
|
||||
expected = numpy.array(
|
||||
expected = [
|
||||
[
|
||||
sum(
|
||||
[
|
||||
s_potential_from_arrays(dipole_array, dot_input)
|
||||
for dipole_array in dipole_config
|
||||
]
|
||||
)
|
||||
for dot_input in dot_inputs
|
||||
]
|
||||
for dipole_config in dipoleses
|
||||
]
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
pdme.util.fast_v_calc.fast_vs_for_dipoleses(dot_inputs, dipoleses),
|
||||
expected,
|
||||
err_msg="Voltages at dot aren't as expected for multidipole calc.",
|
||||
)
|
||||
|
||||
|
||||
def test_fast_electric_field_x_calc():
|
||||
d1 = [1, 2, 3, 4, 5, 6, 7]
|
||||
d2 = [2, 5, 3, 4, -5, -6, 2]
|
||||
|
||||
dipoles = numpy.array([d1, d2])
|
||||
|
||||
dot_inputs = numpy.array([[-1, -1, -1, 11], [2, 3, 1, 5.5]])
|
||||
|
||||
expected_11 = s_electric_field_x_from_arrays(dipoles[0], dot_inputs[0])
|
||||
expected_12 = s_electric_field_x_from_arrays(dipoles[1], dot_inputs[0])
|
||||
expected_21 = s_electric_field_x_from_arrays(dipoles[0], dot_inputs[1])
|
||||
expected_22 = s_electric_field_x_from_arrays(dipoles[1], dot_inputs[1])
|
||||
|
||||
expected = numpy.array([[expected_11, expected_21], [expected_12, expected_22]])
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
pdme.util.fast_v_calc.fast_efieldxs_for_dipoles(dot_inputs, dipoles),
|
||||
expected,
|
||||
err_msg="E x fast calc at dot aren't as expected.",
|
||||
)
|
||||
|
||||
|
||||
def test_fast_electric_field_x_calc_multidipoles():
|
||||
d1 = [1, 2, 3, 4, 5, 6, 7]
|
||||
d2 = [2, 5, 3, 4, -5, -6, 2]
|
||||
|
||||
dipoles = numpy.array([[d1, d2]])
|
||||
|
||||
dot_inputs = numpy.array([[-1, -1, -1, 11], [2, 3, 1, 5.5]])
|
||||
|
||||
expected_11 = s_electric_field_x_from_arrays(dipoles[0][0], dot_inputs[0])
|
||||
expected_12 = s_electric_field_x_from_arrays(dipoles[0][1], dot_inputs[0])
|
||||
expected_21 = s_electric_field_x_from_arrays(dipoles[0][0], dot_inputs[1])
|
||||
expected_22 = s_electric_field_x_from_arrays(dipoles[0][1], dot_inputs[1])
|
||||
|
||||
expected = numpy.array([[expected_11 + expected_12, expected_21 + expected_22]])
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
pdme.util.fast_v_calc.fast_efieldxs_for_dipoleses(dot_inputs, dipoles),
|
||||
expected,
|
||||
err_msg="E x fast calc at dot aren't as expected for multidipole calc.",
|
||||
)
|
||||
|
||||
|
||||
def test_fast_electric_field_x_calc_big_multidipole():
|
||||
|
||||
dipoleses = numpy.array(
|
||||
[
|
||||
[
|
||||
0.0010151687742365581690202135,
|
||||
0.00077627527320628609782627266,
|
||||
0.00043313471258511003340648713,
|
||||
0.000077184305988088453637005111,
|
||||
[1, 1, 5, 6, 3, 1, 1],
|
||||
[5, 3, 2, 13, 1, 1, 2],
|
||||
[-5, -5, -3, -1, -3, 8, 3],
|
||||
],
|
||||
[
|
||||
0.000041099091967966890657097060,
|
||||
0.0019377687238977568792327845,
|
||||
0.0085903193415282984161225029,
|
||||
0.00014557676715208209310911838,
|
||||
[-3, -1, -2, -2, -6, 3, 4],
|
||||
[8, 0, 2, 0, 1, 5, 5],
|
||||
[1, 4, -4, -1, -3, -5, 6],
|
||||
],
|
||||
]
|
||||
)
|
||||
|
||||
dot_inputs = numpy.array(
|
||||
[
|
||||
[1, 1, 0, 1],
|
||||
[2, 5, 6, 2],
|
||||
[3, 1, 3, 3],
|
||||
[0.5, 0.5, 0.5, 4],
|
||||
]
|
||||
)
|
||||
|
||||
expected = [
|
||||
[
|
||||
sum(
|
||||
[
|
||||
s_electric_field_x_from_arrays(dipole_array, dot_input)
|
||||
for dipole_array in dipole_config
|
||||
]
|
||||
)
|
||||
for dot_input in dot_inputs
|
||||
]
|
||||
for dipole_config in dipoleses
|
||||
]
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
pdme.util.fast_v_calc.fast_vs_for_dipoleses(dot_inputs, dipoles),
|
||||
pdme.util.fast_v_calc.fast_efieldxs_for_dipoleses(dot_inputs, dipoleses),
|
||||
expected,
|
||||
err_msg="Voltages at dot aren't as expected for multidipole calc.",
|
||||
err_msg="E x fast calc at dot aren't as expected for multidipole calc.",
|
||||
)
|
||||
|
||||
|
||||
@ -107,28 +217,3 @@ def test_between():
|
||||
expected = numpy.array([False, False, True])
|
||||
|
||||
numpy.testing.assert_array_equal(actual, expected, err_msg="Between calc wrong")
|
||||
|
||||
|
||||
def test_fast_v_calc_asymmetric_multidipoles_but_symmetric():
|
||||
# expected format is [px, py, pz, sx, sy, sz, e1, e2, w]
|
||||
d1 = [1, 2, 3, 4, 5, 6, 1, 1, 7 / 2]
|
||||
d2 = [2, 5, 3, 4, -5, -6, 2, 2, 2 / 2]
|
||||
|
||||
dipoles = numpy.array([[d1, d2]])
|
||||
|
||||
dot_inputs = numpy.array([[-1, -1, -1, 11], [2, 3, 1, 5.5]])
|
||||
# expected_ij is for dot i, dipole j
|
||||
expected_11 = 0.00001421963287022476
|
||||
expected_12 = 0.00001107180225755457
|
||||
expected_21 = 0.000345021108583681380388722
|
||||
expected_22 = 0.0000377061050587914705139781
|
||||
|
||||
expected = numpy.array([[expected_11 + expected_12, expected_21 + expected_22]])
|
||||
|
||||
numpy.testing.assert_allclose(
|
||||
pdme.util.fast_v_calc.fast_vs_for_asymmetric_dipoleses(
|
||||
dot_inputs, dipoles, 1e10
|
||||
),
|
||||
expected,
|
||||
err_msg="Voltages at dot aren't as expected for multidipole calc.",
|
||||
)
|
||||
|
Loading…
x
Reference in New Issue
Block a user