piola_fem_grad_post.icc
Go to the documentation of this file.
#ifndef _RHEOLEF_PIOLA_FEM_GRAD_POST_ICC
#define _RHEOLEF_PIOLA_FEM_GRAD_POST_ICC
#include "piola_fem.h"
namespace rheolef {
// ----------------------------------------------------------------------------
// grad-post treatment: tensor symmetry & surface projector
// ----------------------------------------------------------------------------
// vector-valued grad:
template<class T>
static
void
grad_post (
const piola<T>& p,
const details::differentiate_option& gopt,
const T& u,
point_basic<T>& grad_u)
{
if (gopt.surfacic && p.map_d + 1 == p.d) {
// apply also the tangential projection P=(I-n*n)
// grad_s(u) = P*grad(u)
grad_u = p.P*grad_u; // TODO: DVT_OPTIM_2D
}
}
// tensor-valued grad:
template<class T>
static
void
grad_post (
const piola<T>& p,
const details::differentiate_option& gopt,
const point_basic<T>& u,
tensor_basic<T>& grad_u)
{
if (gopt.symmetrized) {
// compute D(u) from grad(u):
grad_u = (grad_u + trans(grad_u))/2; // TODO: DVT_OPTIM_2D
}
if (gopt.surfacic) {
// will apply also the tangential projection P=(I-n*n)
// surfacic grad, where P=I-nxn is the map projector
// grad_s(u) = grad(u)*P
// or
// Ds(u) = P*D(u)*P
if (! gopt.symmetrized) {
grad_u = grad_u*p.P; // TODO: DVT_OPTIM_2D
} else {
grad_u = p.P*grad_u*p.P; // TODO: DVT_OPTIM_2D
}
}
if (! p.ignore_sys_coord && p.sys_coord != space_constant::cartesian) {
// set the specific grad_u(theta,theta) = vr/r term in the axi case
size_t i_comp_r = (p.sys_coord == space_constant::axisymmetric_rz) ? 0 : 1;
size_t i_comp_theta = 2;
const T& ur = u [i_comp_r];
const T& r = p.F [i_comp_r];
// pb: r==0 when using gauss_lobatto, e.g. for the characteristic method
// e.g. in VF: 2D(u):D(v) dx -> (ur/r)*(vr/r)*r dr*dz = (ur*vr)/r dr*dz
// or for direct interpolation of grad(uh), D(uh), div(uh), etc in the axisym. case
check_macro (1+abs(r) != 1, "grad(): singular axisymmetric (1/r) weight (HINT: avoid interpolate() or change quadrature formulae)");
grad_u (i_comp_theta, i_comp_theta) = ur/r;
}
}
// tensor3-valued grad:
template<class T>
static
void
grad_post (
const piola<T>& p,
const details::differentiate_option& gopt,
const tensor_basic<T>& u,
tensor3_basic<T>& grad_u)
{
// only used by u.grad(tau) = grad(tau)*u
// e.g. in doc/pexamples/transport_tensor_dg.cc : grad_h(sigma)*u
// From BirArmHas-vol1-1987 page 589 : formula for u.grad(tau)
// here, in zr or rz geometries, the u_theta component is zero
// and then there is no additional component to add here
// TODO: DVT_TENSOR3_GRAD_AXI otherwise, for a general grad(tau) tensor3 usage, there are additional components
}
}// namespace rheolef
#endif // _RHEOLEF_PIOLA_FEM_GRAD_POST_ICC
rheolef::tensor3_basic< T >
tensor_basic< T > tensor3_basic< T >
Definition: piola_fem.h:137
rheolef::point_basic< T >
point_basic< T >
Definition: piola_fem.h:135
T
Expr1::float_type T
Definition: field_expr.h:261
check_macro
check_macro(expr1.have_homogeneous_space(Xh1), "dual(expr1,expr2); expr1 should have homogeneous space. HINT: use dual(interpolate(Xh, expr1),expr2)")
rheolef::space_constant::cartesian
@ cartesian
Definition: space_constant.h:122
rheolef::space_constant::axisymmetric_rz
@ axisymmetric_rz
Definition: space_constant.h:123
rheolef
This file is part of Rheolef.
Definition: compiler_eigen.h:37
rheolef::trans
csr< T, sequential > trans(const csr< T, sequential > &a)
trans(a): see the form page for the full documentation
Definition: csr.h:455
piola_fem.h
grad_u
Definition: combustion_exact.icc:34
p
Definition: sphere.icc:25
p::d
size_t d
Definition: sphere.icc:32
u
Definition: leveque.h:25
grad_u
tensor grad_u
Definition: transport_tensor_exact.icc:26