106 template <
class InputVector,
int spacedim>
109 const InputVector &solution,
110 const unsigned int component);
136 template <
class InputVector,
int spacedim>
140 const InputVector &solution,
141 const unsigned int component)
145 std::vector<typename InputVector::value_type>
values(1);
151 std::vector<Vector<typename InputVector::value_type>>
values(
156 return values[0](component);
167 for (
unsigned int i = 0; i < dim; ++i)
217 template <
class InputVector,
int spacedim>
220 const InputVector &solution,
221 const unsigned int component);
251 template <
class InputVector,
int spacedim>
255 const InputVector &solution,
256 const unsigned int component)
260 std::vector<Tensor<1, dim, typename InputVector::value_type>>
values(
268 std::vector<Tensor<1, dim, typename InputVector::value_type>>>
301 const double radicand =
302 ::sqr(
d[0][0] -
d[1][1]) + 4 * ::sqr(
d[0][1]);
304 0.5 * (
d[0][0] +
d[1][1] + std::sqrt(radicand)),
305 0.5 * (
d[0][0] +
d[1][1] - std::sqrt(radicand))};
425 const double am =
trace(
d) / 3.;
429 for (
unsigned int i = 0; i < 3; ++i)
432 const double ss01 = s[0][1] * s[0][1], ss12 = s[1][2] * s[1][2],
433 ss02 = s[0][2] * s[0][2];
435 const double J2 = (s[0][0] * s[0][0] + s[1][1] * s[1][1] +
436 s[2][2] * s[2][2] + 2 * (ss01 + ss02 + ss12)) /
439 (Utilities::fixed_power<3>(s[0][0]) +
440 Utilities::fixed_power<3>(s[1][1]) +
441 Utilities::fixed_power<3>(s[2][2]) + 3. * s[0][0] * (ss01 + ss02) +
442 3. * s[1][1] * (ss01 + ss12) + 3. * s[2][2] * (ss02 + ss12) +
443 6. * s[0][1] * s[0][2] * s[1][2]) /
446 const double R = std::sqrt(4. * J2 / 3.);
448 double EE[3] = {0, 0, 0};
456 EE[0] = EE[1] = EE[2] = am;
461 const double R3 = R * R * R;
462 const double XX = 4. * J3 / R3;
470 const double a = (XX > 0 ? -1. : 1.) * R / 2;
471 EE[0] = EE[1] = am + a;
477 EE[0] = am + R * std::cos(theta);
478 EE[1] = am + R * std::cos(theta + 2. / 3. *
numbers::PI);
479 EE[2] = am + R * std::cos(theta + 4. / 3. *
numbers::PI);
517 for (
unsigned int i = 0; i < dim; ++i)
518 for (
unsigned int j = i + 1; j < dim; ++j)
520 const double s = (
d[i][j] +
d[j][i]) / 2;
521 d[i][j] =
d[j][i] = s;
556 template <
class InputVector,
int spacedim>
559 const InputVector &solution,
560 const unsigned int component);
590 template <
class InputVector,
int spacedim>
594 const InputVector &solution,
595 const unsigned int component)
599 std::vector<Tensor<2, dim, typename InputVector::value_type>>
values(
607 std::vector<Tensor<2, dim, typename InputVector::value_type>>>
647 for (
unsigned int i = 0; i < dim; ++i)
648 for (
unsigned int j = i + 1; j < dim; ++j)
649 for (
unsigned int k = j + 1; k < dim; ++k)
651 const double s = (
d[i][j][k] +
d[i][k][j] +
d[j][i][k] +
652 d[j][k][i] +
d[k][i][j] +
d[k][j][i]) /
654 d[i][j][k] =
d[i][k][j] =
d[j][i][k] =
d[j][k][i] =
d[k][i][j] =
659 for (
unsigned int i = 0; i < dim; ++i)
660 for (
unsigned int j = i + 1; j < dim; ++j)
664 const double s = (
d[i][i][j] +
d[i][j][i] +
d[j][i][i]) / 3;
665 d[i][i][j] =
d[i][j][i] =
d[j][i][i] = s;
669 const double t = (
d[i][j][j] +
d[j][i][j] +
d[j][j][i]) / 3;
670 d[i][j][j] =
d[j][i][j] =
d[j][j][i] = t;
675 template <
int order,
int dim>
740 template <
class DerivativeDescription,
748 const InputVector &solution,
749 const unsigned int component,
751 typename DerivativeDescription::Derivative &derivative)
781 std::vector<typename DoFHandler<dim, spacedim>::active_cell_iterator>
791 typename DerivativeDescription::Derivative projected_derivative;
794 x_fe_midpoint_value.
reinit(cell);
800 const typename DerivativeDescription::ProjectedDerivative
801 this_midpoint_value =
802 DerivativeDescription::get_projected_derivative(fe_midpoint_value,
824 GridTools::get_active_neighbors<DoFHandler<dim, spacedim>>(
825 cell, active_neighbors);
830 auto neighbor_ptr = active_neighbors.begin();
831 for (; neighbor_ptr != active_neighbors.end(); ++neighbor_ptr)
833 const auto neighbor = *neighbor_ptr;
836 x_fe_midpoint_value.
reinit(neighbor);
842 const typename DerivativeDescription::ProjectedDerivative
843 neighbor_midpoint_value =
844 DerivativeDescription::get_projected_derivative(
845 neighbor_fe_midpoint_value, solution, component);
859 const double distance = y.
norm();
873 for (
unsigned int i = 0; i < dim; ++i)
874 for (
unsigned int j = 0; j < dim; ++j)
875 Y[i][j] += y[i] * y[j];
880 typename DerivativeDescription::ProjectedDerivative
881 projected_finite_difference =
882 (neighbor_midpoint_value - this_midpoint_value);
883 projected_finite_difference /= distance;
885 projected_derivative +=
outer_product(y, projected_finite_difference);
902 derivative = Y_inverse * projected_derivative;
915 template <
class DerivativeDescription,
926 const InputVector &solution,
927 const unsigned int component)
930 if (std::get<0>(*cell)->is_locally_owned() ==
false)
931 *std::get<1>(*cell) = 0;
934 typename DerivativeDescription::Derivative derivative;
937 approximate_cell<DerivativeDescription, dim, InputVector, spacedim>(
947 *std::get<1>(*cell) =
963 template <
class DerivativeDescription,
970 const InputVector &solution,
971 const unsigned int component,
982 std::tuple<typename DoFHandler<dim, spacedim>::active_cell_iterator,
994 [&mapping, &dof_handler, &solution, component](
998 approximate<DerivativeDescription, dim, InputVector, spacedim>(
999 cell, mapping, dof_handler, solution, component);
1001 std::function<void(const internal::Assembler::CopyData &)>(),
1015 template <
int dim,
class InputVector,
int spacedim>
1019 const InputVector &solution,
1021 const unsigned int component)
1023 internal::approximate_derivative<internal::Gradient<dim>, dim>(
1028 template <
int dim,
class InputVector,
int spacedim>
1031 const InputVector &solution,
1033 const unsigned int component)
1039 internal::approximate_derivative<internal::Gradient<dim>, dim>(
1040 reference_cell.template get_default_linear_mapping<dim, spacedim>(),
1048 template <
int dim,
class InputVector,
int spacedim>
1052 const InputVector &solution,
1054 const unsigned int component)
1056 internal::approximate_derivative<internal::SecondDerivative<dim>, dim>(
1061 template <
int dim,
class InputVector,
int spacedim>
1064 const InputVector &solution,
1066 const unsigned int component)
1072 internal::approximate_derivative<internal::SecondDerivative<dim>, dim>(
1073 reference_cell.template get_default_linear_mapping<dim, spacedim>(),
1081 template <
int dim,
int spacedim,
class InputVector,
int order>
1086 const InputVector &solution,
1094 const unsigned int component)
1098 mapping, dof, solution, component, cell, derivative);
1103 template <
int dim,
int spacedim,
class InputVector,
int order>
1107 const InputVector &solution,
1115 const unsigned int component)
1119 cell->reference_cell()
1120 .template get_default_linear_mapping<dim, spacedim>(),
1130 template <
int dim,
int order>
1142 #include "derivative_approximation.inst"
static ProjectedDerivative get_projected_derivative(const FEValues< dim, spacedim > &fe_values, const InputVector &solution, const unsigned int component)
static const UpdateFlags update_flags
static void symmetrize(Derivative &derivative_tensor)
static double derivative_norm(const Derivative &d)
static const UpdateFlags update_flags
static ProjectedDerivative get_projected_derivative(const FEValues< dim, spacedim > &fe_values, const InputVector &solution, const unsigned int component)
static double derivative_norm(const Derivative &d)
static void symmetrize(Derivative &derivative_tensor)
static ProjectedDerivative get_projected_derivative(const FEValues< dim, spacedim > &fe_values, const InputVector &solution, const unsigned int component)
static double derivative_norm(const Derivative &d)
static const UpdateFlags update_flags
static void symmetrize(Derivative &derivative_tensor)
cell_iterator end() const
const Triangulation< dim, spacedim > & get_triangulation() const
active_cell_iterator begin_active(const unsigned int level=0) const
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
const hp::FECollection< dim, spacedim > & get_fe_collection() const
void get_function_hessians(const ReadVector< Number > &fe_function, std::vector< Tensor< 2, spacedim, Number >> &hessians) const
const Point< spacedim > & quadrature_point(const unsigned int q_point) const
void get_function_values(const ReadVector< Number > &fe_function, std::vector< Number > &values) const
void get_function_gradients(const ReadVector< Number > &fe_function, std::vector< Tensor< 1, spacedim, Number >> &gradients) const
const FiniteElement< dim, spacedim > & get_fe() const
unsigned int n_components() const
numbers::NumberTraits< Number >::real_type norm() const
unsigned int n_active_cells() const
bool is_mixed_mesh() const
const std::vector< ReferenceCell > & get_reference_cells() const
const FEValuesType & get_present_fe_values() const
void reinit(const TriaIterator< DoFCellAccessor< dim, spacedim, lda >> &cell, const unsigned int q_index=numbers::invalid_unsigned_int, const unsigned int mapping_index=numbers::invalid_unsigned_int, const unsigned int fe_index=numbers::invalid_unsigned_int)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcInsufficientDirections()
static ::ExceptionBase & ExcVectorLengthVsNActiveCells(int arg1, int arg2)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
#define AssertIndexRange(index, range)
#define AssertThrow(cond, exc)
@ update_hessians
Second derivatives of shape functions.
@ update_values
Shape function values.
@ update_gradients
Shape function gradients.
@ update_quadrature_points
Transformed quadrature points.
void approximate_derivative(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof_handler, const InputVector &solution, const unsigned int component, Vector< float > &derivative_norm)
void approximate_cell(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof_handler, const InputVector &solution, const unsigned int component, const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell, typename DerivativeDescription::Derivative &derivative)
void approximate(const SynchronousIterators< std::tuple< typename DoFHandler< dim, spacedim >::active_cell_iterator, Vector< float >::iterator >> &cell, const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof_handler, const InputVector &solution, const unsigned int component)
void approximate_gradient(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const InputVector &solution, Vector< float > &derivative_norm, const unsigned int component=0)
double derivative_norm(const Tensor< order, dim > &derivative)
void approximate_derivative_tensor(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const InputVector &solution, const typename DoFHandler< dim, spacedim >::active_cell_iterator &cell, Tensor< order, dim > &derivative, const unsigned int component=0)
void approximate_second_derivative(const Mapping< dim, spacedim > &mapping, const DoFHandler< dim, spacedim > &dof, const InputVector &solution, Vector< float > &derivative_norm, const unsigned int component=0)
Expression fabs(const Expression &x)
Expression acos(const Expression &x)
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
@ eigenvalues
Eigenvalue vector is filled.
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
VectorType::value_type * begin(VectorType &V)
VectorType::value_type * end(VectorType &V)
void run(const std::vector< std::vector< Iterator >> &colored_iterators, Worker worker, Copier copier, const ScratchData &sample_scratch_data, const CopyData &sample_copy_data, const unsigned int queue_length=2 *MultithreadInfo::n_threads(), const unsigned int chunk_size=8)
static constexpr double PI
constexpr DEAL_II_HOST SymmetricTensor< 4, dim, Number > outer_product(const SymmetricTensor< 2, dim, Number > &t1, const SymmetricTensor< 2, dim, Number > &t2)
constexpr DEAL_II_HOST Number trace(const SymmetricTensor< 2, dim2, Number > &)
constexpr DEAL_II_HOST SymmetricTensor< 2, dim, Number > invert(const SymmetricTensor< 2, dim, Number > &)
constexpr DEAL_II_HOST SymmetricTensor< 2, dim, Number > symmetrize(const Tensor< 2, dim, Number > &t)
constexpr DEAL_II_HOST Number determinant(const SymmetricTensor< 2, dim, Number > &)