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/typetraits.hh>
13 #include<dune/geometry/type.hh>
14 #include<dune/geometry/quadraturerules.hh>
36 template<
typename Eps>
65 template<
typename EG,
typename LFS,
typename X,
typename M>
67 const LFS& lfsv, M& mat)
const
70 typedef typename LFS::Traits::FiniteElementType::
71 Traits::Basis::Traits BasisTraits;
73 typedef typename BasisTraits::DomainField DF;
74 static const unsigned dimDLocal = BasisTraits::dimDomainLocal;
76 typedef typename BasisTraits::RangeField RF;
77 typedef typename BasisTraits::Range Range;
78 static const unsigned dimR = BasisTraits::dimRange;
81 static_assert(dimR == 3 || dimR == 2,
82 "Works only in 2D or 3D");
85 "Grids ctype and Finite Elements DomainFieldType must match");
88 typedef Dune::QuadratureRule<DF,dimDLocal> QR;
89 typedef Dune::QuadratureRules<DF,dimDLocal> QRs;
90 Dune::GeometryType gt = eg.geometry().type();
91 const QR& rule = QRs::rule(gt,qorder);
94 for(
typename QR::const_iterator it=rule.begin();
95 it!=rule.end(); ++it) {
97 std::vector<Range> phi(lfsu.size());
98 lfsu.finiteElement().basis().evaluateFunction(it->position(),phi);
101 typename Eps::Traits::RangeType epsval;
102 eps.evaluate(eg.entity(), it->position(), epsval);
104 RF factor = it->weight()
105 * eg.geometry().integrationElement(it->position()) * epsval;
107 for(
unsigned i = 0; i < lfsu.size(); ++i)
108 for(
unsigned j = 0; j < lfsu.size(); ++j)
109 mat.accumulate(lfsv,i,lfsu,j,factor * (phi[i] * phi[j]));
129 template<
typename Mu>
136 template <
unsigned d>
138 static const unsigned dim =
144 template<
typename RF>
146 jacobianToCurl(FieldVector<RF, 1> &curl,
147 const FieldMatrix<RF, 2, 2> &jacobian)
149 curl[0] = jacobian[1][0] - jacobian[0][1];
151 template<
typename RF>
153 jacobianToCurl(FieldVector<RF, 3> &curl,
154 const FieldMatrix<RF, 3, 3> &jacobian)
156 for(
unsigned i = 0; i < 3; ++i)
157 curl[i] = jacobian[(i+2)%3][(i+1)%3] - jacobian[(i+1)%3][(i+2)%3];
183 template<
typename EG,
typename LFS,
typename X,
typename M>
185 const LFS& lfsv, M& mat)
const
188 typedef typename LFS::Traits::FiniteElementType::Traits::Basis::Traits
191 typedef typename BasisTraits::DomainField DF;
192 static const unsigned dimDLocal = BasisTraits::dimDomainLocal;
194 typedef typename BasisTraits::RangeField RF;
195 static const unsigned dimR = BasisTraits::dimRange;
197 typedef typename BasisTraits::Jacobian Jacobian;
201 static_assert(dimR == 3 || dimR == 2,
202 "Works only in 2D or 3D");
205 "Grids ctype and Finite Elements DomainFieldType must match");
208 typedef Dune::QuadratureRule<DF,dimDLocal> QR;
209 typedef Dune::QuadratureRules<DF,dimDLocal> QRs;
210 Dune::GeometryType gt = eg.geometry().type();
211 const QR& rule = QRs::rule(gt,qorder);
214 for(
typename QR::const_iterator it=rule.begin();
215 it!=rule.end(); ++it) {
217 std::vector<Jacobian> J(lfsu.size());
218 lfsu.finiteElement().basis().evaluateJacobian(it->position(),J);
220 std::vector<Curl> rotphi(lfsu.size());
221 for(
unsigned i = 0; i < lfsu.size(); ++i)
222 jacobianToCurl(rotphi[i], J[i]);
225 typename Mu::Traits::RangeType muval;
226 mu.evaluate(eg.entity(), it->position(), muval);
228 RF factor = it->weight()
229 * eg.geometry().integrationElement(it->position()) / muval;
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]));
247 #endif // DUNE_PDELAB_ELECTRODYNAMIC_HH
Contruct matrix T for the Electrodynamic operator.
Definition: electrodynamic.hh:37
void jacobian_volume(const EG &eg, const LFS &lfsu, const X &x, const LFS &lfsv, M &mat) const
Definition: electrodynamic.hh:184
Implement alpha_volume() based on jacobian_volume()
Definition: defaultimp.hh:611
Definition: electrodynamic.hh:165
Contruct matrix S for the Electrodynamic operator.
Definition: electrodynamic.hh:130
static const unsigned int value
Definition: gridfunctionspace/tags.hh:175
Default flags for all local operators.
Definition: flags.hh:18
sparsity pattern generator
Definition: pattern.hh:13
Definition: electrodynamic.hh:163
Electrodynamic_T(const Eps &eps_, int qorder_=2)
Construct an Electrodynamic_T localoperator.
Definition: electrodynamic.hh:57
Electrodynamic_S(const Mu &mu_, int qorder_=2)
Construct an Electrodynamic_S localoperator.
Definition: electrodynamic.hh:175
Definition: adaptivity.hh:27
static const int dim
Definition: adaptivity.hh:83
Definition: electrodynamic.hh:47
Definition: electrodynamic.hh:45
void jacobian_volume(const EG &eg, const LFS &lfsu, const X &x, const LFS &lfsv, M &mat) const
Definition: electrodynamic.hh:66
const EG & eg
Definition: constraints.hh:280