|
Reference documentation for deal.II version 9.2.0
|
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\)
\(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\)
\(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\)
\(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Go to the documentation of this file.
16 #ifndef dealii_integrators_grad_div_h
17 #define dealii_integrators_grad_div_h
69 const double dx = factor * fe.
JxW(k);
70 for (
unsigned int i = 0; i < n_dofs; ++i)
71 for (
unsigned int j = 0; j < n_dofs; ++j)
78 M(i, j) +=
dx * divu * divv;
92 template <
int dim,
typename number>
97 const double factor = 1.)
106 const double dx = factor * fetest.
JxW(k);
107 for (
unsigned int i = 0; i < n_dofs; ++i)
112 for (
unsigned int d = 0;
d < dim; ++
d)
113 du += input[
d][k][
d];
115 result(i) +=
dx * du * divv;
143 const double dx = factor * fe.
JxW(k);
145 for (
unsigned int i = 0; i < n_dofs; ++i)
146 for (
unsigned int j = 0; j < n_dofs; ++j)
152 double un = 0., vn = 0.;
153 for (
unsigned int d = 0;
d < dim; ++
d)
159 M(i, j) +=
dx * 2. * penalty * un * vn;
160 M(i, j) -=
dx * (divu * vn + divv * un);
187 const ArrayView<
const std::vector<double>> & input,
189 const ArrayView<
const std::vector<double>> & data,
201 const double dx = factor * fe.
JxW(k);
206 for (
unsigned int d = 0;
d < dim; ++
d)
208 umgn += (input[
d][k] - data[
d][k]) * n[
d];
209 divu += Dinput[
d][k][
d];
212 for (
unsigned int i = 0; i < n_dofs; ++i)
217 for (
unsigned int d = 0;
d < dim; ++
d)
221 dx * (2. * penalty * umgn * vn - divv * umgn - divu * vn);
244 double factor2 = -1.)
256 const double fi = factor1;
257 const double fe = (factor2 < 0) ? factor1 : factor2;
258 const double f = .5 * (fi + fe);
262 const double dx = fe1.
JxW(k);
264 for (
unsigned int i = 0; i < n_dofs; ++i)
265 for (
unsigned int j = 0; j < n_dofs; ++j)
280 for (
unsigned int d = 0;
d < dim; ++
d)
288 dx * (-.5 * fi * divvi * uni - .5 * fi * divui * vni +
289 f * penalty * uni * vni);
291 dx * (.5 * fi * divvi * une - .5 * fe * divue * vni -
292 f * penalty * vni * une);
294 dx * (-.5 * fe * divve * uni + .5 * fi * divui * vne -
295 f * penalty * uni * vne);
297 dx * (.5 * fe * divve * une + .5 * fe * divue * vne +
298 f * penalty * une * vne);
323 const ArrayView<
const std::vector<double>> & input1,
325 const ArrayView<
const std::vector<double>> & input2,
328 double int_factor = 1.,
329 double ext_factor = -1.)
339 const double fi = int_factor;
340 const double fe = (ext_factor < 0) ? int_factor : ext_factor;
341 const double penalty = .5 * pen * (fi + fe);
346 const double dx = fe1.
JxW(k);
352 for (
unsigned int d = 0;
d < dim; ++
d)
354 uni += input1[
d][k] * n[
d];
355 une += input2[
d][k] * n[
d];
356 divui += Dinput1[
d][k][
d];
357 divue += Dinput2[
d][k][
d];
360 for (
unsigned int i = 0; i < n1; ++i)
368 for (
unsigned int d = 0;
d < dim; ++
d)
374 result1(i) +=
dx * (-.5 * fi * divvi * uni -
375 .5 * fi * divui * vni + penalty * uni * vni);
376 result1(i) +=
dx * (.5 * fi * divvi * une -
377 .5 * fe * divue * vni - penalty * vni * une);
378 result2(i) +=
dx * (-.5 * fe * divve * uni +
379 .5 * fi * divui * vne - penalty * uni * vne);
380 result2(i) +=
dx * (.5 * fe * divve * une +
381 .5 * fe * divue * vne + penalty * une * vne);
void nitsche_residual(Vector< double > &result, const FEValuesBase< dim > &fe, const ArrayView< const std::vector< double >> &input, const ArrayView< const std::vector< Tensor< 1, dim >>> &Dinput, const ArrayView< const std::vector< double >> &data, double penalty, double factor=1.)
void ip_residual(Vector< double > &result1, Vector< double > &result2, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, const ArrayView< const std::vector< double >> &input1, const ArrayView< const std::vector< Tensor< 1, dim >>> &Dinput1, const ArrayView< const std::vector< double >> &input2, const ArrayView< const std::vector< Tensor< 1, dim >>> &Dinput2, double pen, double int_factor=1., double ext_factor=-1.)
void cell_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, double factor=1.)
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const Tensor< 1, spacedim > & normal_vector(const unsigned int i) const
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
void ip_matrix(FullMatrix< double > &M11, FullMatrix< double > &M12, FullMatrix< double > &M21, FullMatrix< double > &M22, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, double penalty, double factor1=1., double factor2=-1.)
const unsigned int dofs_per_cell
#define AssertVectorVectorDimension(VEC, DIM1, DIM2)
#define DEAL_II_NAMESPACE_OPEN
#define AssertDimension(dim1, dim2)
void cell_residual(Vector< number > &result, const FEValuesBase< dim > &fetest, const ArrayView< const std::vector< Tensor< 1, dim >>> &input, const double factor=1.)
double JxW(const unsigned int quadrature_point) const
Library of integrals over cells and faces.
const FiniteElement< dim, spacedim > & get_fe() const
void nitsche_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, double penalty, double factor=1.)
const unsigned int n_quadrature_points
#define DEAL_II_NAMESPACE_CLOSE