3 #ifndef DUNE_PDELAB_LOCALOPERATOR_ELECTRODYNAMIC_HH 4 #define DUNE_PDELAB_LOCALOPERATOR_ELECTRODYNAMIC_HH 12 #include <dune/common/deprecated.hh> 13 #include <dune/common/fmatrix.hh> 14 #include <dune/common/fvector.hh> 16 #include <dune/geometry/quadraturerules.hh> 28 namespace ElectrodynamicImpl {
35 constexpr std::size_t
dimOfCurl(std::size_t dimOfSpace)
42 throw std::invalid_argument(
"Only applicable for dimensions 1-3");
47 const FieldMatrix<RF, 2, 2> &jacobian)
49 curl[0] = jacobian[1][0] - jacobian[0][1];
53 const FieldMatrix<RF, 3, 3> &jacobian)
55 for(
unsigned i = 0; i < 3; ++i)
56 curl[i] = jacobian[(i+2)%3][(i+1)%3] - jacobian[(i+1)%3][(i+2)%3];
71 template<
typename Eps>
81 static constexpr
bool doPatternVolume =
true;
82 static constexpr
bool doAlphaVolume =
true;
90 eps_(eps), qorder_(qorder)
94 eps_(
std::move(eps)), qorder_(qorder)
100 template<
typename EG,
typename LFS,
typename X,
typename M>
102 const LFS& lfsv, M& mat)
const 105 typename LFS::Traits::FiniteElementType::Traits::Basis::Traits;
108 static constexpr
unsigned dimR = BasisTraits::dimRange;
109 static_assert(dimR == 3 || dimR == 2,
"Works only in 2D or 3D");
111 using ctype =
typename EG::Geometry::ctype;
112 using DF =
typename BasisTraits::DomainField;
114 "Finite Elements DomainFieldType must match");
116 using Range =
typename BasisTraits::Range;
117 std::vector<Range> phi(lfsu.size());
123 lfsu.finiteElement().basis().evaluateFunction(qp.position(),phi);
126 auto factor = qp.weight()
127 * eg.geometry().integrationElement(qp.position())
128 * eps_(eg.entity(), qp.position());
130 for(
unsigned i = 0; i < lfsu.size(); ++i)
131 for(
unsigned j = 0; j < lfsu.size(); ++j)
132 mat.accumulate(lfsv,i,lfsu,j,factor * (phi[i] * phi[j]));
149 return { std::forward<Eps>(eps), qorder };
163 template<
typename Mu>
173 static constexpr
bool doPatternVolume =
true;
174 static constexpr
bool doAlphaVolume =
true;
182 mu_(mu), qorder_(qorder)
186 mu_(
std::move(mu)), qorder_(qorder)
192 template<
typename EG,
typename LFS,
typename X,
typename M>
194 const LFS& lfsv, M& mat)
const 200 typename LFS::Traits::FiniteElementType::Traits::Basis::Traits;
203 static constexpr
unsigned dimR = BasisTraits::dimRange;
204 static_assert(dimR == 3 || dimR == 2,
"Works only in 2D or 3D");
206 using ctype =
typename EG::Geometry::ctype;
207 using DF =
typename BasisTraits::DomainField;
209 "Finite Elements DomainFieldType must match");
211 using Jacobian =
typename BasisTraits::Jacobian;
212 std::vector<Jacobian> J(lfsu.size());
214 using RF =
typename BasisTraits::RangeField;
215 using Curl = FieldVector<RF, dimOfCurl(dimR)>;
216 std::vector<Curl> rotphi(lfsu.size());
222 lfsu.finiteElement().basis().evaluateJacobian(qp.position(),J);
223 for(
unsigned i = 0; i < lfsu.size(); ++i)
227 auto factor = qp.weight()
228 * eg.geometry().integrationElement(qp.position())
229 / mu_(eg.entity(), qp.position());
231 for(
unsigned i = 0; i < lfsu.size(); ++i)
232 for(
unsigned j = 0; j < lfsu.size(); ++j)
233 mat.accumulate(lfsv,i,lfsu,j,factor * (rotphi[i] * rotphi[j]));
251 return { std::forward<Mu>(mu), qorder };
258 #endif // DUNE_PDELAB_LOCALOPERATOR_ELECTRODYNAMIC_HH Default flags for all local operators.
Definition: flags.hh:18
Construct matrix T for the Electrodynamic operator.
Definition: electrodynamic.hh:72
Electrodynamic_T(Eps &&eps, int qorder=2)
Definition: electrodynamic.hh:93
Electrodynamic_T< std::decay_t< Eps > > makeLocalOperatorEdynT(Eps &&eps, int qorder=2)
construct an Electrodynamic_T operator
Definition: electrodynamic.hh:147
sparsity pattern generator
Definition: pattern.hh:13
constexpr std::size_t dimOfCurl(std::size_t dimOfSpace)
Definition: electrodynamic.hh:35
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
Electrodynamic_S(const Mu &mu, int qorder=2)
Construct an Electrodynamic_S localoperator.
Definition: electrodynamic.hh:181
Electrodynamic_S(Mu &&mu, int qorder=2)
Definition: electrodynamic.hh:185
QuadratureRuleWrapper< QuadratureRule< typename Geometry::ctype, Geometry::mydimension > > quadratureRule(const Geometry &geo, std::size_t order, QuadratureType::Enum quadrature_type=QuadratureType::GaussLegendre)
Returns a quadrature rule for the given geometry.
Definition: quadraturerules.hh:117
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Electrodynamic_S< std::decay_t< Mu > > makeLocalOperatorEdynS(Mu &&mu, int qorder=2)
construct an Electrodynamic_S operator
Definition: electrodynamic.hh:249
Implement alpha_volume() based on jacobian_volume()
Definition: numericalresidual.hh:35
Contruct matrix S for the Electrodynamic operator.
Definition: electrodynamic.hh:164
void jacobian_volume(const EG &eg, const LFS &lfsu, const X &x, const LFS &lfsv, M &mat) const
Definition: electrodynamic.hh:101
void jacobian_volume(const EG &eg, const LFS &lfsu, const X &x, const LFS &lfsv, M &mat) const
Definition: electrodynamic.hh:193
void jacobianToCurl(FieldVector< RF, 1 > &curl, const FieldMatrix< RF, 2, 2 > &jacobian)
Definition: electrodynamic.hh:46
Electrodynamic_T(const Eps &eps, int qorder=2)
Construct an Electrodynamic_T localoperator.
Definition: electrodynamic.hh:89