3 #ifndef DUNE_PDELAB_ELECTRODYNAMIC_HH
4 #define DUNE_PDELAB_ELECTRODYNAMIC_HH
8 #include<dune/common/exceptions.hh>
9 #include<dune/common/fmatrix.hh>
10 #include<dune/common/fvector.hh>
11 #include<dune/common/static_assert.hh>
12 #include<dune/common/typetraits.hh>
14 #include<dune/geometry/type.hh>
15 #include<dune/geometry/quadraturerules.hh>
37 template<
typename Eps>
66 template<
typename EG,
typename LFS,
typename X,
typename M>
68 const LFS& lfsv, M& mat)
const
71 typedef typename LFS::Traits::FiniteElementType::
72 Traits::Basis::Traits BasisTraits;
74 typedef typename BasisTraits::DomainField DF;
75 static const unsigned dimDLocal = BasisTraits::dimDomainLocal;
77 typedef typename BasisTraits::RangeField RF;
78 typedef typename BasisTraits::Range Range;
79 static const unsigned dimR = BasisTraits::dimRange;
82 dune_static_assert(dimR == 3 || dimR == 2,
83 "Works only in 2D or 3D");
86 "Grids ctype and Finite Elements DomainFieldType must match");
89 typedef Dune::QuadratureRule<DF,dimDLocal> QR;
90 typedef Dune::QuadratureRules<DF,dimDLocal> QRs;
91 Dune::GeometryType gt = eg.geometry().type();
92 const QR& rule = QRs::rule(gt,qorder);
95 for(
typename QR::const_iterator it=rule.begin();
96 it!=rule.end(); ++it) {
98 std::vector<Range> phi(lfsu.size());
99 lfsu.finiteElement().basis().evaluateFunction(it->position(),phi);
102 typename Eps::Traits::RangeType epsval;
103 eps.evaluate(eg.entity(), it->position(), epsval);
105 RF factor = it->weight()
106 * eg.geometry().integrationElement(it->position()) * epsval;
108 for(
unsigned i = 0; i < lfsu.size(); ++i)
109 for(
unsigned j = 0; j < lfsu.size(); ++j)
110 mat.accumulate(lfsv,i,lfsu,j,factor * (phi[i] * phi[j]));
130 template<
typename Mu>
137 template <
unsigned d>
139 static const unsigned dim =
145 template<
typename RF>
147 jacobianToCurl(FieldVector<RF, 1> &curl,
148 const FieldMatrix<RF, 2, 2> &jacobian)
150 curl[0] = jacobian[1][0] - jacobian[0][1];
152 template<
typename RF>
154 jacobianToCurl(FieldVector<RF, 3> &curl,
155 const FieldMatrix<RF, 3, 3> &jacobian)
157 for(
unsigned i = 0; i < 3; ++i)
158 curl[i] = jacobian[(i+2)%3][(i+1)%3] - jacobian[(i+1)%3][(i+2)%3];
184 template<
typename EG,
typename LFS,
typename X,
typename M>
186 const LFS& lfsv, M& mat)
const
189 typedef typename LFS::Traits::FiniteElementType::Traits::Basis::Traits
192 typedef typename BasisTraits::DomainField DF;
193 static const unsigned dimDLocal = BasisTraits::dimDomainLocal;
195 typedef typename BasisTraits::RangeField RF;
196 static const unsigned dimR = BasisTraits::dimRange;
198 typedef typename BasisTraits::Jacobian Jacobian;
202 dune_static_assert(dimR == 3 || dimR == 2,
203 "Works only in 2D or 3D");
206 "Grids ctype and Finite Elements DomainFieldType must match");
209 typedef Dune::QuadratureRule<DF,dimDLocal> QR;
210 typedef Dune::QuadratureRules<DF,dimDLocal> QRs;
211 Dune::GeometryType gt = eg.geometry().type();
212 const QR& rule = QRs::rule(gt,qorder);
215 for(
typename QR::const_iterator it=rule.begin();
216 it!=rule.end(); ++it) {
218 std::vector<Jacobian> J(lfsu.size());
219 lfsu.finiteElement().basis().evaluateJacobian(it->position(),J);
221 std::vector<Curl> rotphi(lfsu.size());
222 for(
unsigned i = 0; i < lfsu.size(); ++i)
223 jacobianToCurl(rotphi[i], J[i]);
226 typename Mu::Traits::RangeType muval;
227 mu.evaluate(eg.entity(), it->position(), muval);
229 RF factor = it->weight()
230 * eg.geometry().integrationElement(it->position()) / muval;
232 for(
unsigned i = 0; i < lfsu.size(); ++i)
233 for(
unsigned j = 0; j < lfsu.size(); ++j)
234 mat.accumulate(lfsv,i,lfsu,j,factor * (rotphi[i] * rotphi[j]));
248 #endif // DUNE_PDELAB_ELECTRODYNAMIC_HH
void jacobian_volume(const EG &eg, const LFS &lfsu, const X &x, const LFS &lfsv, M &mat) const
Definition: electrodynamic.hh:185
Default flags for all local operators.
Definition: flags.hh:18
void jacobian_volume(const EG &eg, const LFS &lfsu, const X &x, const LFS &lfsv, M &mat) const
Definition: electrodynamic.hh:67
Definition: electrodynamic.hh:48
Contruct matrix T for the Electrodynamic operator.
Definition: electrodynamic.hh:38
static const int dim
Definition: adaptivity.hh:82
Implement alpha_volume() based on jacobian_volume()
Definition: defaultimp.hh:611
Electrodynamic_T(const Eps &eps_, int qorder_=2)
Construct an Electrodynamic_T localoperator.
Definition: electrodynamic.hh:58
Electrodynamic_S(const Mu &mu_, int qorder_=2)
Construct an Electrodynamic_S localoperator.
Definition: electrodynamic.hh:176
Contruct matrix S for the Electrodynamic operator.
Definition: electrodynamic.hh:131
Definition: electrodynamic.hh:166
static const unsigned int value
Definition: gridfunctionspace/tags.hh:175
sparsity pattern generator
Definition: pattern.hh:14
const EG & eg
Definition: common/constraints.hh:277
Definition: electrodynamic.hh:46
Definition: electrodynamic.hh:164