From 0ad38afb0f3b33259e7f1ed45e2351b0e1130deb Mon Sep 17 00:00:00 2001
From: Elias Pipping <elias.pipping@fu-berlin.de>
Date: Mon, 4 Feb 2013 16:18:10 +0100
Subject: [PATCH] Implement the L^2(Gamma_F) norm for the state
---
src/one-body-sample.cc | 31 ++++++++++++++++++++++---------
1 file changed, 22 insertions(+), 9 deletions(-)
diff --git a/src/one-body-sample.cc b/src/one-body-sample.cc
index 676c7317..22ef3bf0 100644
--- a/src/one-body-sample.cc
+++ b/src/one-body-sample.cc
@@ -47,11 +47,13 @@
#include <dune/istl/scaledidmatrix.hh>
#include <dune/fufem/assemblers/functionalassembler.hh>
+#include <dune/fufem/assemblers/localassemblers/boundarymassassembler.hh>
#include <dune/fufem/assemblers/localassemblers/l2functionalassembler.hh>
#include <dune/fufem/assemblers/localassemblers/massassembler.hh>
#include <dune/fufem/assemblers/localassemblers/stvenantkirchhoffassembler.hh>
#include <dune/fufem/assemblers/localassemblers/vonmisesstressassembler.hh>
#include <dune/fufem/assemblers/operatorassembler.hh>
+#include <dune/fufem/boundarypatch.hh>
#include <dune/fufem/dunepython.hh>
#include <dune/fufem/functions/basisgridfunction.hh>
#include <dune/fufem/functions/constantfunction.hh>
@@ -134,13 +136,6 @@ template <class FunctionMap> void initPython(FunctionMap &functions) {
.toC<typename FunctionMap::Base>(functions);
}
-template <class VectorType>
-double diff_two_norm(VectorType const &v1, VectorType const &v2) {
- VectorType tmp = v1;
- tmp -= v2;
- return tmp.two_norm();
-}
-
int main(int argc, char *argv[]) {
try {
Dune::Timer timer;
@@ -152,7 +147,9 @@ int main(int argc, char *argv[]) {
typedef Dune::FieldVector<double, dims> SmallVector;
typedef Dune::FieldMatrix<double, dims, dims> SmallMatrix;
+ typedef Dune::FieldMatrix<double, 1, 1> SmallSingletonMatrix;
typedef Dune::BCRSMatrix<SmallMatrix> MatrixType;
+ typedef Dune::BCRSMatrix<SmallSingletonMatrix> SingletonMatrixType;
typedef Dune::BlockVector<SmallVector> VectorType;
typedef Dune::BlockVector<Dune::FieldVector<double, 1>> SingletonVectorType;
@@ -289,7 +286,23 @@ int main(int argc, char *argv[]) {
localStiffness(E, nu);
OperatorAssembler<P1Basis, P1Basis>(p1Basis, p1Basis)
.assemble(localStiffness, stiffnessMatrix);
- };
+ }
+
+ SingletonMatrixType frictionalBoundaryMassMatrix;
+ EnergyNorm<SingletonMatrixType, SingletonVectorType> const stateEnergyNorm(
+ frictionalBoundaryMassMatrix);
+ {
+ BoundaryPatch<GridView> const frictionalBoundary(leafView,
+ frictionalNodes);
+
+ BoundaryMassAssembler<
+ GridType, BoundaryPatch<GridView>, P1Basis::LocalFiniteElement,
+ P1Basis::LocalFiniteElement, SmallSingletonMatrix> const
+ frictionalBoundaryMassAssembler(frictionalBoundary);
+
+ OperatorAssembler<P1Basis, P1Basis>(p1Basis, p1Basis).assemble(
+ frictionalBoundaryMassAssembler, frictionalBoundaryMassMatrix);
+ }
auto const nodalIntegrals =
assemble_frictional<GridType, GridView, SmallVector, P1Basis>(
@@ -412,7 +425,7 @@ int main(int argc, char *argv[]) {
{
SingletonVectorType computed_state;
stateUpdater->extractState(computed_state);
- double const correction = diff_two_norm(computed_state, alpha);
+ double const correction = stateEnergyNorm.diff(computed_state, alpha);
if (state_fpi <= 2 // Let the first two steps pass through unchanged
|| correction < minimalCorrectionReduction * lastCorrection) {
alpha = computed_state;
--
GitLab