7 #ifndef __advection_polynomial_h 8 #define __advection_polynomial_h 10 #include <advection/parameters.h> 11 #include <amandus/integrator.h> 12 #include <deal.II/base/polynomial.h> 13 #include <deal.II/base/tensor.h> 14 #include <deal.II/integrators/advection.h> 15 #include <deal.II/integrators/l2.h> 40 const std::vector<Polynomials::Polynomial<double>> potentials_1d);
42 virtual void cell(DoFInfo<dim>& dinfo, IntegrationInfo<dim>& info)
const;
43 virtual void boundary(DoFInfo<dim>& dinfo, IntegrationInfo<dim>& info)
const;
46 dealii::SmartPointer<const Parameters, class PolynomialRHS<dim>>
parameters;
82 const std::vector<Polynomials::Polynomial<double>> potentials_1d);
84 virtual void cell(DoFInfo<dim>& dinfo, IntegrationInfo<dim>& info)
const;
85 virtual void boundary(DoFInfo<dim>& dinfo, IntegrationInfo<dim>& info)
const;
86 virtual void face(DoFInfo<dim>& dinfo1, DoFInfo<dim>& dinfo2, IntegrationInfo<dim>& info1,
87 IntegrationInfo<dim>& info2)
const;
90 dealii::SmartPointer<const Parameters, class PolynomialError<dim>>
parameters;
98 const std::vector<Polynomials::Polynomial<double>> potentials_1d)
100 , potentials_1d(potentials_1d)
101 , velocity(dim,
std::vector<double>(1))
103 this->use_boundary =
true;
104 this->use_face =
false;
113 std::vector<double> rhs(info.fe_values(0).n_quadrature_points, 0.);
115 std::vector<double> px(2);
116 std::vector<double> py(2);
117 for (
unsigned int k = 0; k < info.fe_values(0).n_quadrature_points; ++k)
119 const double x = info.fe_values(0).quadrature_point(k)(0);
120 const double y = info.fe_values(0).quadrature_point(k)(1);
127 L2::L2(dinfo.vector(0).block(0), info.fe_values(0), rhs);
134 std::vector<double> rhs(info.fe_values(0).n_quadrature_points, 0.);
135 std::vector<double> null(info.fe_values(0).n_quadrature_points, 0.);
137 std::vector<double> px(2);
138 std::vector<double> py(2);
139 for (
unsigned int k = 0; k < info.fe_values(0).n_quadrature_points; ++k)
141 const double x = info.fe_values(0).quadrature_point(k)(0);
142 const double y = info.fe_values(0).quadrature_point(k)(1);
146 rhs[k] = px[0] * py[0];
149 dealii::LocalIntegrators::Advection::upwind_value_residual(
150 dinfo.vector(0).block(0), info.fe_values(0), null, rhs,
velocity);
157 const Parameters& par,
const std::vector<Polynomials::Polynomial<double>>
potentials_1d)
161 this->use_boundary =
false;
162 this->use_face =
false;
171 Assert(dinfo.n_values() >= 2, ExcDimensionMismatch(dinfo.n_values(), 4));
173 std::vector<double> px(3);
174 std::vector<double> py(3);
175 for (
unsigned int k = 0; k < info.fe_values(0).n_quadrature_points; ++k)
177 const double x = info.fe_values(0).quadrature_point(k)(0);
178 const double y = info.fe_values(0).quadrature_point(k)(1);
181 const double dx = info.fe_values(0).JxW(k);
183 Tensor<1, dim> Du = info.gradients[0][0][k];
184 Du[0] -= px[1] * py[0];
185 Du[1] -= px[0] * py[1];
186 double u = info.values[0][0][k];
190 dinfo.value(0) += (u * u) * dx;
192 dinfo.value(1) += (Du * Du) * dx;
195 for (
unsigned int i = 0; i < 2; ++i)
196 dinfo.value(i) = std::sqrt(dinfo.value(i));
208 IntegrationInfo<dim>&)
const virtual void boundary(DoFInfo< dim > &dinfo, IntegrationInfo< dim > &info) const
Definition: advection/polynomial.h:132
std::vector< Polynomials::Polynomial< double > > potentials_1d
Definition: advection/polynomial.h:91
Definition: advection/polynomial.h:78
virtual void face(DoFInfo< dim > &dinfo1, DoFInfo< dim > &dinfo2, IntegrationInfo< dim > &info1, IntegrationInfo< dim > &info2) const
Definition: advection/polynomial.h:207
Definition: advection/matrix.h:17
dealii::SmartPointer< const Parameters, class PolynomialError< dim > > parameters
Definition: advection/polynomial.h:90
virtual void cell(DoFInfo< dim > &dinfo, IntegrationInfo< dim > &info) const
Definition: advection/polynomial.h:111
Definition: advection/polynomial.h:36
dealii::SmartPointer< const Parameters, class PolynomialRHS< dim > > parameters
Definition: advection/polynomial.h:46
std::vector< unsigned int > error_types
Definition: integrator.h:113
PolynomialError(const Parameters &par, const std::vector< Polynomials::Polynomial< double >> potentials_1d)
Definition: advection/polynomial.h:156
std::vector< std::vector< double > > velocity
Definition: advection/polynomial.h:49
std::vector< Polynomials::Polynomial< double > > potentials_1d
Definition: advection/polynomial.h:47
virtual void cell(DoFInfo< dim > &dinfo, IntegrationInfo< dim > &info) const
Definition: advection/polynomial.h:169
virtual void boundary(DoFInfo< dim > &dinfo, IntegrationInfo< dim > &info) const
Definition: advection/polynomial.h:201
Definition: integrator.h:29