4 #ifndef DUNE_PDELAB_LOCALOPERATOR_VECTORWAVE_HH
5 #define DUNE_PDELAB_LOCALOPERATOR_VECTORWAVE_HH
10 #include <dune/common/fvector.hh>
12 #include <dune/geometry/type.hh>
13 #include <dune/geometry/quadraturerules.hh>
27 namespace VectorWave {
43 template<
class Imp,
class GV,
class RF =
double,
class Time_ =
double>
54 typedef FieldVector<DomainField, dimension>
Domain;
59 typedef FieldVector<RangeField, dimension>
Range;
62 typedef typename GV::template Codim<0>::Entity
Element;
73 RangeField
epsilon(
const Element&
e,
const Domain& xl)
const
74 {
return asImp().epsilonGlobal(e.geometry().global(xl)); }
90 RangeField
mu(
const Element&
e,
const Domain& xl)
const
91 {
return asImp().muGlobal(e.geometry().global(xl)); }
98 bool muChanged(Time t1, Time t2)
const {
return true; }
107 const Imp &asImp()
const {
return *
static_cast<const Imp*
>(
this); }
120 template<
class GV,
class RF =
double,
class Time =
double>
122 public Parameters<ConstantParameters<GV, RF, Time>, GV, RF, Time>
130 template<
typename Domain>
133 template<
typename Domain>
163 template<
class Params>
189 R0(Params ¶ms_, std::size_t qorder_ = 0) :
190 params(params_), qorder(qorder_)
193 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
196 const LFSV& lfsv, M& mat)
const
199 typedef typename LFSU::Traits::FiniteElementType FEU;
200 typedef typename LFSV::Traits::FiniteElementType FEV;
202 const FEU &feU = lfsu.finiteElement();
203 const FEV &feV = lfsv.finiteElement();
205 typedef typename FEU::Traits::Basis BasisU;
206 typedef typename FEV::Traits::Basis BasisV;
208 const BasisU &basisU = feU.basis();
209 const BasisV &basisV = feV.basis();
211 typedef typename Params::RangeField RF;
213 typedef typename BasisU::Traits::DomainField DF;
214 static const std::size_t dimDLocal = BasisU::Traits::dimDomainLocal;
216 typedef typename BasisU::Traits::Jacobian JacobianU;
217 typedef typename BasisV::Traits::Jacobian JacobianV;
222 static const J2CU &j2CU = J2CU();
223 static const J2CV &j2CV = J2CV();
225 typedef typename J2CU::Curl CurlU;
226 typedef typename J2CV::Curl CurlV;
228 static_assert(J2CU::dimCurl == J2CV::dimCurl,
"Curl dimension "
229 "of ansatz and test functions must match in "
233 typedef QuadratureRules<DF,dimDLocal> QRs;
234 typedef QuadratureRule<DF,dimDLocal> QR;
235 typedef typename QR::const_iterator QRIterator;
236 GeometryType gt = eg.geometry().type();
237 const QR& rule = QRs::rule(gt,qorder);
239 std::vector<JacobianU> jacobianU(lfsu.size());
240 std::vector<JacobianV> jacobianV(lfsv.size());
242 std::vector<CurlU> curlU(lfsu.size());
243 std::vector<CurlV> curlV(lfsv.size());
246 for(QRIterator it=rule.begin(); it != rule.end(); ++it) {
248 basisU.evaluateJacobian(it->position(), jacobianU);
249 basisV.evaluateJacobian(it->position(), jacobianV);
251 for(std::size_t j = 0; j < lfsu.size(); ++j)
252 j2CU(jacobianU[j], curlU[j]);
253 for(std::size_t i = 0; i < lfsv.size(); ++i)
254 j2CV(jacobianV[i], curlV[i]);
256 RF factor = it->weight()
257 * eg.geometry().integrationElement(it->position())
258 / params.mu(eg.entity(), it->position());
260 for(std::size_t j = 0; j < lfsu.size(); ++j)
261 for(std::size_t i = 0; i < lfsv.size(); ++i)
262 mat(i,j) += factor * (curlU[j] * curlV[i]);
268 params.setTime(time);
292 template<
class Params>
306 R1(Params ¶ms_, std::size_t qorder_ = 0) {}
332 template<
class Params>
358 R2(Params ¶ms_, std::size_t qorder_ = 2)
359 : params(params_), qorder(qorder_)
362 template<
typename EG,
typename LFSU,
typename X,
typename LFSV,
365 const LFSV& lfsv, M& mat)
const
368 typedef typename LFSU::Traits::FiniteElementType FEU;
369 typedef typename LFSV::Traits::FiniteElementType FEV;
371 const FEU &feU = lfsu.finiteElement();
372 const FEV &feV = lfsv.finiteElement();
374 typedef typename FEU::Traits::Basis BasisU;
375 typedef typename FEV::Traits::Basis BasisV;
377 const BasisU &basisU = feU.basis();
378 const BasisV &basisV = feV.basis();
380 typedef typename Params::RangeField RF;
382 typedef typename BasisU::Traits::DomainField DF;
383 static const std::size_t dimDLocal = BasisU::Traits::dimDomainLocal;
385 typedef typename BasisU::Traits::Range RangeU;
386 typedef typename BasisV::Traits::Range RangeV;
389 typedef QuadratureRules<DF,dimDLocal> QRs;
390 typedef QuadratureRule<DF,dimDLocal> QR;
391 typedef typename QR::const_iterator QRIterator;
392 GeometryType gt = eg.geometry().type();
393 const QR& rule = QRs::rule(gt,qorder);
395 std::vector<RangeU> phiU(lfsu.size());
396 std::vector<RangeV> phiV(lfsv.size());
399 for(QRIterator it=rule.begin(); it != rule.end(); ++it) {
401 basisU.evaluateFunction(it->position(), phiU);
402 basisV.evaluateFunction(it->position(), phiV);
404 RF factor = it->weight()
405 * eg.geometry().integrationElement(it->position())
406 * params.epsilon(eg.entity(), it->position());
408 for(std::size_t j = 0; j < lfsu.size(); ++j)
409 for(std::size_t i = 0; i < lfsv.size(); ++i)
410 mat(i,j) += factor * (phiU[j] * phiV[i]);
416 params.setTime(time);
427 #endif // DUNE_PDELAB_LOCALOPERATOR_VECTORWAVE_HH
Local operator for the vector wave problem, no-temporal-derivatives part.
Definition: vectorwave.hh:164
RangeField mu(const Element &e, const Domain &xl) const
evaluate magnetic permeability
Definition: vectorwave.hh:90
Implement alpha_volume() based on jacobian_volume()
Definition: defaultimp.hh:611
static const std::size_t dimension
export dimension (both domain and range)
Definition: vectorwave.hh:49
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: vectorwave.hh:195
void setTime(R t_)
set time for subsequent evaluation
Definition: idefault.hh:104
bool muChanged(Time t1, Time t2) const
Definition: vectorwave.hh:135
const E & e
Definition: interpolate.hh:172
RangeField epsilon(const Element &e, const Domain &xl) const
evaluate dielectric permittivity
Definition: vectorwave.hh:73
RF muGlobal(const Domain &) const
Definition: vectorwave.hh:134
void setTime(typename Params::Time time)
set time on the function object
Definition: vectorwave.hh:415
void setTime(typename Params::Time time)
set time on the parameter object
Definition: vectorwave.hh:267
GV::ctype DomainField
field type of domain
Definition: vectorwave.hh:52
bool muChanged(Time t1, Time t2) const
Whether mu changes between the given time steps.
Definition: vectorwave.hh:98
Default flags for all local operators.
Definition: flags.hh:18
R2(Params ¶ms_, std::size_t qorder_=2)
Construct a local operator object.
Definition: vectorwave.hh:358
sparsity pattern generator
Definition: pattern.hh:13
FieldVector< RangeField, dimension > Range
vector type of range
Definition: vectorwave.hh:59
R1(Params ¶ms_, std::size_t qorder_=0)
Construct a local operator object.
Definition: vectorwave.hh:306
void jacobian_volume(const EG &eg, const LFSU &lfsu, const X &x, const LFSV &lfsv, M &mat) const
Definition: vectorwave.hh:364
Local operator for the vector wave problem, second-temporal-derivatives part.
Definition: vectorwave.hh:333
Time_ Time
export type of temporal values
Definition: vectorwave.hh:46
ConstantParameters(RF epsilon, RF mu)
Definition: vectorwave.hh:128
Definition: vectorwave.hh:348
Definition: vectorwave.hh:179
RF epsilonGlobal(const Domain &) const
Definition: vectorwave.hh:131
void setTime(const Time &time)
set the time for subsequent evaluation
Definition: vectorwave.hh:104
Definition: vectorwave.hh:178
bool epsilonChanged(Time t1, Time t2) const
Definition: vectorwave.hh:132
Definition: adaptivity.hh:27
FieldVector< DomainField, dimension > Domain
vector type of domain
Definition: vectorwave.hh:54
RF RangeField
field type of range
Definition: vectorwave.hh:57
R0(Params ¶ms_, std::size_t qorder_=0)
Construct a local operator object.
Definition: vectorwave.hh:189
Local operator for the vector wave problem, one-temporal-derivative part.
Definition: vectorwave.hh:293
extract the curl of a function from the jacobian of that function
Definition: jacobiantocurl.hh:27
Parameter class interface for the vector wave local operators.
Definition: vectorwave.hh:44
GV::template Codim< 0 >::Entity Element
element type of GridView
Definition: vectorwave.hh:62
bool epsilonChanged(Time t1, Time t2) const
Whether epsilon changes between the given time steps.
Definition: vectorwave.hh:81
const EG & eg
Definition: constraints.hh:280
Default class for additional methods in instationary local operators.
Definition: idefault.hh:89
Definition: vectorwave.hh:347
Homogenous parameter class for the vector wave local operators.
Definition: vectorwave.hh:121