Skip to content

Commit

Permalink
MAINT: PEP8 compliancy
Browse files Browse the repository at this point in the history
  • Loading branch information
@Fraser-Birks
Fraser-Birks committed Nov 6, 2023
1 parent 3931d00 commit 6d6463c
Show file tree
Hide file tree
Showing 7 changed files with 1,328 additions and 1,233 deletions.
94 changes: 46 additions & 48 deletions matscipy/cauchy_born.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,8 @@
from scipy.stats.qmc import LatinHypercube, scale

import ase


class RegressionModel:
"""Simple regression model object wrapper for np.linalg.lstsq for
predicting shifts."""
Expand All @@ -27,7 +29,7 @@ def fit(self, phi, values):
values : array_like
value of function to be fitted evaluated at each data point
"""
self.model = np.linalg.lstsq(phi, values,rcond=None)[0]
self.model = np.linalg.lstsq(phi, values, rcond=None)[0]

def predict(self, phi):
"""Evaluate the simple least-squares regression model.
Expand All @@ -46,8 +48,8 @@ def predict(self, phi):
"""

return phi @ self.model
def predict_gradient(self,grad_phi):

def predict_gradient(self, grad_phi):
"""Evaluate the simple least-squares regression model using
the derivative of the basis functions to get the gradient.
Expand All @@ -64,11 +66,11 @@ def predict_gradient(self,grad_phi):
1D vector of least squares model gradient predictions for simple model.
"""
return grad_phi @ self.model

def save(self):
"""Saves the regression model to file: CB_model.npy"""
np.save('CB_model.npy',self.model)
np.save('CB_model.npy', self.model)

def load(self):
"""Loads the regression model from file: CB_model.npy"""
self.model = np.load('CB_model.npy')
Expand Down Expand Up @@ -104,7 +106,7 @@ def __init__(self, el, a0, calc, lattice=Diamond):
self.lattice1mask = None # mask for the lattice 1 atoms
self.lattice2mask = None # mask for the lattice 2 atoms

def set_sublattices(self, atoms, A,read_from_atoms=False):
def set_sublattices(self, atoms, A, read_from_atoms=False):
"""Apply a small strain to all atoms in the supplied atoms structure
and determine which atoms belong to which sublattice using forces. NOTE
as this method is based on forces, it does not work in cases where atom
Expand Down Expand Up @@ -181,11 +183,13 @@ def set_sublattices(self, atoms, A,read_from_atoms=False):
if all(element == lattice1mask[0] for element in lattice1mask):
lattice1mask = force_diff_lattice[:, 2] > 0
lattice2mask = force_diff_lattice[:, 2] < 0

#set new array in atoms
try: #make new arrays if they didn't already exist
atoms.new_array('lattice1mask',np.zeros(len(atoms),dtype=bool))
atoms.new_array('lattice2mask',np.zeros(len(atoms),dtype=bool))

# set new array in atoms
try: # make new arrays if they didn't already exist
atoms.new_array('lattice1mask', np.zeros(
len(atoms), dtype=bool))
atoms.new_array('lattice2mask', np.zeros(
len(atoms), dtype=bool))
except RuntimeError:
pass
lattice1maskatoms = atoms.arrays['lattice1mask']
Expand All @@ -197,7 +201,7 @@ def set_sublattices(self, atoms, A,read_from_atoms=False):
self.lattice1mask = lattice1mask
self.lattice2mask = lattice2mask

def switch_sublattices(self,atoms):
def switch_sublattices(self, atoms):
"""Switch the sublattice masks on a set of atoms around
such that lattice1mask = lattice2mask and vice versa
Expand Down Expand Up @@ -276,12 +280,14 @@ def f_gl(E_vec, calc, unitcell):
cell = unitcell_copy.get_cell()
cell_rescale = np.transpose(x[:, :] @ cell)
unitcell_copy.set_cell(cell_rescale, scale_atoms=True)
initial_shift[:] = unitcell_copy.get_positions()[1] - unitcell_copy.get_positions()[0]
initial_shift[:] = unitcell_copy.get_positions(
)[1] - unitcell_copy.get_positions()[0]

# relax cell
opt = LBFGS(unitcell_copy, logfile=None)
opt.run(fmax=1e-10)
relaxed_shift[:] = unitcell_copy.get_positions()[1] - unitcell_copy.get_positions()[0]
relaxed_shift[:] = unitcell_copy.get_positions(
)[1] - unitcell_copy.get_positions()[0]

# get shift
shift_diff = relaxed_shift - initial_shift
Expand Down Expand Up @@ -800,9 +806,8 @@ def predict_shifts(
for i in range(natoms):
shifts[i, :] = np.transpose(A) @ R[i, :, :] @ shifts_no_rr[i, :]
return shifts


def permutation(self,strain_vec, perm_shift):

def permutation(self, strain_vec, perm_shift):
"""Rotate a Green-Lagrange strain tensor in Voigt notation about
the 111 axis a number of times given by perm_shift. This
corresponds to permuting the top three components and bottom three
Expand Down Expand Up @@ -859,8 +864,6 @@ def evaluate_shift_model(self, E, method='regression'):
# takes E in the form of a list of matrices [natoms, :, :]
# define function for evaluating model



