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 {
37 template<
class Params>
40 const Params *params_;
43 DUNE_DEPRECATED_MSG(
"You are using Dune::PDELab::Electrodynamic_T or " 44 "Dune::PDELab::Electrodynamic_S with an old-style " 45 "parameter class for eps/mu. Please use a " 46 "callable (function object, (generic) lambda, or " 47 "a function pointer) instead.")
48 Functor(const Params ¶ms) : params_(¶ms) {}
50 template<
class Element,
class Pos>
51 typename Params::Traits::RangeType
54 typename Params::Traits::RangeType result;
55 params_->evaluate(e, xl, result);
60 template<
class Params,
class =
void>
63 template<
class Params>
66 std::enable_if_t<sizeof(typename Params::Traits::RangeType)> > :
70 template<
class Params>
79 constexpr std::size_t
dimOfCurl(std::size_t dimOfSpace)
86 throw std::invalid_argument(
"Only applicable for dimensions 1-3");
91 const FieldMatrix<RF, 2, 2> &jacobian)
93 curl[0] = jacobian[1][0] - jacobian[0][1];
97 const FieldMatrix<RF, 3, 3> &jacobian)
99 for(
unsigned i = 0; i < 3; ++i)
100 curl[i] = jacobian[(i+2)%3][(i+1)%3] - jacobian[(i+1)%3][(i+2)%3];
115 template<
typename Eps>
121 static constexpr
bool oldstyle =
127 static constexpr
bool doPatternVolume =
true;
128 static constexpr
bool doAlphaVolume =
true;
144 eps_(eps), qorder_(qorder)
148 eps_(
std::move(eps)), qorder_(qorder)
154 template<
typename EG,
typename LFS,
typename X,
typename M>
156 const LFS& lfsv, M& mat)
const 159 typename LFS::Traits::FiniteElementType::Traits::Basis::Traits;
162 static constexpr
unsigned dimR = BasisTraits::dimRange;
163 static_assert(dimR == 3 || dimR == 2,
"Works only in 2D or 3D");
165 using ctype =
typename EG::Geometry::ctype;
166 using DF =
typename BasisTraits::DomainField;
168 "Finite Elements DomainFieldType must match");
170 using Range =
typename BasisTraits::Range;
171 std::vector<Range> phi(lfsu.size());
177 lfsu.finiteElement().basis().evaluateFunction(qp.position(),phi);
180 auto factor = qp.weight()
181 * eg.geometry().integrationElement(qp.position())
182 * eps_(eg.entity(), qp.position());
184 for(
unsigned i = 0; i < lfsu.size(); ++i)
185 for(
unsigned j = 0; j < lfsu.size(); ++j)
186 mat.accumulate(lfsv,i,lfsu,j,factor * (phi[i] * phi[j]));
203 return { std::forward<Eps>(eps), qorder };
217 template<
typename Mu>
223 static constexpr
bool oldstyle =
229 static constexpr
bool doPatternVolume =
true;
230 static constexpr
bool doAlphaVolume =
true;
246 mu_(mu), qorder_(qorder)
250 mu_(
std::move(mu)), qorder_(qorder)
256 template<
typename EG,
typename LFS,
typename X,
typename M>
258 const LFS& lfsv, M& mat)
const 264 typename LFS::Traits::FiniteElementType::Traits::Basis::Traits;
267 static constexpr
unsigned dimR = BasisTraits::dimRange;
268 static_assert(dimR == 3 || dimR == 2,
"Works only in 2D or 3D");
270 using ctype =
typename EG::Geometry::ctype;
271 using DF =
typename BasisTraits::DomainField;
273 "Finite Elements DomainFieldType must match");
275 using Jacobian =
typename BasisTraits::Jacobian;
276 std::vector<Jacobian> J(lfsu.size());
278 using RF =
typename BasisTraits::RangeField;
279 using Curl = FieldVector<RF, dimOfCurl(dimR)>;
280 std::vector<Curl> rotphi(lfsu.size());
286 lfsu.finiteElement().basis().evaluateJacobian(qp.position(),J);
287 for(
unsigned i = 0; i < lfsu.size(); ++i)
291 auto factor = qp.weight()
292 * eg.geometry().integrationElement(qp.position())
293 / mu_(eg.entity(), qp.position());
295 for(
unsigned i = 0; i < lfsu.size(); ++i)
296 for(
unsigned j = 0; j < lfsu.size(); ++j)
297 mat.accumulate(lfsv,i,lfsu,j,factor * (rotphi[i] * rotphi[j]));
315 return { std::forward<Mu>(mu), qorder };
322 #endif // DUNE_PDELAB_LOCALOPERATOR_ELECTRODYNAMIC_HH const Entity & e
Definition: localfunctionspace.hh:111
Electrodynamic_S< std::decay_t< Mu > > makeLocalOperatorEdynS(Mu &&mu, int qorder=2)
construct an Electrodynamic_S operator
Definition: electrodynamic.hh:313
void jacobian_volume(const EG &eg, const LFS &lfsu, const X &x, const LFS &lfsv, M &mat) const
Definition: electrodynamic.hh:257
void jacobian_volume(const EG &eg, const LFS &lfsu, const X &x, const LFS &lfsv, M &mat) const
Definition: electrodynamic.hh:155
Electrodynamic_S(Mu &&mu, int qorder=2)
Definition: electrodynamic.hh:249
Construct matrix T for the Electrodynamic operator.
Definition: electrodynamic.hh:116
Params::Traits::RangeType operator()(const Element &e, const Pos &xl) const
Definition: electrodynamic.hh:52
Electrodynamic_S(const Mu &mu, int qorder=2)
Construct an Electrodynamic_S localoperator.
Definition: electrodynamic.hh:245
For backward compatibility – Do not use this!
Definition: adaptivity.hh:27
void jacobianToCurl(FieldVector< RF, 1 > &curl, const FieldMatrix< RF, 2, 2 > &jacobian)
Definition: electrodynamic.hh:90
Definition: electrodynamic.hh:38
Electrodynamic_T< std::decay_t< Eps > > makeLocalOperatorEdynT(Eps &&eps, int qorder=2)
construct an Electrodynamic_T operator
Definition: electrodynamic.hh:201
Electrodynamic_T(Eps &&eps, int qorder=2)
Definition: electrodynamic.hh:147
Contruct matrix S for the Electrodynamic operator.
Definition: electrodynamic.hh:218
sparsity pattern generator
Definition: pattern.hh:13
Definition: electrodynamic.hh:61
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
Implement alpha_volume() based on jacobian_volume()
Definition: numericalresidual.hh:35
Electrodynamic_T(const Eps &eps, int qorder=2)
Construct an Electrodynamic_T localoperator.
Definition: electrodynamic.hh:143
constexpr std::size_t dimOfCurl(std::size_t dimOfSpace)
Definition: electrodynamic.hh:79
Default flags for all local operators.
Definition: flags.hh:18
static const unsigned int value
Definition: gridfunctionspace/tags.hh:139
std::conditional_t< IsOldstyleParams< Params >::value, Functor< Params >, Params > ConstRefOrFunctor
Definition: electrodynamic.hh:73