7 #ifndef __matrix_allen_cahn_h 8 #define __matrix_allen_cahn_h 10 #include <amandus/integrator.h> 11 #include <deal.II/integrators/divergence.h> 12 #include <deal.II/integrators/l2.h> 13 #include <deal.II/integrators/laplace.h> 14 #include <deal.II/meshworker/integration_info.h> 30 virtual void cell(MeshWorker::DoFInfo<dim>& dinfo, MeshWorker::IntegrationInfo<dim>& info)
const;
31 virtual void boundary(MeshWorker::DoFInfo<dim>& dinfo,
32 MeshWorker::IntegrationInfo<dim>& info)
const;
33 virtual void face(MeshWorker::DoFInfo<dim>& dinfo1, MeshWorker::DoFInfo<dim>& dinfo2,
34 MeshWorker::IntegrationInfo<dim>& info1,
35 MeshWorker::IntegrationInfo<dim>& info2)
const;
45 this->use_boundary =
false;
50 Matrix<dim>::cell(MeshWorker::DoFInfo<dim>& dinfo, MeshWorker::IntegrationInfo<dim>& info)
const 52 AssertDimension(dinfo.n_matrices(), 1);
55 if (info.values.size() > 0)
57 AssertDimension(info.values[0][0].size(), info.fe_values(0).n_quadrature_points);
58 std::vector<double> fu(info.fe_values(0).n_quadrature_points);
59 for (
unsigned int k = 0; k < fu.size(); ++k)
61 const double u = info.values[0][0][k];
62 fu[k] = (3. * u * u - 1.);
77 MeshWorker::IntegrationInfo<dim>& info1,
78 MeshWorker::IntegrationInfo<dim>& info2)
const 80 const unsigned int deg = info1.fe_values(0).get_fe().tensor_degree();
81 Laplace::ip_matrix(dinfo1.matrix(0,
false).matrix,
82 dinfo1.matrix(0,
true).matrix,
83 dinfo2.matrix(0,
true).matrix,
84 dinfo2.matrix(0,
false).matrix,
87 Laplace::compute_penalty(dinfo1, dinfo2, deg, deg),
void weighted_mass_matrix(dealii::FullMatrix< double > &M, const dealii::FEValuesBase< dim > &fe, const dealii::TensorFunction< 2, dim > &weight)
Definition: darcy/integrators.h:41
void cell_matrix(dealii::FullMatrix< double > &M, 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: matrix_integrators.h:23
virtual void cell(MeshWorker::DoFInfo< dim > &dinfo, MeshWorker::IntegrationInfo< dim > &info) const
Definition: allen_cahn/matrix.h:50
virtual void face(MeshWorker::DoFInfo< dim > &dinfo1, MeshWorker::DoFInfo< dim > &dinfo2, MeshWorker::IntegrationInfo< dim > &info1, MeshWorker::IntegrationInfo< dim > &info2) const
Definition: allen_cahn/matrix.h:76
virtual void boundary(MeshWorker::DoFInfo< dim > &dinfo, MeshWorker::IntegrationInfo< dim > &info) const
Definition: allen_cahn/matrix.h:70
Matrix(double diffusion)
Definition: allen_cahn/matrix.h:42
double D
Definition: allen_cahn/matrix.h:38
Definition: allen_cahn/matrix.h:25
Definition: allen_cahn/matrix.h:19
Definition: integrator.h:29