16#ifndef dealii_integrators_elasticity_h
17#define dealii_integrators_elasticity_h
53 const double factor = 1.)
63 const double dx = factor * fe.
JxW(k);
64 for (
unsigned int i = 0; i < n_dofs; ++i)
65 for (
unsigned int j = 0; j < n_dofs; ++j)
66 for (
unsigned int d1 = 0; d1 < dim; ++d1)
67 for (
unsigned int d2 = 0; d2 < dim; ++d2)
82 template <
int dim,
typename number>
97 for (
unsigned int k = 0; k < nq; ++k)
99 const double dx = factor * fe.
JxW(k);
100 for (
unsigned int i = 0; i < n_dofs; ++i)
101 for (
unsigned int d1 = 0; d1 < dim; ++d1)
102 for (
unsigned int d2 = 0; d2 < dim; ++d2)
104 result(i) +=
dx * .25 *
105 (input[d1][k][d2] + input[d2][k][d1]) *
136 const double dx = factor * fe.
JxW(k);
138 for (
unsigned int i = 0; i < n_dofs; ++i)
139 for (
unsigned int j = 0; j < n_dofs; ++j)
140 for (
unsigned int d1 = 0; d1 < dim; ++d1)
144 M(i, j) +=
dx * 2. * penalty * u * v;
145 for (
unsigned int d2 = 0; d2 < dim; ++d2)
191 const double dx = factor * fe.
JxW(k);
193 for (
unsigned int i = 0; i < n_dofs; ++i)
194 for (
unsigned int j = 0; j < n_dofs; ++j)
201 for (
unsigned int d = 0;
d < dim; ++
d)
208 for (
unsigned int d1 = 0; d1 < dim; ++d1)
214 M(i, j) +=
dx * 2. * penalty * u * v;
216 M(i, j) +=
dx * (ngradun * v + ngradvn * u);
217 for (
unsigned int d2 = 0; d2 < dim; ++d2)
255 template <
int dim,
typename number>
259 const ArrayView<
const std::vector<double>> & input,
261 const ArrayView<
const std::vector<double>> & data,
272 const double dx = factor * fe.
JxW(k);
274 for (
unsigned int i = 0; i < n_dofs; ++i)
275 for (
unsigned int d1 = 0; d1 < dim; ++d1)
277 const double u = input[d1][k];
279 const double g = data[d1][k];
280 result(i) +=
dx * 2. * penalty * (u - g) * v;
282 for (
unsigned int d2 = 0; d2 < dim; ++d2)
285 result(i) -= .5 *
dx * v * Dinput[d1][k][d2] * n[d2];
287 result(i) -= .5 *
dx * v * Dinput[d2][k][d1] * n[d2];
289 result(i) -= .5 *
dx * (u - g) *
292 result(i) -= .5 *
dx * (u - g) *
307 template <
int dim,
typename number>
312 const ArrayView<
const std::vector<double>> & input,
314 const ArrayView<
const std::vector<double>> & data,
325 const double dx = factor * fe.
JxW(k);
327 for (
unsigned int i = 0; i < n_dofs; ++i)
335 for (
unsigned int d = 0;
d < dim; ++
d)
337 udotn += n[
d] * input[
d][k];
338 gdotn += n[
d] * data[
d][k];
340 ngradun += n * Dinput[
d][k] * n[
d];
343 for (
unsigned int d1 = 0; d1 < dim; ++d1)
345 const double u = input[d1][k] - udotn * n[d1];
348 const double g = data[d1][k] - gdotn * n[d1];
349 result(i) +=
dx * 2. * penalty * (u - g) * v;
351 result(i) +=
dx * (ngradun * v + ngradvn * (u - g));
352 for (
unsigned int d2 = 0; d2 < dim; ++d2)
355 result(i) -= .5 *
dx * Dinput[d1][k][d2] * n[d2] * v;
357 result(i) -= .5 *
dx * Dinput[d2][k][d1] * n[d2] * v;
359 result(i) -= .5 *
dx * (u - g) *
363 result(i) -= .5 *
dx * (u - g) *
385 template <
int dim,
typename number>
390 const ArrayView<
const std::vector<double>> & input,
401 const double dx = factor * fe.
JxW(k);
403 for (
unsigned int i = 0; i < n_dofs; ++i)
404 for (
unsigned int d1 = 0; d1 < dim; ++d1)
406 const double u = input[d1][k];
408 result(i) +=
dx * 2. * penalty * u * v;
410 for (
unsigned int d2 = 0; d2 < dim; ++d2)
413 result(i) -= .5 *
dx * v * Dinput[d1][k][d2] * n[d2];
415 result(i) -= .5 *
dx * v * Dinput[d2][k][d1] * n[d2];
417 result(i) -= .5 *
dx * u *
420 result(i) -= .5 *
dx * u *
439 const double int_factor = 1.,
440 const double ext_factor = -1.)
455 const double nu1 = int_factor;
456 const double nu2 = (ext_factor < 0) ? int_factor : ext_factor;
457 const double penalty = .5 * pen * (nu1 + nu2);
461 const double dx = fe1.
JxW(k);
463 for (
unsigned int i = 0; i < n_dofs; ++i)
464 for (
unsigned int j = 0; j < n_dofs; ++j)
465 for (
unsigned int d1 = 0; d1 < dim; ++d1)
472 M11(i, j) +=
dx * penalty * u1 *
v1;
473 M12(i, j) -=
dx * penalty * u2 *
v1;
474 M21(i, j) -=
dx * penalty * u1 * v2;
475 M22(i, j) +=
dx * penalty * u2 * v2;
477 for (
unsigned int d2 = 0; d2 < dim; ++d2)
480 M11(i, j) -= .25 *
dx * nu1 *
483 M12(i, j) -= .25 *
dx * nu2 *
486 M21(i, j) += .25 *
dx * nu1 *
489 M22(i, j) += .25 *
dx * nu2 *
493 M11(i, j) -= .25 *
dx * nu1 *
496 M12(i, j) -= .25 *
dx * nu2 *
499 M21(i, j) += .25 *
dx * nu1 *
502 M22(i, j) += .25 *
dx * nu2 *
506 M11(i, j) -= .25 *
dx * nu1 *
509 M12(i, j) += .25 *
dx * nu1 *
512 M21(i, j) -= .25 *
dx * nu2 *
515 M22(i, j) += .25 *
dx * nu2 *
519 M11(i, j) -= .25 *
dx * nu1 *
522 M12(i, j) += .25 *
dx * nu1 *
525 M21(i, j) -= .25 *
dx * nu2 *
528 M22(i, j) += .25 *
dx * nu2 *
538 template <
int dim,
typename number>
544 const ArrayView<
const std::vector<double>> & input1,
546 const ArrayView<
const std::vector<double>> & input2,
549 double int_factor = 1.,
550 double ext_factor = -1.)
561 const double nu1 = int_factor;
562 const double nu2 = (ext_factor < 0) ? int_factor : ext_factor;
563 const double penalty = .5 * pen * (nu1 + nu2);
568 const double dx = fe1.
JxW(k);
571 for (
unsigned int i = 0; i < n1; ++i)
572 for (
unsigned int d1 = 0; d1 < dim; ++d1)
576 const double u1 = input1[d1][k];
577 const double u2 = input2[d1][k];
579 result1(i) +=
dx * penalty * u1 *
v1;
580 result1(i) -=
dx * penalty * u2 *
v1;
581 result2(i) -=
dx * penalty * u1 * v2;
582 result2(i) +=
dx * penalty * u2 * v2;
584 for (
unsigned int d2 = 0; d2 < dim; ++d2)
589 (nu1 * Dinput1[d1][k][d2] + nu2 * Dinput2[d1][k][d2]) *
593 (nu1 * Dinput1[d1][k][d2] + nu2 * Dinput2[d1][k][d2]) *
598 (nu1 * Dinput1[d2][k][d1] + nu2 * Dinput2[d2][k][d1]) *
602 (nu1 * Dinput1[d2][k][d1] + nu2 * Dinput2[d2][k][d1]) *
605 result1(i) -= .25 *
dx * nu1 *
608 result2(i) -= .25 *
dx * nu2 *
612 result1(i) -= .25 *
dx * nu1 *
615 result2(i) -= .25 *
dx * nu2 *
const unsigned int dofs_per_cell
const unsigned int n_quadrature_points
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
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
double JxW(const unsigned int quadrature_point) const
const FiniteElement< dim, spacedim > & get_fe() const
unsigned int n_components() const
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define Assert(cond, exc)
#define AssertVectorVectorDimension(VEC, DIM1, DIM2)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
void nitsche_tangential_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, double penalty, double factor=1.)
void ip_matrix(FullMatrix< double > &M11, FullMatrix< double > &M12, FullMatrix< double > &M21, FullMatrix< double > &M22, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, const double pen, const double int_factor=1., const double ext_factor=-1.)
void nitsche_residual(Vector< number > &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 nitsche_residual_homogeneous(Vector< number > &result, const FEValuesBase< dim > &fe, const ArrayView< const std::vector< double > > &input, const ArrayView< const std::vector< Tensor< 1, dim > > > &Dinput, double penalty, double factor=1.)
void ip_residual(Vector< number > &result1, Vector< number > &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 nitsche_tangential_residual(Vector< number > &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 nitsche_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, double penalty, double factor=1.)
void cell_residual(Vector< number > &result, const FEValuesBase< dim > &fe, const ArrayView< const std::vector< Tensor< 1, dim > > > &input, double factor=1.)
void cell_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const double factor=1.)
Library of integrals over cells and faces.
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)