7 #ifndef schroedinger_coulomb_h 8 #define schroedinger_coulomb_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> 15 #include <schroedinger/parameters.h> 32 virtual void cell(MeshWorker::DoFInfo<dim>& dinfo, MeshWorker::IntegrationInfo<dim>& info)
const;
33 virtual void boundary(MeshWorker::DoFInfo<dim>& dinfo,
34 MeshWorker::IntegrationInfo<dim>& info)
const;
35 virtual void face(MeshWorker::DoFInfo<dim>& dinfo1, MeshWorker::DoFInfo<dim>& dinfo2,
36 MeshWorker::IntegrationInfo<dim>& info1,
37 MeshWorker::IntegrationInfo<dim>& info2)
const;
40 SmartPointer<const Parameters, class Coulomb<dim>>
parameters;
47 this->use_boundary =
true;
48 this->use_face =
false;
54 MeshWorker::IntegrationInfo<dim>& info)
const 56 const FEValuesBase<dim>& fe = info.fe_values(0);
57 FullMatrix<double>& M = dinfo.matrix(0,
false).matrix;
60 const unsigned int n_dofs = fe.dofs_per_cell;
61 const unsigned int n_components = fe.get_fe().n_components();
63 for (
unsigned int k = 0; k < fe.n_quadrature_points; ++k)
65 const double potential = -
parameters->depth / fe.quadrature_point(k).norm();
66 const double dx = fe.JxW(k);
67 for (
unsigned int i = 0; i < n_dofs; ++i)
70 for (
unsigned int d = 0; d < n_components; ++d)
72 dx * potential * fe.shape_value_component(i, k, d) * fe.shape_value_component(i, k, d);
76 for (
unsigned int j = i + 1; j < n_dofs; ++j)
79 for (
unsigned int d = 0; d < n_components; ++d)
81 dx * potential * fe.shape_value_component(j, k, d) * fe.shape_value_component(i, k, d);
89 if (dinfo.n_matrices() == 2)
90 L2::mass_matrix(dinfo.matrix(1,
false).matrix, info.fe_values(0));
96 typename MeshWorker::IntegrationInfo<dim>& info)
const 98 const unsigned int deg = info.fe_values(0).get_fe().tensor_degree();
99 Laplace::nitsche_matrix(dinfo.matrix(0,
false).matrix,
101 Laplace::compute_penalty(dinfo, dinfo, deg, deg));
107 MeshWorker::IntegrationInfo<dim>& info1,
108 MeshWorker::IntegrationInfo<dim>& info2)
const 110 const unsigned int deg = info1.fe_values(0).get_fe().tensor_degree();
111 Laplace::ip_matrix(dinfo1.matrix(0,
false).matrix,
112 dinfo1.matrix(0,
true).matrix,
113 dinfo2.matrix(0,
true).matrix,
114 dinfo2.matrix(0,
false).matrix,
117 Laplace::compute_penalty(dinfo1, dinfo2, deg, deg));
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 boundary(MeshWorker::DoFInfo< dim > &dinfo, MeshWorker::IntegrationInfo< dim > &info) const
Definition: coulomb.h:95
virtual void face(MeshWorker::DoFInfo< dim > &dinfo1, MeshWorker::DoFInfo< dim > &dinfo2, MeshWorker::IntegrationInfo< dim > &info1, MeshWorker::IntegrationInfo< dim > &info2) const
Definition: coulomb.h:106
virtual void cell(MeshWorker::DoFInfo< dim > &dinfo, MeshWorker::IntegrationInfo< dim > &info) const
Definition: coulomb.h:53
SmartPointer< const Parameters, class Coulomb< dim > > parameters
Definition: coulomb.h:40
Definition: integrator.h:29