diff --git a/dune/mmdg/mmdg.hh b/dune/mmdg/mmdg.hh index b6d5c222ba74ef10883cadd072bf576aa6ff228e..67251829207bc23cea09c440c40d8a363cde93cb 100644 --- a/dune/mmdg/mmdg.hh +++ b/dune/mmdg/mmdg.hh @@ -354,6 +354,8 @@ private: { const Scalar interfaceUpdate2 = interfaceUpdate1 * (iFrame[j] * qpGlobal); + const Scalar couplingUpdate4 = + 0.5 * weight * betaEvaluation * qpGlobal[j] * qpI; //quadrature for // int_iElem beta * phi_iElem,i * phi_iElem,j ds @@ -366,13 +368,13 @@ private: // -int_iElem beta * phi_iElem,i * avg(phi_elem,j) ds //for elem in {elemIn, elemOut} (*Base::A)[elemInIdxSLE + j + 1][iElemIdxSLE + i + 1] -= - couplingUpdate3; + couplingUpdate4; (*Base::A)[iElemIdxSLE + i + 1][elemInIdxSLE + j + 1] -= - couplingUpdate3; + couplingUpdate4; (*Base::A)[elemOutIdxSLE + j + 1][iElemIdxSLE + i + 1] -= - couplingUpdate3; + couplingUpdate4; (*Base::A)[iElemIdxSLE + i + 1][elemOutIdxSLE + j + 1] -= - couplingUpdate3; + couplingUpdate4; } //quadrature for @@ -440,13 +442,15 @@ private: const auto lambda_Kpar = [&](const Coordinate& qpGlobal, const Scalar weight) { + CoordinateMatrix KEvaluation = Base::problem_.Kparallel(qpGlobal); Scalar dEvaluation2 = Base::problem_.aperture(qpGlobal); dEvaluation2 *= dEvaluation2; + Coordinate KparDotTau_i; //storage for Kparallel * iFrame[i] + for (int i = 0; i < dim-1 ; i++) { - Coordinate KparDotTau_i; - Base::problem_.Kparallel(qpGlobal).mv(iFrame[i], KparDotTau_i); + KEvaluation.mv(iFrame[i], KparDotTau_i); (*Base::A)[iElemIdxSLE + i + 1][iElemIdxSLE + i + 1] += weight * dEvaluation2 * ( KparDotTau_i * iFrame[i] ); @@ -705,7 +709,7 @@ private: weight * qpI * dEvaluation2 * (interfaceUpdate1 - 2.0 * interfaceUpdate2); - for (int j = 0; j < i; i++) + for (int j = 0; j < i; j++) { KEvaluation.mv(iFrame[j], KparDotTau_j);