def evaluate_shift_regression(eps):
"""Evaluate the regression model a given set of Green-Lagrange
strain vectors.
Expand Down Expand Up @@ -1399,11 +1402,10 @@ def basis_function_evaluation(self, E_vecs):

# return the design matrix

#print('design matrix size',np.shape(np.concatenate((E_vecs, phi_2nd_term), axis=1)))
# print('design matrix size',np.shape(np.concatenate((E_vecs, phi_2nd_term), axis=1)))
return np.concatenate((E_vecs, phi_2nd_term), axis=1)


def grad_basis_function_evaluation(self,E_vecs,dE_vecs):
def grad_basis_function_evaluation(self, E_vecs, dE_vecs):
"""Function that evaluates the derivative of the polynomial
basis functions used in the regression model on a set of
strain vectors to build the gradient design matrix grad_phi.
Expand All @@ -1422,25 +1424,23 @@ def grad_basis_function_evaluation(self,E_vecs,dE_vecs):
Design matrix composed of all of the gradient basis functions evaluated at
every data point given in E_vecs.
"""
#take a list of component wise gradients strain with respect
#to some other variable (for example, crack tip position, alpha)
#and propagate this through the regression model to get a set of
#derivatives of the shift with respect to the external variable
# take a list of component wise gradients strain with respect
# to some other variable (for example, crack tip position, alpha)
# and propagate this through the regression model to get a set of
# derivatives of the shift with respect to the external variable

#take dE_vecs in the form of [natoms,6]
# take dE_vecs in the form of [natoms,6]
triu_indices = np.triu_indices(6)
grad_phi_2nd_term = np.zeros(
[np.shape(E_vecs)[0], np.shape(triu_indices)[1]])
for i, dE in enumerate(dE_vecs):
e_de_outer_p = np.outer(E_vecs[i,:],dE)
e_de_outer_p = np.outer(E_vecs[i, :], dE)
de_dalpha = e_de_outer_p + np.transpose(e_de_outer_p)
grad_phi_2nd_term[i,:] = de_dalpha[triu_indices]

return np.concatenate((dE_vecs, grad_phi_2nd_term),axis=1)

grad_phi_2nd_term[i, :] = de_dalpha[triu_indices]

return np.concatenate((dE_vecs, grad_phi_2nd_term), axis=1)

def evaluate_shift_gradient_regression(self,E,dE):
def evaluate_shift_gradient_regression(self, E, dE):
"""Evaluate the gradient of the regression model a given set of Green-Lagrange
strain tensors and their gradients.
Expand All @@ -1449,7 +1449,7 @@ def evaluate_shift_gradient_regression(self,E,dE):
E : array_like
2D array of size [natoms,6], where each row is a Green-Lagrange
strain tensor written in Voigt notation.
dE : array_like
2D array of size [natoms, 6], where each row is the derivative of the
corresponding row in E with respect to some parameter
Expand All @@ -1460,7 +1460,7 @@ def evaluate_shift_gradient_regression(self,E,dE):
prediction of gradient of the Cauchy-Born shift corrector model
at each of the given E points.
"""

# turn E and dE into voigt vectors
E_voigt = np.zeros([np.shape(E)[0], 6])
E_voigt[:, 0] = E[:, 0, 0]
Expand All @@ -1469,7 +1469,7 @@ def evaluate_shift_gradient_regression(self,E,dE):
E_voigt[:, 3] = E[:, 1, 2]
E_voigt[:, 4] = E[:, 0, 2]
E_voigt[:, 5] = E[:, 0, 1]

dE_voigt = np.zeros([np.shape(dE)[0], 6])
dE_voigt[:, 0] = dE[:, 0, 0]
dE_voigt[:, 1] = dE[:, 1, 1]
Expand All @@ -1490,14 +1490,13 @@ def evaluate_shift_gradient_regression(self,E,dE):

# Predict the shift gradients using the fitted regression model
predictions[:, 0] = self.RM.predict_gradient(
self.grad_basis_function_evaluation(epsx,depsx))
self.grad_basis_function_evaluation(epsx, depsx))
predictions[:, 1] = self.RM.predict_gradient(
self.grad_basis_function_evaluation(epsy,depsy))
self.grad_basis_function_evaluation(epsy, depsy))
predictions[:, 2] = self.RM.predict_gradient(
self.grad_basis_function_evaluation(epsz,depsz))
self.grad_basis_function_evaluation(epsz, depsz))
return predictions


def get_shift_gradients(
self,
A,
Expand Down Expand Up @@ -1570,9 +1569,9 @@ def get_shift_gradients(
A, atoms, F_func=F_func, E_func=E_func,
coordinates=coordinates, de=de, *args, **kwargs)

#print(E_higher,E_lower)
# print(E_higher,E_lower)
dE = (E_higher-E_lower)/(2*de)
#print(dE)
# print(dE)
natoms = len(atoms)
# get the cauchy born shifts unrotated
dshifts_no_rr = self.evaluate_shift_gradient_regression(E, dE)
Expand All @@ -1582,14 +1581,13 @@ def get_shift_gradients(
# and to get them back into the lab frame
for i in range(natoms):
dshifts[i, :] = np.transpose(A) @ dshifts_no_rr[i, :]
#dshift_2[i, :] = np.transpose(A) @ dshift_2_no_rr[i, :]
# dshift_2[i, :] = np.transpose(A) @ dshift_2_no_rr[i, :]

# need to adjust gradients for different lattices
dshifts[self.lattice1mask] = -0.5*(dshifts[self.lattice1mask])
dshifts[self.lattice2mask] = 0.5*(dshifts[self.lattice2mask])

return dshifts

return dshifts

def save_regression_model(self):
"""Saves the regression model to file: CB_model.npy"""
Expand All @@ -1602,6 +1600,6 @@ def load_regression_model(self):

def get_model(self):
return self.RM.model

def set_model(self, model):
self.RM.model = model
self.RM.model = model
Loading

0 comments on commit 6d6463c

Please sign in to comment.