Commit 8311fd7b authored by Artem Pulkin's avatar Artem Pulkin

potentials: add epsilon to 3p descriptors

parent a81e35a0
Pipeline #43483 passed with stage
in 2 minutes and 29 seconds
......@@ -802,7 +802,7 @@ def kernel_mlsf_g5(
double[::1] r_data,
double[:, ::1] cartesian_row,
double[:, ::1] cartesian_col,
double a, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles,
double a, double epsilon, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles,
int[::1] species_row,
int[::1] species_mask,
double[::1] out,
......@@ -820,7 +820,7 @@ def kernel_mlsf_g5(
# --- preamble ---
cdef double _fn_cutoff1, _fn_cutoff2, _fn_exponent, _fn_pw
cdef double _prefactor = pow(2, 1 - zeta)
cdef double _prefactor = pow(2, 1 - zeta) * epsilon
cdef int pre_compute_r_cutoff_fn_handle = pre_compute_r_handles[0], pre_compute_r_exp_fn_handle = pre_compute_r_handles[2]
# ----------------
......@@ -870,7 +870,7 @@ def kernel_g_mlsf_g5(
double[::1] r_data,
double[:, ::1] cartesian_row,
double[:, ::1] cartesian_col,
double a, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles,
double a, double epsilon, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles,
int[::1] species_row,
int[::1] species_mask,
double[:, :, ::1] out,
......@@ -889,7 +889,7 @@ def kernel_g_mlsf_g5(
# --- preamble ---
cdef double _fn_cutoff1, _fn_cutoff2, _fn_exponent, _fn_pw
cdef double _prefactor = pow(2, 1 - zeta)
cdef double _prefactor = pow(2, 1 - zeta) * epsilon
cdef int pre_compute_r_cutoff_fn_handle = pre_compute_r_handles[0], pre_compute_r_exp_fn_handle = pre_compute_r_handles[2]
cdef double _fp_cutoff1, _fp_cutoff2
cdef int pre_compute_r_cutoff_fp_handle = pre_compute_r_handles[1]
......@@ -968,7 +968,7 @@ def kernel_mlsf_g4(
double[::1] r_data,
double[:, ::1] cartesian_row,
double[:, ::1] cartesian_col,
double a, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles,
double a, double epsilon, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles,
int[::1] species_row,
int[::1] species_mask,
double[::1] out,
......@@ -986,7 +986,7 @@ def kernel_mlsf_g4(
# --- preamble ---
cdef double _fn_cutoff1, _fn_cutoff2, _fn_cutoff3, _fn_exponent, _fn_pw
cdef double _prefactor = pow(2, 1 - zeta)
cdef double _prefactor = pow(2, 1 - zeta) * epsilon
cdef int pre_compute_r_cutoff_fn_handle = pre_compute_r_handles[0], pre_compute_r_exp_fn_handle = pre_compute_r_handles[2]
cdef double _r3, _r3_factor, g5_fun
# ----------------
......@@ -1039,7 +1039,7 @@ def kernel_g_mlsf_g4(
double[::1] r_data,
double[:, ::1] cartesian_row,
double[:, ::1] cartesian_col,
double a, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles,
double a, double epsilon, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles,
int[::1] species_row,
int[::1] species_mask,
double[:, :, ::1] out,
......@@ -1058,7 +1058,7 @@ def kernel_g_mlsf_g4(
# --- preamble ---
cdef double _fn_cutoff1, _fn_cutoff2, _fn_cutoff3, _fn_exponent, _fn_pw
cdef double _prefactor = pow(2, 1 - zeta)
cdef double _prefactor = pow(2, 1 - zeta) * epsilon
cdef int pre_compute_r_cutoff_fn_handle = pre_compute_r_handles[0], pre_compute_r_exp_fn_handle = pre_compute_r_handles[2]
cdef double _r3, _r3_factor, g5_fun
cdef double _fp_cutoff1, _fp_cutoff2
......
......@@ -117,7 +117,7 @@ known_families = {} # This dict is here for restoring json representations of p
class LocalPotentialFamily:
def __init__(self, coordination_number, parameters, cutoff, kernel, kernel_gradient=None, parameter_units=None,
tag=None, proto=None, pre_compute_r=None):
parameter_defaults=None, tag=None, proto=None, pre_compute_r=None):
"""
Represents a local potential family with the same potential shape and
arbitrary parameter values.
......@@ -136,6 +136,8 @@ class LocalPotentialFamily:
A callable for evaluating potential gradients.
parameter_units : ParameterUnits
A dictionary with units per parameter.
parameter_defaults : dict
A dictionary with parameter default values.
tag : str
Optional tag
proto : class
......@@ -153,6 +155,7 @@ class LocalPotentialFamily:
self.coordination_number = coordination_number
self.parameters = parameters
self.parameter_units = parameter_units
self.parameter_defaults = parameter_defaults if parameter_defaults is not None else dict()
self.__cutoff_handle__ = cutoff
self.kernels = dict(kernel=kernel, kernel_gradient=kernel_gradient)
self.tag = tag
......@@ -198,17 +201,18 @@ class LocalPotentialFamily:
result : LocalPotential
The potential.
"""
kwargs = {k: (np.array(v) if isinstance(v, (list, tuple)) else v) for k, v in kwargs.items()}
parameters = self.parameter_defaults.copy()
parameters.update({k: (np.array(v) if isinstance(v, (list, tuple)) else v) for k, v in kwargs.items()})
for k, r in self.parameters.items():
if k not in kwargs:
if k not in parameters:
raise ValueError(f"Parameter '{k}' is missing from defined parameters")
if r is not None:
l, u = r
p = kwargs[k]
p = parameters[k]
if not np.all((l <= p) <= u):
raise ValueError(f"Parameter '{k}' is out of bounds: {l} <= {p} <= {u}")
return self.proto(self.coordination_number, kwargs, self.cutoff(kwargs), self.kernels, family=self, tag=tag)
return self.proto(self.coordination_number, parameters, self.cutoff(parameters), self.kernels, family=self, tag=tag)
def __call__(self, *args, **kwargs):
return self.instantiate(*args, **kwargs)
......@@ -790,6 +794,7 @@ behler5_descriptor_family = LocalPotentialFamily(
coordination_number=3,
parameters=dict(a=(0, np.inf), eta=(0, np.inf), l=(-1, 1), zeta=(0, np.inf)),
parameter_units=ParameterUnits(dict(a="angstrom", eta="1/angstrom**2")),
parameter_defaults=dict(epsilon=1),
cutoff="a",
kernel=kernel_mlsf_g5,
kernel_gradient=kernel_g_mlsf_g5,
......@@ -800,6 +805,7 @@ behler4_descriptor_family = LocalPotentialFamily(
coordination_number=3,
parameters=dict(a=(0, np.inf), eta=(0, np.inf), l=(-1, 1), zeta=(0, np.inf)),
parameter_units=ParameterUnits(dict(a="angstrom", eta="1/angstrom**2")),
parameter_defaults=dict(epsilon=1),
cutoff="a",
kernel=kernel_mlsf_g4,
kernel_gradient=kernel_g_mlsf_g4,
......
......@@ -541,25 +541,25 @@ class TestBehlerSF5(AnglePotentialTestMixin, SerializationMixin, TestCase):
cls.degenerate = False
@classmethod
def __raw_fun__(cls, r1, r2, r12_cos, a, eta, l, zeta):
return 2 ** (1 - zeta) * (1 + l * r12_cos) ** zeta * np.exp(- eta * (r1 ** 2 + r2 ** 2)) * (.5 + np.cos(np.pi * r1 / a) / 2) * (.5 + np.cos(np.pi * r2 / a) / 2)
def __raw_fun__(cls, r1, r2, r12_cos, a, eta, l, zeta, epsilon):
return epsilon * 2 ** (1 - zeta) * (1 + l * r12_cos) ** zeta * np.exp(- eta * (r1 ** 2 + r2 ** 2)) * (.5 + np.cos(np.pi * r1 / a) / 2) * (.5 + np.cos(np.pi * r2 / a) / 2)
@classmethod
def __raw_fun_prime_r1__(cls, r1, r2, r12_cos, a, eta, l, zeta):
function_value = cls.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta)
return - function_value * eta * 2 * r1 - 2 ** (1 - zeta) * (1 + l * r12_cos) ** zeta * np.exp(- eta * (r1 ** 2 + r2 ** 2)) * np.sin(np.pi * r1 / a) / 2 * (.5 + np.cos(np.pi * r2 / a) / 2) * np.pi / a
def __raw_fun_prime_r1__(cls, r1, r2, r12_cos, a, eta, l, zeta, epsilon):
function_value = cls.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta, epsilon)
return - function_value * eta * 2 * r1 - epsilon * 2 ** (1 - zeta) * (1 + l * r12_cos) ** zeta * np.exp(- eta * (r1 ** 2 + r2 ** 2)) * np.sin(np.pi * r1 / a) / 2 * (.5 + np.cos(np.pi * r2 / a) / 2) * np.pi / a
@classmethod
def __raw_fun_prime_r2__(cls, r1, r2, r12_cos, a, eta, l, zeta):
function_value = cls.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta)
return - function_value * eta * 2 * r2 - 2 ** (1 - zeta) * (1 + l * r12_cos) ** zeta * np.exp(- eta * (r1 ** 2 + r2 ** 2)) * np.sin(np.pi * r2 / a) / 2 * (.5 + np.cos(np.pi * r1 / a) / 2) * np.pi / a
def __raw_fun_prime_r2__(cls, r1, r2, r12_cos, a, eta, l, zeta, epsilon):
function_value = cls.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta, epsilon)
return - function_value * eta * 2 * r2 - epsilon * 2 ** (1 - zeta) * (1 + l * r12_cos) ** zeta * np.exp(- eta * (r1 ** 2 + r2 ** 2)) * np.sin(np.pi * r2 / a) / 2 * (.5 + np.cos(np.pi * r1 / a) / 2) * np.pi / a
@classmethod
def __raw_fun_prime_cos__(cls, r1, r2, r12_cos, a, eta, l, zeta):
return l * zeta * 2 ** (1 - zeta) * (1 + l * r12_cos) ** (zeta - 1) * np.exp(- eta * (r1 ** 2 + r2 ** 2)) * (.5 + np.cos(np.pi * r1 / a) / 2) * (.5 + np.cos(np.pi * r2 / a) / 2)
def __raw_fun_prime_cos__(cls, r1, r2, r12_cos, a, eta, l, zeta, epsilon):
return epsilon * l * zeta * 2 ** (1 - zeta) * (1 + l * r12_cos) ** (zeta - 1) * np.exp(- eta * (r1 ** 2 + r2 ** 2)) * (.5 + np.cos(np.pi * r1 / a) / 2) * (.5 + np.cos(np.pi * r2 / a) / 2)
@classmethod
def __cutoff__(cls, a, eta, l, zeta):
def __cutoff__(cls, a, eta, l, zeta, epsilon):
return a
......@@ -584,32 +584,32 @@ class TestBehlerSF4(AnglePotentialTestMixin, SerializationMixin, TestCase):
return - np.exp(- eta * r3 ** 2) * (2 * eta * r3 * (.5 + np.cos(np.pi * r3 / a) / 2) + np.pi / a * np.sin(np.pi * r3 / a) / 2) * (r3 < a)
@classmethod
def __raw_fun__(cls, r1, r2, r12_cos, a, eta, l, zeta):
return TestBehlerSF5.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta) * cls.__r3_factor__(r1, r2, r12_cos, a, eta)
def __raw_fun__(cls, r1, r2, r12_cos, a, eta, l, zeta, epsilon):
return TestBehlerSF5.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta, epsilon) * cls.__r3_factor__(r1, r2, r12_cos, a, eta)
@classmethod
def __raw_fun_prime_r1__(cls, r1, r2, r12_cos, a, eta, l, zeta):
def __raw_fun_prime_r1__(cls, r1, r2, r12_cos, a, eta, l, zeta, epsilon):
r3 = cls.__r3__(r1, r2, r12_cos)
dr3_dr1 = (r1 - r2 * r12_cos) / r3
return TestBehlerSF5.__raw_fun_prime_r1__(r1, r2, r12_cos, a, eta, l, zeta) * cls.__r3_factor__(r1, r2, r12_cos, a, eta) + \
TestBehlerSF5.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta) * cls.__r3_factor_prime__(r1, r2, r12_cos, a, eta) * dr3_dr1
return TestBehlerSF5.__raw_fun_prime_r1__(r1, r2, r12_cos, a, eta, l, zeta, epsilon) * cls.__r3_factor__(r1, r2, r12_cos, a, eta) + \
TestBehlerSF5.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta, epsilon) * cls.__r3_factor_prime__(r1, r2, r12_cos, a, eta) * dr3_dr1
@classmethod
def __raw_fun_prime_r2__(cls, r1, r2, r12_cos, a, eta, l, zeta):
def __raw_fun_prime_r2__(cls, r1, r2, r12_cos, a, eta, l, zeta, epsilon):
r3 = cls.__r3__(r1, r2, r12_cos)
dr3_dr2 = (r2 - r1 * r12_cos) / r3
return TestBehlerSF5.__raw_fun_prime_r2__(r1, r2, r12_cos, a, eta, l, zeta) * cls.__r3_factor__(r1, r2, r12_cos, a, eta) + \
TestBehlerSF5.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta) * cls.__r3_factor_prime__(r1, r2, r12_cos, a, eta) * dr3_dr2
return TestBehlerSF5.__raw_fun_prime_r2__(r1, r2, r12_cos, a, eta, l, zeta, epsilon) * cls.__r3_factor__(r1, r2, r12_cos, a, eta) + \
TestBehlerSF5.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta, epsilon) * cls.__r3_factor_prime__(r1, r2, r12_cos, a, eta) * dr3_dr2
@classmethod
def __raw_fun_prime_cos__(cls, r1, r2, r12_cos, a, eta, l, zeta):
def __raw_fun_prime_cos__(cls, r1, r2, r12_cos, a, eta, l, zeta, epsilon):
r3 = cls.__r3__(r1, r2, r12_cos)
dr3_dcos = -r1 * r2 / r3
return TestBehlerSF5.__raw_fun_prime_cos__(r1, r2, r12_cos, a, eta, l, zeta) * cls.__r3_factor__(r1, r2, r12_cos, a, eta) + \
TestBehlerSF5.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta) * cls.__r3_factor_prime__(r1, r2, r12_cos, a, eta) * dr3_dcos
return TestBehlerSF5.__raw_fun_prime_cos__(r1, r2, r12_cos, a, eta, l, zeta, epsilon) * cls.__r3_factor__(r1, r2, r12_cos, a, eta) + \
TestBehlerSF5.__raw_fun__(r1, r2, r12_cos, a, eta, l, zeta, epsilon) * cls.__r3_factor_prime__(r1, r2, r12_cos, a, eta) * dr3_dcos
@classmethod
def __cutoff__(cls, a, eta, l, zeta):
def __cutoff__(cls, a, eta, l, zeta, epsilon):
return a
......
......@@ -56,10 +56,10 @@
{
"name": "mlsf_g5",
"type": "potential-3",
"parameters": "double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles",
"parameters": "double epsilon, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles",
"preamble": [
"cdef double _fn_cutoff1, _fn_cutoff2, _fn_exponent, _fn_pw",
"cdef double _prefactor = pow(2, 1 - zeta)",
"cdef double _prefactor = pow(2, 1 - zeta) * epsilon",
"cdef int pre_compute_r_cutoff_fn_handle = pre_compute_r_handles[0], pre_compute_r_exp_fn_handle = pre_compute_r_handles[2]"
],
"preamble_grad": [
......@@ -83,10 +83,10 @@
{
"name": "mlsf_g4",
"type": "potential-3",
"parameters": "double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles",
"parameters": "double epsilon, double eta, double l, double zeta, double[:, ::1] pre_compute_r, int[::1] pre_compute_r_handles",
"preamble": [
"cdef double _fn_cutoff1, _fn_cutoff2, _fn_cutoff3, _fn_exponent, _fn_pw",
"cdef double _prefactor = pow(2, 1 - zeta)",
"cdef double _prefactor = pow(2, 1 - zeta) * epsilon",
"cdef int pre_compute_r_cutoff_fn_handle = pre_compute_r_handles[0], pre_compute_r_exp_fn_handle = pre_compute_r_handles[2]",
"cdef double _r3, _r3_factor, g5_fun"
],
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment