dominik/nmreval
1
0
Fork 0
forked from IPKM/nmreval
Code Pull Requests Activity

added f-omega option to bds fit functions; closes #133

Browse Source
This commit is contained in:
Dominik Demuth 2023-11-09 17:26:57 +01:00
parent 070509c691
commit dbb35cdba4
1 changed files with 91 additions and 31 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

122
src/nmreval/models/bds.py
View File

@ -10,14 +10,20 @@ class _AbstractBDS:
equation = ''
params = [r'\Delta\epsilon', r'\tau_{0}']
bounds = [(0, None), (0, None)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
susceptibility = None
iscomplex = True
@classmethod
def func(cls, x, *args, complex_mode: int = 0, **kwargs):
def func(cls, x, *args, complex_mode: int = 0, xaxis: str = 'freq', **kwargs):
# args[0] : Delta epsilon
# args[1:] : every other parameter
chi = args[0] * cls.susceptibility(2*np.pi*x, *args[1:], **kwargs)
_w = _convert_x_to_omega(x, xaxis=xaxis)
chi = args[0] * cls.susceptibility(_w, *args[1:], **kwargs)
if complex_mode == 0:
return chi
elif complex_mode == 1:
@ -64,11 +70,16 @@ class HavriliakNegamiAlphaGammaBDS:
equation = r'\Delta\epsilon / [1-(i\omega\tau)^{\gamma}]^{\alpha}'
params = [r'\Delta\epsilon', r'\tau_{0}', r'\alpha', r'\alpha\gamma']
bounds = [(0, None), (0, None), (0, 1), (0, 1)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
iscomplex = True
@staticmethod
def func(x, deps, tau, alpha, alphagamma, complex_mode: int = 0, **kwargs):
chi = deps * HavriliakNegami.susceptibility(2*np.pi*x, tau, alpha, alphagamma/alpha, **kwargs)
def func(x, deps, tau, alpha, alphagamma, complex_mode: int = 0, xaxis: str = 'freq', **kwargs):
_w = _convert_x_to_omega(x, xaxis=xaxis)
chi = deps * HavriliakNegami.susceptibility(_w, tau, alpha, alphagamma/alpha, **kwargs)
if complex_mode == 0:
return chi
elif complex_mode == 1:
@ -115,13 +126,17 @@ class HNWithHF:
equation = r'\Delta\epsilon HN(\omega, \tau, \alpha, \gamma) / CD(\omega, \tau_{c}, \alpha\gamma-\delta)'
params = [r'\Delta\epsilon', r'\tau', r'\alpha', r'\gamma', r'\tau_{c}', '\delta']
bounds = [(0, None), (0, None), (0, 1), (0, 1), (0, None), (0, 1)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
iscomplex = True
@staticmethod
def func(x, deps, tau, alpha, gamma, tauc, delta, complex_mode: int = 0):
w = 2 * np.pi * x
numerator = (1 + 1j*w*tauc) ** (gamma-delta)
denominator = (1 + (1j*w*tau)**alpha) ** gamma
def func(x, deps, tau, alpha, gamma, tauc, delta, complex_mode: int = 0, xaxis: str = 'freq'):
_w = _convert_x_to_omega(x, xaxis=xaxis)
numerator = (1 + 1j*_w*tauc) ** (gamma-delta)
denominator = (1 + (1j*_w*tau)**alpha) ** gamma
epsilon = deps * np.conjugate(numerator / denominator)
@ -144,8 +159,8 @@ class _CCWithHF:
iscomplex = True
@staticmethod
def func(x, deps, tau, alpha, tauc, delta, complex_mode: int = 0):
return HNWithHF.func(x, deps, tau, alpha, 1, tauc, delta, complex_mode=complex_mode)
def func(x, deps, tau, alpha, tauc, delta, complex_mode: int = 0, xaxis: str = 'freq'):
return HNWithHF.func(x, deps, tau, alpha, 1, tauc, delta, complex_mode=complex_mode, xaxis=xaxis)
class _CDWithHF:
@ -157,8 +172,8 @@ class _CDWithHF:
iscomplex = True
@staticmethod
def func(x, deps, tau, gamma, tauc, delta, complex_mode: int = 0):
return HNWithHF.func(x, deps, tau, 1, gamma, tauc, delta, complex_mode=complex_mode)
def func(x, deps, tau, gamma, tauc, delta, complex_mode: int = 0, xaxis: str = 'freq'):
return HNWithHF.func(x, deps, tau, 1, gamma, tauc, delta, complex_mode=complex_mode, xaxis=xaxis)
class PowerLawBDS:
@ -167,16 +182,21 @@ class PowerLawBDS:
equation = r'A / (i\omega)^{n}'
params = ['A', 'n']
bounds = [(None, None), (None, None)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
iscomplex = True
@staticmethod
def func(x, a, n, complex_mode: int = 0):
def func(x, a, n, complex_mode: int = 0, xaxis: str = 'freq'):
_w = _convert_x_to_omega(x, xaxis=xaxis)
if complex_mode == 0:
ret_val = np.exp(1j*n*np.pi/2) * a / x**n
ret_val = np.exp(1j*n*np.pi/2) * a / _w**n
elif complex_mode == 1:
ret_val = np.cos(n*np.pi/2) * a / x**n
ret_val = np.cos(n*np.pi/2) * a / _w**n
elif complex_mode == 2:
ret_val = np.sin(n*np.pi/2) * a / x**n
ret_val = np.sin(n*np.pi/2) * a / _w**n
else:
raise ValueError(f'{complex_mode!r} is not 0, 1, 2')
@ -189,17 +209,22 @@ class DCCondBDS:
equation = r'i\sigma_{dc}/\epsilon_{0}\omega'
params = [r'\sigma_{dc}']
bounds = [(0, None)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
iscomplex = True
@staticmethod
def func(x, sigma, complex_mode: int = 0):
def func(x, sigma, complex_mode: int = 0, xaxis: str = 'freq'):
_w = _convert_x_to_omega(x, xaxis=xaxis)
if complex_mode == 0:
ret_val = np.zeros(x.shape, dtype=complex)
ret_val += 1j * sigma / x / epsilon0
ret_val += 1j * sigma / _w / epsilon0
elif complex_mode == 1:
ret_val = np.zeros(x.shape)
elif complex_mode == 2:
ret_val = sigma / x / epsilon0
ret_val = sigma / _w / epsilon0
else:
raise ValueError(f'{complex_mode!r} is not 0, 1, 2')
@ -210,11 +235,16 @@ class DerivativeHavriliakNegami:
name = 'Derivative HN'
type = 'Dielectric Spectroscopy'
params = [r'\Delta\epsilon', r'\tau', r'\alpha', r'\gamma']
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
bounds = [(0, None), (0, None), (0, 1), (0, 1)]
@staticmethod
def func(x, eps, tau, a, g):
omtau = 2*np.pi*x * tau
def func(x, eps, tau, a, g, xaxis: str = 'freq'):
w = _convert_x_to_omega(x, xaxis=xaxis)
omtau = w * tau
theta = np.arctan2(np.sin(np.pi*a/2.), omtau**(-a) + np.cos(np.pi*a/2.))
numer = a*g * omtau**a * np.cos(np.pi*a/2. - (1+g)*theta)
@ -228,10 +258,15 @@ class DerivativeColeCole:
type = 'Dielectric Spectroscopy'
params = [r'\Delta\epsilon', r'\tau', r'\alpha']
bounds = [(0, None), (0, None), (0, 1)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
@staticmethod
def func(x, eps, tau, alpha):
omtau = 2*np.pi*x * tau
def func(x, eps, tau, alpha, xaxis: str = 'freq'):
w = _convert_x_to_omega(x, xaxis=xaxis)
omtau = w * tau
theta = np.arctan2(np.sin(np.pi*alpha/2.), omtau**(-alpha) + np.cos(np.pi*alpha/2.))
numer = alpha * omtau**alpha * np.cos(np.pi*alpha/2. - 2*theta)
@ -245,10 +280,14 @@ class DerivativeColeDavidson:
type = 'Dielectric Spectroscopy'
params = [r'\Delta\epsilon', r'\tau', r'\gamma']
bounds = [(0, None), (0, None), (0, 1)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
@staticmethod
def func(x, eps, tau, g):
omtau = 2*np.pi*x * tau
def func(x, eps, tau, g, xaxis: str = 'freq'):
w = _convert_x_to_omega(x, xaxis=xaxis)
omtau = w * tau
numer = g * omtau * np.sin((1+g)*np.arctan(omtau))
denom = (1 + omtau**2)**((1+g)/2.)
@ -260,10 +299,13 @@ class _DerivativeHNWithHF:
type = 'Dielectric Spectroscopy'
params = [r'\Delta\epsilon', r'\tau', r'\alpha', r'\gamma', r'\tau_{c}', r'\delta']
bounds = [(0, None), (0, None), (0, 1), (0, 1), (0, None), (0, 1)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
@staticmethod
def func(x, deps, tau, alpha, gamma, tauc, delta):
w = 2*np.pi*x
def func(x, deps, tau, alpha, gamma, tauc, delta, xaxis: str = 'freq'):
w = _convert_x_to_omega(x, xaxis=xaxis)
a_pi2 = alpha*np.pi/2
w_lamb = (w*tau)**alpha
@ -287,10 +329,13 @@ class DerivativeCCWithHF:
type = 'Dielectric Spectroscopy'
params = [r'\Delta\epsilon', r'\tau', r'\alpha', r'\tau_{c}', r'\delta']
bounds = [(0, None), (0, None), (0, 1), (0, None), (0, 1)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
@staticmethod
def func(x, deps, tau, alpha, tauc, delta):
return _DerivativeHNWithHF.func(x, deps, tau, alpha, 1, tauc, delta)
def func(x, deps, tau, alpha, tauc, delta, xaxis: str = 'freq'):
return _DerivativeHNWithHF.func(x, deps, tau, alpha, 1, tauc, delta, xaxis=xaxis)
class DerivativeCDWithHF:
@ -298,7 +343,22 @@ class DerivativeCDWithHF:
type = 'Dielectric Spectroscopy'
params = [r'\Delta\epsilon', r'\tau', r'\gamma', r'\tau_{c}', r'\delta']
bounds = [(0, None), (0, None), (0, 1), (0, None), (0, 1)]
choices = [
('x axis', 'xaxis', {'Frequency': 'freq', 'Omega': 'omega'})
]
@staticmethod
def func(x, deps, tau, gamma, tauc, delta):
return _DerivativeHNWithHF.func(x, deps, tau, 1, gamma, tauc, delta)
def func(x, deps, tau, gamma, tauc, delta, xaxis: str = 'freq'):
return _DerivativeHNWithHF.func(x, deps, tau, 1, gamma, tauc, delta, xaxis=xaxis)
def _convert_x_to_omega(x, xaxis: str = 'freq'):
if xaxis not in ['freq', 'omega']:
raise ValueError(f'Argument `xaxis` is `freq` or `omega`, given is {xaxis!r}')
if xaxis == 'freq':
_w = 2 * np.pi * x
else:
_w = x
return _w