From e434ab9a016547ee69db0a3c59f2e8d5083f38ba Mon Sep 17 00:00:00 2001
From: Oliver Sander <sander@igpm.rwth-aachen.de>
Date: Wed, 26 Oct 2011 15:23:24 +0000
Subject: [PATCH] implement frame-invariance test for the energy
[[Imported from SVN: r8055]]
---
test/cosseratenergytest.cc | 121 ++++++++++++++++++++++++++++++++++++-
1 file changed, 119 insertions(+), 2 deletions(-)
diff --git a/test/cosseratenergytest.cc b/test/cosseratenergytest.cc
index f6b4ae3c..0e0bd8d4 100644
--- a/test/cosseratenergytest.cc
+++ b/test/cosseratenergytest.cc
@@ -4,6 +4,8 @@
#include <dune/localfunctions/lagrange/pqkfactory.hh>
+#include <dune/fufem/functions/constantfunction.hh>
+
#include <dune/gfe/rigidbodymotion.hh>
#include <dune/gfe/cosseratenergystiffness.hh>
@@ -102,10 +104,115 @@ void testDerivativeOfRotationMatrix(const std::vector<TargetSpace>& corners)
}
}
+//////////////////////////////////////////////////////////////////////////////////////
+// Test invariance of the energy functional under rotations
+//////////////////////////////////////////////////////////////////////////////////////
+
+template <class GridType>
+void testEnergy(const GridType* grid, const std::vector<TargetSpace>& coefficients) {
+
+ PQkLocalFiniteElementCache<double,double,GridType::dimension,1> feCache;
+ typedef typename PQkLocalFiniteElementCache<double,double,GridType::dimension,1>::FiniteElementType LocalFiniteElement;
+
+ //LocalGeodesicFEFunction<domainDim,double,LocalFiniteElement,TargetSpace> f(feCache.get(element),corners);
+
+ ParameterTree materialParameters;
+ materialParameters["thickness"] = "0.1";
+ materialParameters["mu"] = "3.8462e+05";
+ materialParameters["lambda"] = "2.7149e+05";
+ materialParameters["mu_c"] = "3.8462e+05";
+ materialParameters["L_c"] = "0.1";
+ materialParameters["q"] = "2.5";
+
+ ConstantFunction<Dune::FieldVector<double,GridType::dimension>, Dune::FieldVector<double,3> > zeroFunction(Dune::FieldVector<double,3>(0));
+
+ CosseratEnergyLocalStiffness<typename GridType::LeafGridView,LocalFiniteElement,3> assembler(materialParameters,
+ NULL,
+ &zeroFunction);
+ std::vector<TargetSpace> rotatedCoefficients(coefficients.size());
+
+ std::vector<Rotation<3> > testRotations;
+ ValueFactory<Rotation<3> >::get(testRotations);
+ for (size_t i=0; i<testRotations.size(); i++) {
+ /////////////////////////////////////////////////////////////////////////
+ // Multiply the given configuration by the test rotation.
+ // The energy should remain unchanged.
+ /////////////////////////////////////////////////////////////////////////
+ FieldMatrix<double,3,3> matrix;
+ testRotations[i].matrix(matrix);
+
+ for (size_t j=0; j<coefficients.size(); j++) {
+ FieldVector<double,3> tmp;
+ matrix.mv(coefficients[j].r, tmp);
+ rotatedCoefficients[j].r = tmp;
+
+ rotatedCoefficients[j].q = testRotations[i].mult(rotatedCoefficients[j].q);
+ }
+
+ std::cout << "energy: " << assembler.energy(*grid->template leafbegin<0>(),
+ feCache.get(grid->template leafbegin<0>()->type()),
+ rotatedCoefficients) << std::endl;
-int main(int argc, char** argv)
+ }
+
+}
+
+
+template <int domainDim>
+void testFrameInvariance()
+{
+ // ////////////////////////////////////////////////////////
+ // Make a test grid consisting of a single simplex
+ // ////////////////////////////////////////////////////////
+
+ typedef UGGrid<domainDim> GridType;
+
+ GridFactory<GridType> factory;
+
+ FieldVector<double,dim> pos(0);
+ factory.insertVertex(pos);
+
+ for (int i=0; i<domainDim+1; i++) {
+ pos = 0;
+ pos[i] = 1;
+ factory.insertVertex(pos);
+ }
+
+ std::vector<unsigned int> v(domainDim+1);
+ for (int i=0; i<domainDim+1; i++)
+ v[i] = i;
+ factory.insertElement(GeometryType(GeometryType::simplex,dim), v);
+
+ const GridType* grid = factory.createGrid();
+
+ // //////////////////////////////////////////////////////////
+ // Test whether the energy is invariant under isometries
+ // //////////////////////////////////////////////////////////
+
+ std::vector<TargetSpace> testPoints;
+ ValueFactory<TargetSpace>::get(testPoints);
+
+ // Set up elements of SE(3)
+ std::vector<TargetSpace> coefficients(dim+1);
+
+ MultiIndex index(dim+1, testPoints.size());
+ int numIndices = index.cycle();
+
+ for (int i=0; i<numIndices; i++, ++index) {
+
+ for (int j=0; j<dim+1; j++)
+ coefficients[j] = testPoints[index[j]];
+
+ testEnergy<GridType>(grid, coefficients);
+
+ }
+
+}
+
+
+int main(int argc, char** argv) try
{
const int domainDim = 2;
std::cout << " --- Testing Rotation<3>, domain dimension: " << domainDim << " ---" << std::endl;
@@ -131,4 +238,14 @@ int main(int argc, char** argv)
}
-}
+ //////////////////////////////////////////////////////////////////////////////////////
+ // Test invariance of the energy functional under rotations
+ //////////////////////////////////////////////////////////////////////////////////////
+
+ testFrameInvariance<2>();
+
+} catch (Exception e) {
+
+ std::cout << e << std::endl;
+
+ }
--
GitLab