Source code for sfepy.terms.terms_fibres

from __future__ import absolute_import
import numpy as nm

from sfepy.base.base import Struct
from sfepy.terms.terms_hyperelastic_tl import HyperElasticTLBase
from sfepy.mechanics.tensors import dim2sym
from sfepy.homogenization.utils import iter_sym
from six.moves import range

[docs] def create_omega(fdir): r""" Create the fibre direction tensor :math:`\omega_{ij} = d_i d_j`. """ n_el, n_qp, dim, _ = fdir.shape sym = dim2sym(dim) omega = nm.empty((n_el, n_qp, sym, 1), dtype=nm.float64) for ii, (ir, ic) in enumerate(iter_sym(dim)): omega[..., ii, 0] = fdir[..., ir, 0] * fdir[..., ic, 0] return omega
[docs] def compute_fibre_strain(green_strain, omega): """ Compute the Green strain projected to the fibre direction. """ eps = nm.zeros_like(omega[..., :1, :]) for ii in range(omega.shape[2]): eps[..., 0, 0] += omega[..., ii, 0] * green_strain[..., ii, 0] return eps
def _setdefault_fibre_data(self, state): """ Returns `fibre_data` :class:``Struct`` for storing/caching fibre-related data. """ cache = self.set_default('fibre_cache', {}) _, _, key = self.get_mapping(state, return_key=True) data_key = key + ('fibre_data',) if data_key in cache: fibre_data = cache[data_key] else: fibre_data = Struct() cache[data_key] = fibre_data return fibre_data
[docs] class FibresActiveTLTerm(HyperElasticTLBase): r""" Hyperelastic active fibres term. Effective stress :math:`S_{ij} = A f_{\rm max} \exp{\left\{-(\frac{\epsilon - \varepsilon_{\rm opt}}{s})^2\right\}} d_i d_j`, where :math:`\epsilon = E_{ij} d_i d_j` is the Green strain :math:`\ull{E}` projected to the fibre direction :math:`\ul{d}`. :Definition: .. math:: \int_{\Omega} S_{ij}(\ul{u}) \delta E_{ij}(\ul{u};\ul{v}) :Arguments: - material_1 : :math:`f_{\rm max}` - material_2 : :math:`\varepsilon_{\rm opt}` - material_3 : :math:`s` - material_4 : :math:`\ul{d}` - material_5 : :math:`A` - virtual : :math:`\ul{v}` - state : :math:`\ul{u}` """ name = 'dw_tl_fib_a' arg_types = ('material_1', 'material_2', 'material_3', 'material_4', 'material_5', 'virtual', 'state') arg_shapes = {'material_1' : '1, 1', 'material_2' : '1, 1', 'material_3' : '1, 1', 'material_4' : 'D, 1', 'material_5' : '1, 1', 'virtual' : ('D', 'state'), 'state' : 'D'} family_data_names = ['green_strain']
[docs] def get_fargs(self, mat1, mat2, mat3, mat4, mat5, virtual, state, mode=None, term_mode=None, diff_var=None, **kwargs): fibre_data = _setdefault_fibre_data(self, state) n_el, _, _, _, _ = self.get_data_shape(state) mat1 = HyperElasticTLBase.tile_mat(mat1, n_el) mat2 = HyperElasticTLBase.tile_mat(mat2, n_el) mat3 = HyperElasticTLBase.tile_mat(mat3, n_el) mat4 = HyperElasticTLBase.tile_mat(mat4, n_el) mat5 = HyperElasticTLBase.tile_mat(mat5, n_el) fargs = HyperElasticTLBase.get_fargs(self, (mat1, mat2, mat3, mat4, mat5), virtual, state, mode, term_mode, diff_var, fibre_data=fibre_data, **kwargs) return fargs
[docs] @staticmethod def stress_function(out, pars, green_strain, fibre_data=None): fmax, eps_opt, s, fdir, act = pars omega = fibre_data.get('omega', None) if omega is None: omega = fibre_data.omega = create_omega(fdir) eps = fibre_data.eps = compute_fibre_strain(green_strain, omega) tau = fibre_data.tau = act * fmax * nm.exp(-((eps - eps_opt) / s)**2.0) out[:] = omega * tau
[docs] @staticmethod def tan_mod_function(out, pars, green_strain, fibre_data=None): fmax, eps_opt, s, fdir, act = pars omega, eps, tau = fibre_data.omega, fibre_data.eps, fibre_data.tau for ir in range(omega.shape[2]): for ic in range(omega.shape[2]): out[..., ir, ic] = omega[..., ir, 0] * omega[..., ic, 0] out[:] *= -2.0 * ((eps - eps_opt) / (s**2.0)) * tau
[docs] def get_eval_shape(self, mat1, mat2, mat3, mat4, mat5, virtual, state, mode=None, term_mode=None, diff_var=None, **kwargs): n_el, n_qp, dim, n_en, n_c = self.get_data_shape(state) sym = dim2sym(dim) if mode != 'qp': n_qp = 1 return (n_el, n_qp, sym, 1), state.dtype
[docs] class FibresExponentialTLTerm(HyperElasticTLBase): r""" Hyperelastic fibres term with an exponential response. Effective stress :math:`S_{ij} = \max\left(0, \sigma \left[ \exp{\left\{k (\epsilon - \epsilon_0)\right\}} - 1 \right]\right) d_i d_j`, where :math:`\epsilon = E_{ij} d_i d_j` is the Green strain :math:`\ull{E}` projected to the fibre direction :math:`\ul{d}`. :Definition: .. math:: \int_{\Omega} S_{ij}(\ul{u}) \delta E_{ij}(\ul{u};\ul{v}) :Arguments: - material_1 : :math:`\sigma` - material_3 : :math:`k` - material_3 : :math:`\epsilon_{0}` - material_4 : :math:`\ul{d}` - virtual : :math:`\ul{v}` - state : :math:`\ul{u}` """ name = 'dw_tl_fib_e' arg_types = ('material_1', 'material_2', 'material_3', 'material_4', 'virtual', 'state') arg_shapes = {'material_1' : '1, 1', 'material_2' : '1, 1', 'material_3' : '1, 1', 'material_4' : 'D, 1', 'virtual' : ('D', 'state'), 'state' : 'D'} family_data_names = ['green_strain']
[docs] def get_fargs(self, mat1, mat2, mat3, mat4, virtual, state, mode=None, term_mode=None, diff_var=None, **kwargs): fibre_data = _setdefault_fibre_data(self, state) n_el, _, _, _, _ = self.get_data_shape(state) mat1 = HyperElasticTLBase.tile_mat(mat1, n_el) mat2 = HyperElasticTLBase.tile_mat(mat2, n_el) mat3 = HyperElasticTLBase.tile_mat(mat3, n_el) mat4 = HyperElasticTLBase.tile_mat(mat4, n_el) fargs = HyperElasticTLBase.get_fargs(self, (mat1, mat2, mat3, mat4), virtual, state, mode, term_mode, diff_var, fibre_data=fibre_data, **kwargs) return fargs
[docs] @staticmethod def stress_function(out, pars, green_strain, fibre_data=None): sigma, k, eps0, fdir = pars omega = fibre_data.get('omega', None) if omega is None: omega = fibre_data.omega = create_omega(fdir) eps = fibre_data.eps = compute_fibre_strain(green_strain, omega) tau = fibre_data.tau = nm.maximum( 0, sigma * (nm.exp(k * (eps - eps0)) - 1.0), ) out[:] = omega * tau
[docs] @staticmethod def tan_mod_function(out, pars, green_strain, fibre_data=None): sigma, k, eps0, fdir = pars omega, eps, tau = fibre_data.omega, fibre_data.eps, fibre_data.tau for ir in range(omega.shape[2]): for ic in range(omega.shape[2]): out[..., ir, ic] = omega[..., ir, 0] * omega[..., ic, 0] tan_mod = nm.where( tau > 0.0, sigma * k * nm.exp(k * (eps - eps0)), 0.0, ) out[:] *= tan_mod
[docs] def get_eval_shape(self, mat1, mat2, mat3, mat4, virtual, state, mode=None, term_mode=None, diff_var=None, **kwargs): n_el, n_qp, dim, n_en, n_c = self.get_data_shape(state) sym = dim2sym(dim) if mode != 'qp': n_qp = 1 return (n_el, n_qp, sym, 1), state.dtype