Commit 752fd85b authored by Malthe Kjær Bisbo's avatar Malthe Kjær Bisbo
Browse files

in GPR.predict_forces changed eval_std to eval_with_energy_std

parent 05fabfb7
from candidate_operations.candidate_generation import CandidateGenerator, StartGenerator
from candidate_operations.basic_mutations import RattleMutation, PermutationMutation
__all__ = ['CandidateGenerator', 'StartGenerator', 'RattleMutation', 'PermutationMutation']
from candidates.candidate_generation import CandidateGenerator, StartGenerator
from candidates.basic_mutations import RattleMutation, PermutationMutation
__all__ = ['CandidateGenerator', 'StartGenerator', 'RattleMutation', 'RattleMutation2', 'PermutationMutation']
......@@ -209,61 +209,6 @@ def random_pos(box):
pos += vspan[i] * r[0, i]
return pos
class OperationConstraint():
""" Class used to enforce constraints on the positions of
atoms in mutation and crossover operations.
Parameters:
box: Box in which atoms are allowed to be placed. It should
have the form [] [p0, vspan] where 'p0' is the position of
the box corner and 'vspan' is a matrix containing the three
spanning vectors.
xlim: On the form [xmin, xmax], specifying, in the x-direction,
the lower and upper limit of the region atoms can be moved
within.
ylim, zlim: See xlim.
"""
def __init__(self, box=None, xlim=None, ylim=None, zlim=None):
self.box = box
self.xlim = xlim
self.ylim = ylim
self.zlim = zlim
def check_if_valid(self, positions):
""" Returns whether positions are valid under the
constraints or not.
"""
if np.ndim(positions) == 1:
pos = positions.reshape(-1,3)
else:
pos = positions
if self.box is not None:
p0, V = self.box
p_rel = pos - p0 # positions relative to box anchor.
V_inv = np.linalg.inv(V)
p_box = p_rel @ V_inv # positions in box-vector basis.
if (np.any(p_box < 0) or np.any(p_box > 1)):
return False
if self.xlim is not None:
if (np.any(pos[:,0] < self.xlim[0]) or
np.any(pos[:,0] > self.xlim[1])):
return False
if self.ylim is not None:
if (np.any(pos[:,1] < self.ylim[0]) or
np.any(pos[:,1] > self.ylim[1])):
return False
if self.zlim is not None:
if (np.any(pos[:,2] < self.zlim[0]) or
np.any(pos[:,2] > self.zlim[1])):
return False
return True
class StartGenerator(OffspringOperation):
""" Class used to generate random initial candidates.
......
......@@ -105,6 +105,18 @@ class GPR():
self.memory = gpr_memory(self.descriptor, self.prior)
def predict_energy(self, a, eval_std=False):
"""Evaluate the energy predicted by the GPR-model.
parameters:
a: Atoms object
The structure to evaluate.
eval_std: bool
In addition to the force, predict also force contribution
arrising from including the standard deviation of the
predicted energy.
"""
x = self.descriptor.get_feature(a)
k = self.kernel.kernel_vector(x, self.X)
......@@ -118,7 +130,19 @@ class GPR():
else:
return E
def predict_forces(self, a, eval_std=False):
def predict_forces(self, a, eval_with_energy_std=False):
"""Evaluate the force predicted by the GPR-model.
parameters:
a: Atoms object
The structure to evaluate.
eval_with_energy_std: bool
In addition to the force, predict also force contribution
arrising from including the standard deviation of the
predicted energy.
"""
# Calculate descriptor and its gradient
x = self.descriptor.get_feature(a)
......@@ -130,7 +154,7 @@ class GPR():
F = -np.dot(k_ddr.T, self.alpha) + self.prior.forces(a)
if eval_std:
if eval_with_energy_std:
k = self.kernel.kernel_vector(x, self.X)
vk = np.dot(self.K_inv, k)
g = self.K0 - np.dot(k.T, vk)
......
from ase.calculators.calculator import Calculator
from ase.calculators.singlepoint import SinglePointCalculator
class gpr_calculator(Calculator):
......@@ -26,6 +25,6 @@ class gpr_calculator(Calculator):
if self.kappa is None:
F = self.gpr.predict_forces(atoms)
else:
F, Fstd = self.gpr.predict_forces(atoms, eval_std=True)
F, Fstd = self.gpr.predict_forces(atoms, eval_with_energy_std=True)
F = F - self.kappa*Fstd
self.results['forces'] = F
......@@ -89,4 +89,58 @@ def array_to_string(arr, unit='', format='0.4f', max_line_length=80):
msg += '\n'
line_length_counter = len(string)
msg += string
return msg
\ No newline at end of file
return msg
class OperationConstraint():
""" Class used to enforce constraints on the positions of
atoms in mutation and crossover operations.
Parameters:
box: Box in which atoms are allowed to be placed. It should
have the form [] [p0, vspan] where 'p0' is the position of
the box corner and 'vspan' is a matrix containing the three
spanning vectors.
xlim: On the form [xmin, xmax], specifying, in the x-direction,
the lower and upper limit of the region atoms can be moved
within.
ylim, zlim: See xlim.
"""
def __init__(self, box=None, xlim=None, ylim=None, zlim=None):
self.box = box
self.xlim = xlim
self.ylim = ylim
self.zlim = zlim
def check_if_valid(self, positions):
""" Returns whether positions are valid under the
constraints or not.
"""
if np.ndim(positions) == 1:
pos = positions.reshape(-1,3)
else:
pos = positions
if self.box is not None:
p0, V = self.box
p_rel = pos - p0 # positions relative to box anchor.
V_inv = np.linalg.inv(V)
p_box = p_rel @ V_inv # positions in box-vector basis.
if (np.any(p_box < 0) or np.any(p_box > 1)):
return False
if self.xlim is not None:
if (np.any(pos[:,0] < self.xlim[0]) or
np.any(pos[:,0] > self.xlim[1])):
return False
if self.ylim is not None:
if (np.any(pos[:,1] < self.ylim[0]) or
np.any(pos[:,1] > self.ylim[1])):
return False
if self.zlim is not None:
if (np.any(pos[:,2] < self.zlim[0]) or
np.any(pos[:,2] > self.zlim[1])):
return False
return True
\ No newline at end of file
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