mirror of
https://github.com/nqrduck/nqr-blochsimulator.git
synced 2024-10-02 01:50:39 +00:00
Compare commits
No commits in common. "83a488974834661a73d53563c080e45fe9a77628" and "f845135dc2eaa9b3a61d12730816a71b8c20035c" have entirely different histories.
83a4889748
...
f845135dc2
3 changed files with 23 additions and 113 deletions
|
@ -1,7 +0,0 @@
|
||||||
"""The nqr_blochsimulator package contains the classes necessary to simulate NQR experiments using the Bloch equations."""
|
|
||||||
|
|
||||||
from .classes.sample import Sample
|
|
||||||
from .classes.simulation import Simulation
|
|
||||||
from .classes.pulse import PulseArray
|
|
||||||
|
|
||||||
__all__ = ["Sample", "Simulation", "PulseArray"]
|
|
|
@ -1,7 +1,5 @@
|
||||||
import logging
|
from math import pi, sqrt
|
||||||
import numpy as np
|
|
||||||
|
|
||||||
logger = logging.getLogger(__name__)
|
|
||||||
|
|
||||||
class Sample:
|
class Sample:
|
||||||
"""
|
"""
|
||||||
|
@ -12,18 +10,18 @@ class Sample:
|
||||||
|
|
||||||
def __init__(
|
def __init__(
|
||||||
self,
|
self,
|
||||||
name : str,
|
name,
|
||||||
density : float,
|
density,
|
||||||
molar_mass : float,
|
molar_mass,
|
||||||
resonant_frequency : float,
|
resonant_frequency,
|
||||||
gamma : float,
|
gamma,
|
||||||
nuclear_spin : float,
|
nuclear_spin,
|
||||||
spin_factor : float,
|
spin_factor,
|
||||||
powder_factor : float,
|
powder_factor,
|
||||||
filling_factor : float,
|
filling_factor,
|
||||||
T1 : float,
|
T1,
|
||||||
T2 : float,
|
T2,
|
||||||
T2_star : float,
|
T2_star,
|
||||||
atom_density=None,
|
atom_density=None,
|
||||||
sample_volume=None,
|
sample_volume=None,
|
||||||
sample_length=None,
|
sample_length=None,
|
||||||
|
@ -47,7 +45,7 @@ class Sample:
|
||||||
nuclear_spin : float
|
nuclear_spin : float
|
||||||
The nuclear spin quantum number of the sample.
|
The nuclear spin quantum number of the sample.
|
||||||
spin_factor : float
|
spin_factor : float
|
||||||
The spin transition factor of the sample.
|
The spin factor of the sample.
|
||||||
powder_factor : float
|
powder_factor : float
|
||||||
The powder factor of the sample.
|
The powder factor of the sample.
|
||||||
filling_factor : float
|
filling_factor : float
|
||||||
|
@ -102,69 +100,3 @@ class Sample:
|
||||||
)
|
)
|
||||||
else:
|
else:
|
||||||
self.atoms = self.avogadro * self.density / self.molar_mass
|
self.atoms = self.avogadro * self.density / self.molar_mass
|
||||||
|
|
||||||
def pauli_spin_matrices(self, spin):
|
|
||||||
"""
|
|
||||||
Generate the spin matrices for a given spin value.
|
|
||||||
|
|
||||||
Parameters:
|
|
||||||
spin (float): The spin value, which can be a half-integer or integer.
|
|
||||||
|
|
||||||
Returns:
|
|
||||||
tuple: A tuple containing the following elements:
|
|
||||||
Jx (np.ndarray): The x-component of the spin matrix.
|
|
||||||
Jy (np.ndarray): The y-component of the spin matrix.
|
|
||||||
Jz (np.ndarray): The z-component of the spin matrix.
|
|
||||||
J_minus (np.ndarray): The lowering operator matrix.
|
|
||||||
J_plus (np.ndarray): The raising operator matrix.
|
|
||||||
m (np.ndarray): The array of magnetic quantum numbers.
|
|
||||||
"""
|
|
||||||
|
|
||||||
m = np.arange(spin, -spin-1, -1)
|
|
||||||
paulirowlength = int(spin * 2 + 1)
|
|
||||||
|
|
||||||
pauli_z = np.diag(m)
|
|
||||||
pauli_plus = np.zeros((paulirowlength, paulirowlength))
|
|
||||||
pauli_minus = np.zeros((paulirowlength, paulirowlength))
|
|
||||||
|
|
||||||
for row_index in range(paulirowlength - 1):
|
|
||||||
col_index = row_index + 1
|
|
||||||
pauli_plus[row_index, col_index] = np.sqrt(spin * (spin + 1) - m[col_index] * (m[col_index] + 1))
|
|
||||||
|
|
||||||
for row_index in range(1, paulirowlength):
|
|
||||||
col_index = row_index - 1
|
|
||||||
pauli_minus[row_index, col_index] = np.sqrt(spin * (spin + 1) - m[col_index] * (m[col_index] - 1))
|
|
||||||
|
|
||||||
Jx = 0.5 * (pauli_plus + pauli_minus)
|
|
||||||
Jy = -0.5j * (pauli_plus - pauli_minus)
|
|
||||||
Jz = pauli_z
|
|
||||||
|
|
||||||
return Jx, Jy, Jz, pauli_minus, pauli_plus, m
|
|
||||||
|
|
||||||
def calculate_spin_transition_factor(self, I, transition):
|
|
||||||
"""
|
|
||||||
Calculate the prefactor for the envisaged spin transition for a given nuclear spin.
|
|
||||||
|
|
||||||
Parameters:
|
|
||||||
I (float): The nuclear spin value, which can be a half-integer or integer.
|
|
||||||
transition (int): The index of the transition.
|
|
||||||
The transition indices represent the shifts between magnetic quantum numbers m:
|
|
||||||
- 0 represents -1/2 --> 1/2
|
|
||||||
- 1 represents 1/2 --> 3/2
|
|
||||||
- 2 represents 3/2 --> 5/2
|
|
||||||
(only valid transitions based on spin value I are allowed)
|
|
||||||
|
|
||||||
Returns:
|
|
||||||
float: The prefactor for the envisaged spin transition.
|
|
||||||
"""
|
|
||||||
m_values = np.arange(I, -I-1, -1)
|
|
||||||
if transition < 0 or transition >= len(m_values) - 1:
|
|
||||||
raise ValueError(f"Invalid transition for spin {I}. Valid range is 0 to {len(m_values) - 2}")
|
|
||||||
|
|
||||||
Jx, Jy, Jz, J_minus, J_plus, m = self.pauli_spin_matrices(I)
|
|
||||||
trindex = int(len(Jx) / 2 - transition)
|
|
||||||
spinfactor = Jx[trindex - 1, trindex]
|
|
||||||
|
|
||||||
logger.debug(f"Spin transition factor for I={I}, transition={transition}: {np.real(spinfactor)}")
|
|
||||||
logger.info(f"Jx is {Jx}")
|
|
||||||
return np.real(spinfactor)
|
|
||||||
|
|
|
@ -1,13 +1,10 @@
|
||||||
import unittest
|
import unittest
|
||||||
import numpy as np
|
import numpy as np
|
||||||
import logging
|
|
||||||
import matplotlib.pyplot as plt
|
import matplotlib.pyplot as plt
|
||||||
from nqr_blochsimulator.classes.sample import Sample
|
from nqr_blochsimulator.classes.sample import Sample
|
||||||
from nqr_blochsimulator.classes.simulation import Simulation
|
from nqr_blochsimulator.classes.simulation import Simulation
|
||||||
from nqr_blochsimulator.classes.pulse import PulseArray
|
from nqr_blochsimulator.classes.pulse import PulseArray
|
||||||
|
|
||||||
logger = logging.getLogger(__name__)
|
|
||||||
logging.basicConfig(level=logging.INFO)
|
|
||||||
|
|
||||||
class TestSimulation(unittest.TestCase):
|
class TestSimulation(unittest.TestCase):
|
||||||
def setUp(self):
|
def setUp(self):
|
||||||
|
@ -16,14 +13,14 @@ class TestSimulation(unittest.TestCase):
|
||||||
density=1.585e6, # g/m^3
|
density=1.585e6, # g/m^3
|
||||||
molar_mass=440.3, # g/mol
|
molar_mass=440.3, # g/mol
|
||||||
resonant_frequency=83.56e6, # Hz
|
resonant_frequency=83.56e6, # Hz
|
||||||
gamma=43.42, # MHz/T
|
gamma=4.342e7, # Hz/T
|
||||||
nuclear_spin=9 / 2,
|
nuclear_spin=9 / 2,
|
||||||
spin_factor=1.94
|
spin_factor=2,
|
||||||
powder_factor=0.75,
|
powder_factor=0.75,
|
||||||
filling_factor=0.7,
|
filling_factor=0.7,
|
||||||
T1=83.5, # µs
|
T1=83.5e-5, # s
|
||||||
T2=396, # µs
|
T2=396e-6, # s
|
||||||
T2_star=50, # µs
|
T2_star=50e-6, # s
|
||||||
)
|
)
|
||||||
|
|
||||||
simulation_length = 300e-6
|
simulation_length = 300e-6
|
||||||
|
@ -48,8 +45,8 @@ class TestSimulation(unittest.TestCase):
|
||||||
initial_magnetization=1,
|
initial_magnetization=1,
|
||||||
gradient=1,
|
gradient=1,
|
||||||
noise=0.5,
|
noise=0.5,
|
||||||
length_coil=6, # mm
|
length_coil=6e-3,
|
||||||
diameter_coil=3, # mm
|
diameter_coil=3e-3,
|
||||||
number_turns=9,
|
number_turns=9,
|
||||||
q_factor_transmit=100,
|
q_factor_transmit=100,
|
||||||
q_factor_receive=100,
|
q_factor_receive=100,
|
||||||
|
@ -60,7 +57,7 @@ class TestSimulation(unittest.TestCase):
|
||||||
temperature=300,
|
temperature=300,
|
||||||
loss_TX=12,
|
loss_TX=12,
|
||||||
loss_RX=12,
|
loss_RX=12,
|
||||||
conversion_factor=2884, # This is for the LimeSDR based spectrometer
|
conversion_factor=2884 # This is for the LimeSDR based spectrometer
|
||||||
)
|
)
|
||||||
|
|
||||||
def test_simulation(self):
|
def test_simulation(self):
|
||||||
|
@ -72,15 +69,3 @@ class TestSimulation(unittest.TestCase):
|
||||||
plt.ylabel("Magnetization (a.u.)")
|
plt.ylabel("Magnetization (a.u.)")
|
||||||
plt.title("FID of BiPh3")
|
plt.title("FID of BiPh3")
|
||||||
plt.show()
|
plt.show()
|
||||||
|
|
||||||
def test_spin_factor_calculation(self):
|
|
||||||
|
|
||||||
spin = self.sample.nuclear_spin
|
|
||||||
transition = 1
|
|
||||||
|
|
||||||
spin_factor = self.sample.calculate_spin_transition_factor(spin, transition)
|
|
||||||
logger.info("Spin factor: " + str(spin_factor))
|
|
||||||
|
|
||||||
|
|
||||||
if __name__ == "__main__":
|
|
||||||
unittest.main()
|
|
||||||
|
|
Loading…
Reference in a new issue