9 #ifndef __allen_cahn_polynomial_h 10 #define __allen_cahn_polynomial_h 12 #include <amandus/integrator.h> 13 #include <deal.II/base/polynomial.h> 14 #include <deal.II/base/tensor.h> 15 #include <deal.II/integrators/divergence.h> 16 #include <deal.II/integrators/l2.h> 17 #include <deal.II/integrators/laplace.h> 36 PolynomialResidual(
double diffusion,
const Polynomials::Polynomial<double> solution_1d);
38 virtual void cell(DoFInfo<dim>& dinfo, IntegrationInfo<dim>& info)
const;
39 virtual void boundary(DoFInfo<dim>& dinfo, IntegrationInfo<dim>& info)
const;
40 virtual void face(DoFInfo<dim>& dinfo1, DoFInfo<dim>& dinfo2, IntegrationInfo<dim>& info1,
41 IntegrationInfo<dim>& info2)
const;
54 virtual void cell(DoFInfo<dim>& dinfo, IntegrationInfo<dim>& info)
const;
55 virtual void boundary(DoFInfo<dim>& dinfo, IntegrationInfo<dim>& info)
const;
56 virtual void face(DoFInfo<dim>& dinfo1, DoFInfo<dim>& dinfo2, IntegrationInfo<dim>& info1,
57 IntegrationInfo<dim>& info2)
const;
67 const Polynomials::Polynomial<double> solution_1d)
69 , solution_1d(solution_1d)
71 this->use_boundary =
false;
72 this->use_face =
true;
73 this->input_vector_names.push_back(
"Newton iterate");
80 Assert(info.values.size() >= 1, ExcDimensionMismatch(info.values.size(), 1));
81 Assert(info.gradients.size() >= 1, ExcDimensionMismatch(info.values.size(), 1));
83 std::vector<double> rhs(info.fe_values(0).n_quadrature_points, 0.);
85 std::vector<double> px(3);
86 std::vector<double> py(3);
87 for (
unsigned int k = 0; k < info.fe_values(0).n_quadrature_points; ++k)
89 const double x = info.fe_values(0).quadrature_point(k)(0);
90 const double y = info.fe_values(0).quadrature_point(k)(1);
95 rhs[k] =
D * px[2] * py[0] +
D * px[0] * py[2];
97 rhs[k] -= px[0] * py[0] * (px[0] * py[0] * px[0] * py[0] - 1.);
99 const double u = info.values[0][0][k];
100 rhs[k] += u * (u * u - 1.);
103 L2::L2(dinfo.vector(0).block(0), info.fe_values(0), rhs);
111 std::vector<double> null(info.fe_values(0).n_quadrature_points, 0.);
113 const unsigned int deg = info.fe_values(0).get_fe().tensor_degree();
114 Laplace::nitsche_residual(dinfo.vector(0).block(0),
117 info.gradients[0][0],
119 Laplace::compute_penalty(dinfo, dinfo, deg, deg),
126 IntegrationInfo<dim>& info1, IntegrationInfo<dim>& info2)
const 128 const unsigned int deg = info1.fe_values(0).get_fe().tensor_degree();
129 Laplace::ip_residual(dinfo1.vector(0).block(0),
130 dinfo2.vector(0).block(0),
134 info1.gradients[0][0],
136 info2.gradients[0][0],
137 Laplace::compute_penalty(dinfo1, dinfo2, deg, deg),
145 : solution_1d(solution_1d)
147 this->use_boundary =
false;
148 this->use_face =
false;
157 std::vector<double> px(2);
158 std::vector<double> py(2);
159 for (
unsigned int k = 0; k < info.fe_values(0).n_quadrature_points; ++k)
161 const double x = info.fe_values(0).quadrature_point(k)(0);
162 const double y = info.fe_values(0).quadrature_point(k)(1);
165 const double dx = info.fe_values(0).JxW(k);
167 Tensor<1, dim> Du = info.gradients[0][0][k];
168 Du[0] -= px[1] * py[0];
169 Du[1] -= px[0] * py[1];
170 double u = info.values[0][0][k];
174 dinfo.value(i++) += u * u * dx;
176 dinfo.value(i++) += (Du * Du) * dx;
189 IntegrationInfo<dim>&)
const virtual void face(DoFInfo< dim > &dinfo1, DoFInfo< dim > &dinfo2, IntegrationInfo< dim > &info1, IntegrationInfo< dim > &info2) const
Definition: allen_cahn/polynomial.h:188
virtual void boundary(DoFInfo< dim > &dinfo, IntegrationInfo< dim > &info) const
Definition: allen_cahn/polynomial.h:182
Polynomials::Polynomial< double > solution_1d
Definition: allen_cahn/polynomial.h:45
virtual void cell(DoFInfo< dim > &dinfo, IntegrationInfo< dim > &info) const
Definition: allen_cahn/polynomial.h:78
virtual void face(DoFInfo< dim > &dinfo1, DoFInfo< dim > &dinfo2, IntegrationInfo< dim > &info1, IntegrationInfo< dim > &info2) const
Definition: allen_cahn/polynomial.h:125
PolynomialResidual(double diffusion, const Polynomials::Polynomial< double > solution_1d)
Definition: allen_cahn/polynomial.h:66
Definition: allen_cahn/polynomial.h:49
Polynomials::Polynomial< double > solution_1d
Definition: allen_cahn/polynomial.h:60
virtual void boundary(DoFInfo< dim > &dinfo, IntegrationInfo< dim > &info) const
Definition: allen_cahn/polynomial.h:109
double D
Definition: allen_cahn/polynomial.h:44
void cell_residual(dealii::Vector< number > &result, const dealii::FEValuesBase< dim > &fe, const dealii::VectorSlice< const std::vector< std::vector< dealii::Tensor< 1, dim >>>> &input, double lambda=0., double mu=1.)
Definition: elasticity/integrators.h:23
Definition: allen_cahn/polynomial.h:33
PolynomialError(const Polynomials::Polynomial< double > solution_1d)
Definition: allen_cahn/polynomial.h:144
virtual void cell(DoFInfo< dim > &dinfo, IntegrationInfo< dim > &info) const
Definition: allen_cahn/polynomial.h:153
Definition: allen_cahn/matrix.h:19
Definition: integrator.h:29