sfepy.terms.terms_multilinear module¶
- class sfepy.terms.terms_multilinear.ECauchyStressTerm(*args, **kwargs)[source]¶
Evaluate Cauchy stress tensor.
It is given in the usual vector form exploiting symmetry: in 3D it has 6 components with the indices ordered as [11, 22, 33, 12, 13, 23], in 2D it has 3 components with the indices ordered as [11, 22, 12].
- Definition:
\int_{\Omega} D_{ijkl} e_{kl}(\ul{w})
- Call signature:
de_cauchy_stress
(material, parameter)
- Arguments:
material : D_{ijkl}
parameter : \ul{w}
- arg_shapes = {'material': 'S, S', 'parameter': 'D'}¶
- arg_types = ('material', 'parameter')¶
- name = 'de_cauchy_stress'¶
- class sfepy.terms.terms_multilinear.EConvectTerm(*args, **kwargs)[source]¶
Nonlinear convective term.
- Definition:
\int_{\Omega} ((\ul{u} \cdot \nabla) \ul{u}) \cdot \ul{v} \mbox{ , } \int_{\Omega} ((\ul{w} \cdot \nabla) \ul{w}) \cdot \bar{\ul{u}}
- Call signature:
de_convect
(virtual, state)
(parameter_1, parameter_2)
- Arguments 1:
virtual : \ul{v}
state : \ul{u}
- Arguments 2:
parameter_1 : \bar{\ul{u}}
parameter_2 : \ul{w}
- arg_shapes = {'parameter_1': 'D', 'parameter_2': 'D', 'state': 'D', 'virtual': ('D', 'state')}¶
- arg_types = (('virtual', 'state'), ('parameter_1', 'parameter_2'))¶
- modes = ('weak', 'eval')¶
- name = 'de_convect'¶
- class sfepy.terms.terms_multilinear.EDivGradTerm(*args, **kwargs)[source]¶
Vector field diffusion term.
- Definition:
\int_{\Omega} \nu\ \nabla \ul{v} : \nabla \ul{u} \mbox{ , } \int_{\Omega} \nu\ \nabla \ul{u} : \nabla \ul{w} \\ \int_{\Omega} \nabla \ul{v} : \nabla \ul{u} \mbox{ , } \int_{\Omega} \nabla \ul{u} : \nabla \ul{w}
- Call signature:
de_div_grad
(opt_material, virtual, state)
(opt_material, parameter_1, parameter_2)
- Arguments 1:
material : \nu (viscosity, optional)
virtual : \ul{v}
state : \ul{u}
- Arguments 2:
material : \nu (viscosity, optional)
parameter_1 : \ul{u}
parameter_2 : \ul{w}
- arg_shapes = [{'opt_material': '1, 1', 'parameter_1': 'D', 'parameter_2': 'D', 'state': 'D', 'virtual': ('D', 'state')}, {'opt_material': None}]¶
- arg_types = (('opt_material', 'virtual', 'state'), ('opt_material', 'parameter_1', 'parameter_2'))¶
- modes = ('weak', 'eval')¶
- name = 'de_div_grad'¶
- class sfepy.terms.terms_multilinear.EDivTerm(*args, **kwargs)[source]¶
Weighted divergence term.
- Definition:
\int_{\Omega} \nabla \cdot \ul{v} \mbox { , } \int_{\Omega} \nabla \cdot \ul{u} \\ \int_{\Omega} c \nabla \cdot \ul{v} \mbox { , } \int_{\Omega} c \nabla \cdot \ul{u}
- Call signature:
de_div
(opt_material, virtual)
(opt_material, parameter)
- Arguments 1:
material : c (optional)
virtual : \ul{v}
- Arguments 2:
material : c (optional)
parameter : \ul{u}
- arg_shapes = [{'opt_material': '1, 1', 'parameter': 'D', 'virtual': ('D', None)}, {'opt_material': None}]¶
- arg_types = (('opt_material', 'virtual'), ('opt_material', 'parameter'))¶
- modes = ('weak', 'eval')¶
- name = 'de_div'¶
- class sfepy.terms.terms_multilinear.EDotTerm(*args, **kwargs)[source]¶
Volume and surface L^2(\Omega) weighted dot product for both scalar and vector fields. Can be evaluated. Can use derivatives.
- Definition:
\int_{\cal{D}} q p \mbox{ , } \int_{\cal{D}} \ul{v} \cdot \ul{u} \mbox{ , } \int_{\cal{D}} p r \mbox{ , } \int_{\cal{D}} \ul{u} \cdot \ul{w} \\ \int_{\cal{D}} c q p \mbox{ , } \int_{\cal{D}} c \ul{v} \cdot \ul{u} \mbox{ , } \int_{\cal{D}} c p r \mbox{ , } \int_{\cal{D}} c \ul{u} \cdot \ul{w} \\ \int_{\cal{D}} \ul{v} \cdot \ull{M} \cdot \ul{u} \mbox{ , } \int_{\cal{D}} \ul{u} \cdot \ull{M} \cdot \ul{w}
- Call signature:
de_dot
(opt_material, virtual, state)
(opt_material, parameter_1, parameter_2)
- Arguments 1:
material : c or \ull{M} (optional)
virtual : q or \ul{v}
state : p or \ul{u}
- Arguments 2:
material : c or \ull{M} (optional)
parameter_1 : p or \ul{u}
parameter_2 : r or \ul{w}
- arg_shapes = [{'opt_material': '1, 1', 'parameter_1': 1, 'parameter_2': 1, 'state': 1, 'virtual': (1, 'state')}, {'opt_material': None}, {'opt_material': '1, 1', 'parameter_1': 'D', 'parameter_2': 'D', 'state': 'D', 'virtual': ('D', 'state')}, {'opt_material': 'D, D'}, {'opt_material': None}]¶
- arg_types = (('opt_material', 'virtual', 'state'), ('opt_material', 'parameter_1', 'parameter_2'))¶
- integration = 'by_region'¶
- modes = ('weak', 'eval')¶
- name = 'de_dot'¶
- class sfepy.terms.terms_multilinear.EIntegrateOperatorTerm(*args, **kwargs)[source]¶
Volume and surface integral of a test function weighted by a scalar function c.
- Definition:
\int_\Omega q \mbox{ or } \int_\Omega c q
- Call signature:
de_integrate
(opt_material, virtual)
- Arguments:
material : c (optional)
virtual : q
- arg_shapes = [{'opt_material': '1, 1', 'virtual': (1, None)}, {'opt_material': None}]¶
- arg_types = ('opt_material', 'virtual')¶
- integration = 'by_region'¶
- name = 'de_integrate'¶
- class sfepy.terms.terms_multilinear.ELaplaceTerm(*args, **kwargs)[source]¶
Laplace term with c coefficient. Can be evaluated. Can use derivatives.
- Definition:
\int_{\Omega} c \nabla q \cdot \nabla p \mbox{ , } \int_{\Omega} c \nabla \bar{p} \cdot \nabla r
- Call signature:
de_laplace
(opt_material, virtual, state)
(opt_material, parameter_1, parameter_2)
- Arguments 1:
material : c
virtual : q
state : p
- Arguments 2:
material : c
parameter_1 : \bar{p}
parameter_2 : r
- arg_shapes = [{'opt_material': '1, 1', 'parameter_1': 1, 'parameter_2': 1, 'state': 1, 'virtual': (1, 'state')}, {'opt_material': None}]¶
- arg_types = (('opt_material', 'virtual', 'state'), ('opt_material', 'parameter_1', 'parameter_2'))¶
- modes = ('weak', 'eval')¶
- name = 'de_laplace'¶
- class sfepy.terms.terms_multilinear.ELinearElasticTerm(*args, **kwargs)[source]¶
General linear elasticity term, with D_{ijkl} given in the usual matrix form exploiting symmetry: in 3D it is 6\times6 with the indices ordered as [11, 22, 33, 12, 13, 23], in 2D it is 3\times3 with the indices ordered as [11, 22, 12].
- Definition:
\int_{\Omega} D_{ijkl}\ e_{ij}(\ul{v}) e_{kl}(\ul{u}) \mbox{ , } \int_{\Omega} D_{ijkl}\ e_{ij}(\ul{w}) e_{kl}(\ul{u})
- Call signature:
de_lin_elastic
(material, virtual, state)
(material, parameter_1, parameter_2)
- Arguments 1:
material : D_{ijkl}
virtual : \ul{v}
state : \ul{u}
- Arguments 2:
material : D_{ijkl}
parameter_1 : \ul{w}
parameter_2 : \ul{u}
- arg_shapes = {'material': 'S, S', 'parameter_1': 'D', 'parameter_2': 'D', 'state': 'D', 'virtual': ('D', 'state')}¶
- arg_types = (('material', 'virtual', 'state'), ('material', 'parameter_1', 'parameter_2'))¶
- modes = ('weak', 'eval')¶
- name = 'de_lin_elastic'¶
- class sfepy.terms.terms_multilinear.ENonPenetrationPenaltyTerm(*args, **kwargs)[source]¶
Non-penetration condition in the weak sense using a penalty.
- Definition:
\int_{\Gamma} c (\ul{n} \cdot \ul{v}) (\ul{n} \cdot \ul{u})
- Call signature:
de_non_penetration_p
(material, virtual, state)
- Arguments:
material : c
virtual : \ul{v}
state : \ul{u}
- arg_shapes = {'material': '1, 1', 'state': 'D', 'virtual': ('D', 'state')}¶
- arg_types = ('material', 'virtual', 'state')¶
- integration = 'surface'¶
- name = 'de_non_penetration_p'¶
- class sfepy.terms.terms_multilinear.EScalarDotMGradScalarTerm(*args, **kwargs)[source]¶
Volume dot product of a scalar gradient dotted with a material vector with a scalar.
- Definition:
\int_{\Omega} q \ul{y} \cdot \nabla p \mbox{ , } \int_{\Omega} p \ul{y} \cdot \nabla q
- Call signature:
de_s_dot_mgrad_s
(material, virtual, state)
(material, state, virtual)
- Arguments 1:
material : \ul{y}
virtual : q
state : p
- Arguments 2:
material : \ul{y}
state : p
virtual : q
- arg_shapes = [{'material': 'D, 1', 'state/grad_state': 1, 'state/grad_virtual': 1, 'virtual/grad_state': (1, None), 'virtual/grad_virtual': (1, None)}]¶
- arg_types = (('material', 'virtual', 'state'), ('material', 'state', 'virtual'))¶
- modes = ('grad_state', 'grad_virtual')¶
- name = 'de_s_dot_mgrad_s'¶
- class sfepy.terms.terms_multilinear.EStokesTerm(*args, **kwargs)[source]¶
Stokes problem coupling term. Corresponds to weak forms of gradient and divergence terms.
- Definition:
\int_{\Omega} p\ \nabla \cdot \ul{v} \mbox{ , } \int_{\Omega} q\ \nabla \cdot \ul{u} \mbox{ or } \int_{\Omega} c\ p\ \nabla \cdot \ul{v} \mbox{ , } \int_{\Omega} c\ q\ \nabla \cdot \ul{u} \\ \int_{\Omega} r\ \nabla \cdot \ul{w} \mbox{ , } \int_{\Omega} c r\ \nabla \cdot \ul{w}
- Call signature:
de_stokes
(opt_material, virtual, state)
(opt_material, state, virtual)
(opt_material, parameter_v, parameter_s)
- Arguments 1:
material : c (optional)
virtual : \ul{v}
state : p
- Arguments 2:
material : c (optional)
state : \ul{u}
virtual : q
- Arguments 3:
material : c (optional)
parameter_v : \ul{u}
parameter_s : p
- arg_shapes = [{'opt_material': '1, 1', 'parameter_s': 1, 'parameter_v': 'D', 'state/div': 'D', 'state/grad': 1, 'virtual/div': (1, None), 'virtual/grad': ('D', None)}, {'opt_material': None}]¶
- arg_types = (('opt_material', 'virtual', 'state'), ('opt_material', 'state', 'virtual'), ('opt_material', 'parameter_v', 'parameter_s'))¶
- modes = ('grad', 'div', 'eval')¶
- name = 'de_stokes'¶
- class sfepy.terms.terms_multilinear.ETermBase(*args, **kwargs)[source]¶
Reserved letters:
c .. cells q .. quadrature points d-h .. DOFs axes r-z .. auxiliary axes
Layout specification letters:
c .. cells q .. quadrature points v .. variable component - matrix form (v, d) -> vector v*d g .. gradient component d .. local DOF (basis, node) 0 .. all material axes
- can_backend = {'dask_single': None, 'dask_threads': None, 'jax': None, 'jax_vmap': None, 'numpy': <module 'numpy' from '/usr/lib/python3/dist-packages/numpy/__init__.py'>, 'numpy_loop': <module 'numpy' from '/usr/lib/python3/dist-packages/numpy/__init__.py'>, 'numpy_qloop': <module 'numpy' from '/usr/lib/python3/dist-packages/numpy/__init__.py'>, 'opt_einsum': None, 'opt_einsum_dask_single': None, 'opt_einsum_dask_threads': None, 'opt_einsum_loop': None, 'opt_einsum_qloop': None}¶
- layout_letters = 'cqgvd0'¶
- verbosity = 0¶