From 47af137d8face9f669e010f12829d1cb6b445fe2 Mon Sep 17 00:00:00 2001 From: jpsferreira Date: Fri, 22 Mar 2024 11:36:34 +0000 Subject: [PATCH] Update Kinematics class and add new tests --- hyper_surrogate/kinematics.py | 25 ++++---- hyper_surrogate/symbolic.py | 44 ++++++++++---- tests/test_kinematics.py | 32 ++++++++++- tests/test_symbolic_rules.py | 105 +++++++++++++++++++--------------- 4 files changed, 134 insertions(+), 72 deletions(-) diff --git a/hyper_surrogate/kinematics.py b/hyper_surrogate/kinematics.py index f94e4ae..ce0443d 100644 --- a/hyper_surrogate/kinematics.py +++ b/hyper_surrogate/kinematics.py @@ -68,41 +68,42 @@ def left_cauchy_green(f: np.ndarray) -> Any: # and 'jk' are the indices for the second matrix (transposed to 'kj' for multiplication). return np.einsum("nij,nkj->nik", f, f) - def strain_tensor(self, f: np.ndarray) -> Any: + @staticmethod + def strain_tensor(f: np.ndarray) -> Any: """ Compute the strain tensor. Args: - f (np.ndarray): The deformation gradient. + f (np.ndarray): The deformation gradients. batched with shape (N, 3, 3). Returns: - np.ndarray: The strain tensor. + np.ndarray: The strain tensors. batched with shape (N, 3, 3). """ - return 0.5 * (f.T @ f - np.eye(3)) + return 0.5 * (np.einsum("nji,njk->nik", f, f) - np.eye(3)) def stretch_tensor(self, f: np.ndarray) -> Any: """ Compute the stretch tensor. Args: - f (np.ndarray): The deformation gradient. + f (np.ndarray): The deformation gradients. batched with shape (N, 3, 3). Returns: - np.ndarray: The stretch tensor. + np.ndarray: The stretch tensor. batched with shape (N, 3, 3). """ - return np.sqrt(self.right_cauchy_green(f)) + return np.sqrt(np.einsum("nji,njk->nik", f, f)) - def rotation_tensor(self, f: np.ndarray) -> np.ndarray: + def rotation_tensor(self, f: np.ndarray) -> Any: """ - Compute the rotation tensor. + Compute the rotation tensors. Args: - f (np.ndarray): The deformation gradient. + f (np.ndarray): The deformation gradients. batched with shape (N, 3, 3). Returns: - np.ndarray: The rotation tensor. + np.ndarray: The rotation tensors. batched with shape (N, 3, 3). """ - return f @ np.linalg.inv(self.stretch_tensor(f)) + return np.einsum("nij,njk->nik", f, np.linalg.inv(f)) def principal_stretches(self, f: np.ndarray) -> np.ndarray: """ diff --git a/hyper_surrogate/symbolic.py b/hyper_surrogate/symbolic.py index 6fe37e5..efd7182 100644 --- a/hyper_surrogate/symbolic.py +++ b/hyper_surrogate/symbolic.py @@ -25,7 +25,7 @@ def _c_tensor(self) -> sym.Matrix: return sym.Matrix(3, 3, lambda i, j: sym.Symbol(f"C_{i+1}{j+1}")) # multuply c_tensor by itself - def _c_tensor_squared(self) -> sym.Matrix: + def _c_tensor_squared(self) -> Any: """ Compute the square of the c_tensor. @@ -68,7 +68,7 @@ def invariant3(self) -> Any: """ return self.c_tensor.det() - def pk2_tensor(self, sef: sym.Expr) -> sym.MutableDenseNDimArray: + def pk2_tensor(self, sef: sym.Expr) -> sym.Matrix: """ Compute the pk2 tensor. @@ -78,9 +78,9 @@ def pk2_tensor(self, sef: sym.Expr) -> sym.MutableDenseNDimArray: Returns: sym.Matrix: The pk2 tensor. """ - return sym.MutableDenseNDimArray([[sym.diff(sef, self.c_tensor[i, j]) for j in range(3)] for i in range(3)]) + return sym.Matrix([[sym.diff(sef, self.c_tensor[i, j]) for j in range(3)] for i in range(3)]) - def cmat_tensor(self, pk2: sym.Matrix) -> sym.MutableDenseNDimArray: + def cmat_tensor(self, pk2: sym.Matrix) -> sym.ImmutableDenseNDimArray: """ Compute the cmat tensor. @@ -90,27 +90,51 @@ def cmat_tensor(self, pk2: sym.Matrix) -> sym.MutableDenseNDimArray: Returns: sym.MutableDenseNDimArray: The cmat tensor. """ - return sym.MutableDenseNDimArray( + return sym.ImmutableDenseNDimArray( [ [[[sym.diff(pk2[i, j], self.c_tensor[k, ll]) for ll in range(3)] for k in range(3)] for j in range(3)] for i in range(3) ] ) - def substitute(self, symbolic_tensor: sym.Expr, numerical_c_tensor: np.ndarray) -> Any: + def substitute(self, symbolic_tensor: sym.Expr, numerical_c_tensor: np.ndarray, *args: dict) -> Any: """ Automatically substitute numerical values from a given 3x3 numerical matrix into c_tensor. Args: symbolic_tensor (sym.Matrix): A symbolic tensor to substitute numerical values into. numerical_matrix (np.ndarray): A 3x3 numerical matrix to substitute into c_tensor. + args (dict): Additional substitution dictionaries. Returns: - sym.Matrix: The c_tensor with numerical values substituted. + sym.Matrix: The symbolic_tensor with numerical values substituted. Raises: ValueError: If numerical_tensor is not a 3x3 matrix. """ - return symbolic_tensor.subs( - {self.c_tensor[i, j]: np.array(numerical_c_tensor)[i, j] for i in range(3) for j in range(3)} - ) + if not isinstance(numerical_c_tensor, np.ndarray) or numerical_c_tensor.shape != (3, 3): + raise ValueError("c_tensor.shape") + + # Start with substitutions for c_tensor elements + substitutions = {self.c_tensor[i, j]: numerical_c_tensor[i, j] for i in range(3) for j in range(3)} + # Merge additional substitution dictionaries from *args + substitutions.update(*args) + return symbolic_tensor.subs(substitutions) + + def substitute_iterator(self, symbolic_tensor: sym.Expr, numerical_c_tensors: np.ndarray, *args: dict) -> Any: + """ + Automatically substitute numerical values from a given 3x3 numerical matrix into c_tensor. + + Args: + symbolic_tensor (sym.Matrix): A symbolic tensor to substitute numerical values into. + numerical_matrix (np.ndarray): A 3x3 numerical matrix to substitute into c_tensor. + args (dict): Additional substitution dictionaries. + + Returns: + sym.Matrix: The symbolic_tensor with numerical values substituted. + + Raises: + ValueError: If numerical_tensor is not a 3x3 matrix. + """ + for numerical_c_tensor in numerical_c_tensors: + yield self.substitute(symbolic_tensor, numerical_c_tensor, *args) diff --git a/tests/test_kinematics.py b/tests/test_kinematics.py index e6f0e9c..ea2717c 100644 --- a/tests/test_kinematics.py +++ b/tests/test_kinematics.py @@ -2,9 +2,14 @@ import pytest from hyper_surrogate.deformation_gradient import DeformationGradientGenerator as FGen -from hyper_surrogate.kinematics import Kinematics as K +from hyper_surrogate.kinematics import Kinematics -SIZE = 2 +SIZE = 1 + + +@pytest.fixture +def K(): + return Kinematics() @pytest.fixture @@ -18,6 +23,27 @@ def right_cauchys(def_gradients): return np.array([np.matmul(f, f.T) for f in def_gradients]) -def test_right_cauchys(def_gradients): +@pytest.fixture +def left_cauchys(def_gradients): + return np.array([np.matmul(f.T, f) for f in def_gradients]) + + +def test_right_cauchys(def_gradients, K): right_cauchys = K.right_cauchy_green(def_gradients) assert np.allclose(right_cauchys, np.array([np.matmul(f.T, f) for f in def_gradients])) + + +def test_left_cauchys(def_gradients, K): + left_cauchys = K.left_cauchy_green(def_gradients) + assert np.allclose(left_cauchys, np.array([np.matmul(f, f.T) for f in def_gradients])) + + +def test_strain_tensor(def_gradients, K): + strain_tensors = K.strain_tensor(def_gradients) + assert np.allclose(strain_tensors, 0.5 * (np.array([f.T @ f for f in def_gradients]) - np.eye(3))) + + +# def test_stretch_tensor(def_gradients, K): +# stretch_tensors = K.stretch_tensor(def_gradients) +# logging.info(stretch_tensors) +# assert np.allclose(stretch_tensors, np.sqrt(np.array([f.T @ f for f in def_gradients]))) diff --git a/tests/test_symbolic_rules.py b/tests/test_symbolic_rules.py index 56d930b..e690c56 100644 --- a/tests/test_symbolic_rules.py +++ b/tests/test_symbolic_rules.py @@ -1,5 +1,3 @@ -import logging - import numpy as np import pytest import sympy as sym @@ -11,41 +9,59 @@ SIZE = 2 +# numeric fixtures @pytest.fixture -def def_gradients(): +def def_gradients() -> np.ndarray: return FGen(seed=42, size=SIZE).generate() @pytest.fixture -def handler(): +def right_cauchys(def_gradients) -> np.ndarray: + return K.right_cauchy_green(def_gradients) + + +# symbolic fixtures +@pytest.fixture +def handler() -> SymbolicHandler: return SymbolicHandler() +# mooney_rivlin @pytest.fixture -def sef(): - symbolic_handler = SymbolicHandler() - return (symbolic_handler.invariant1 - 3) * sym.Symbol("C10") + (symbolic_handler.invariant2 - 3) * sym.Symbol("C01") +def sef(handler) -> sym.Expr: + return (handler.invariant1 - 3) * sym.Symbol("C10") + (handler.invariant2 - 3) * sym.Symbol("C01") -# kinematics testing @pytest.fixture -def right_cauchys(def_gradients): - return np.array([np.matmul(f, f.T) for f in def_gradients]) +def sef_args() -> dict: + return {"C10": 1, "C01": 1} + + +@pytest.fixture +def pk2(handler, sef) -> sym.Matrix: + return handler.pk2_tensor(sef) + + +@pytest.fixture +def cmat(handler, pk2) -> sym.ImmutableDenseNDimArray: + return handler.cmat_tensor(pk2) + +# testing -def test_symbolic_pk2_cmat(handler, sef): + +def test_symbolic_pk2_cmat(pk2, cmat): # derivative of sef in order to c_tensor - pk2 = handler.pk2_tensor(sef) - cmat = handler.cmat_tensor(pk2) + # assert instance - assert isinstance(pk2, sym.MutableDenseNDimArray) - assert isinstance(cmat, sym.MutableDenseNDimArray) + assert isinstance(pk2, sym.Matrix) + assert isinstance(cmat, sym.ImmutableDenseNDimArray) # assert shape assert pk2.shape == (3, 3) assert cmat.shape == (3, 3, 3, 3) -def test_symbolic_subs(handler): +def test_symbolic_subs_in_c(handler): # substitute the c_tensor with a 3x3 matrix of ones c_tensor = handler.c_tensor # dummy c values @@ -56,34 +72,29 @@ def test_symbolic_subs(handler): assert all(c_tensor_subs[i, j] == 1 for i in range(3) for j in range(3)) -def test_symbolic_subs_in_batch(handler, def_gradients): - # - c_tensor = handler.c_tensor - # get first c from def_gradients - c = K.right_cauchy_green(def_gradients) - # subs c_tensor with c values - c_all = [handler.substitute(c_tensor, cc) for cc in c] - logging.info(c_all) - - -# voigt notation -# def test_voigt_notation(handler): -# c_tensor = handler.c_tensor -# c_voigt = sym.Matrix([c_tensor[0, 0], c_tensor[1, 1], c_tensor[2, 2], c_tensor[0, 1], c_tensor[0, 2], c_tensor[1, 2]]) - -# invariant3 = symbolic_handler.invariant3_symbolic() - - -# def test_derivative(def_gradients,c_tensor_symbolic): -# logging.info(def_gradients) - -# C_10_sym=sym.Symbol("C_10_sym") -# #get symbols from c_tensor_symbolic -# C11,C12,C13,C21,C22,C23,C31,C32,C33=c_tensor_symbolic.args -# I3=C11*C22*C33 - C11*C23*C32 - C12*C21*C33 + C12*C23*C31 + C13*C21*C32 - C13*C22*C31 -# trace=C11+C22+C33 -# I1=trace*(I3)**(-1/3) - -# #Symbolic Strain Energy Function -# SEF=(I1-3)*C_10_sym -# logging.info(SEF) +def test_symbolic_subs_in_pk2(handler, pk2, right_cauchys, sef_args): + # right_cauchys # (N, 3, 3) + # for each c_tensor in c, substitute the pk2 tensor with c values and material parameters values. + assert all( + isinstance(subs, sym.Matrix) + and subs.shape == (3, 3) + and all(isinstance(subs[i, j], sym.Expr) for i in range(3) for j in range(3)) + for subs in handler.substitute_iterator(pk2, right_cauchys, sef_args) + ) + + +def test_symbolic_subs_in_cmat(handler, cmat, right_cauchys, sef_args): + # right_cauchys # (N, 3, 3) + # for each c_tensor in c, substitute the cmat tensor with c values and material parameters values. + assert all( + isinstance(subs, sym.ImmutableDenseNDimArray) + and subs.shape == (3, 3, 3, 3) + and all( + isinstance(subs[i, j, k, ll], sym.Expr) + for i in range(3) + for j in range(3) + for k in range(3) + for ll in range(3) + ) + for subs in handler.substitute_iterator(cmat, right_cauchys, sef_args) + )