diff --git a/dune/mmdg/dg.hh b/dune/mmdg/dg.hh index 395ca14d83faff76c9fe07300860235256bea5b9..a35207028fbe96d35b8325d24ae374b7ccd6bd31 100644 --- a/dune/mmdg/dg.hh +++ b/dune/mmdg/dg.hh @@ -73,11 +73,11 @@ private: //and the indices elemIdxSLE + i + 1, i = 1,...,dim, refer to the //basis function phi_elem,i const int elemIdxSLE = (dim + 1)*elemIdx; - -/* //TODO: can be done outside of the loop? +/* + //TODO: can be done outside of the loop? const Dune::QuadratureRule<double,dim>& rule = Dune::QuadratureRules<double,dim>::rule(geo.type(), problem_.quadratureOrder()); - Dune::FieldVector<double, dim+1> update(0.0); + // Dune::FieldVector<double, dim+1> update(0.0); //NOTE: how are quadrature rules in Dune applied correctly? for (const auto& ip : rule) @@ -88,12 +88,12 @@ private: const double weight = ip.weight(); //quadrature for int_elem q*phi_elem,0 dV - b[elemIdxSLE] += weight * problem_.q(qp); + b[elemIdxSLE] += weight * problem_.q(geo.global(qp)) * geo.integrationElement(qp); //quadrature for int_elem q*phi_elem,i dV for (int i = 0; i < dim; i++) { - b[elemIdxSLE + i + 1] += weight * qp[i] * problem_.q(qp); + b[elemIdxSLE + i + 1] += weight * qp[i] * problem_.q(geo.global(qp)) * geo.integrationElement(qp); } } */