Skip to content
Snippets Groups Projects
Commit 4a41903a authored by AlexanderMüller's avatar AlexanderMüller
Browse files

augmented existing tests to also test for ProductManifold

parent 9a3bf4dd
No related branches found
No related tags found
1 merge request!70added ProductManifold<...>
This commit is part of merge request !70. Comments created here will be created in the context of that merge request.
Hide whitespace changes
Inline Side-by-side
test/averagedistanceassemblertest.cc
+ 21
− 0
View file @ 4a41903a
...@@ -9,6 +9,7 @@ ...@@ -9,6 +9,7 @@
#include <dune/gfe/rotation.hh> #include <dune/gfe/rotation.hh>
#include <dune/gfe/realtuple.hh> #include <dune/gfe/realtuple.hh>
#include <dune/gfe/unitvector.hh> #include <dune/gfe/unitvector.hh>
#include <dune/gfe/productmanifold.hh>
#include <dune/gfe/localgeodesicfefunction.hh> #include <dune/gfe/localgeodesicfefunction.hh>
...@@ -229,9 +230,29 @@ void testRotations() ...@@ -229,9 +230,29 @@ void testRotations()
} }
void testProductManifold()
{
typedef Dune::GFE::ProductManifold<RealTuple<double,5>,UnitVector<double,3>, Rotation<double,3>> TargetSpace;
std::vector<TargetSpace> corners(dim+1);
std::generate(corners.begin(), corners.end(), []() {
return Dune::GFE::randomFieldVector<typename TargetSpace::field_type,TargetSpace::CoordinateType::dimension>(0.9,1.1);
});
TargetSpace argument = corners[0];
testWeightSet(corners, argument);
argument = corners[1];
testWeightSet(corners, argument);
argument = corners[2];
testWeightSet(corners, argument);
}
int main() int main()
{ {
testRealTuples(); testRealTuples();
testUnitVectors(); testUnitVectors();
testRotations(); testRotations();
testProductManifold();
} }
test/localgeodesicfefunctiontest.cc
+ 7
− 0
View file @ 4a41903a
...@@ -14,6 +14,7 @@ ...@@ -14,6 +14,7 @@
#include <dune/gfe/rotation.hh> #include <dune/gfe/rotation.hh>
#include <dune/gfe/realtuple.hh> #include <dune/gfe/realtuple.hh>
#include <dune/gfe/unitvector.hh> #include <dune/gfe/unitvector.hh>
#include <dune/gfe/productmanifold.hh>
#include <dune/gfe/localgeodesicfefunction.hh> #include <dune/gfe/localgeodesicfefunction.hh>
#include "multiindex.hh" #include "multiindex.hh"
...@@ -303,6 +304,8 @@ int main() ...@@ -303,6 +304,8 @@ int main()
test<UnitVector<double,3>,1>(GeometryTypes::simplex(1)); test<UnitVector<double,3>,1>(GeometryTypes::simplex(1));
test<Rotation<double,3>,1>(GeometryTypes::simplex(1)); test<Rotation<double,3>,1>(GeometryTypes::simplex(1));
test<RigidBodyMotion<double,3>,1>(GeometryTypes::simplex(1)); test<RigidBodyMotion<double,3>,1>(GeometryTypes::simplex(1));
typedef Dune::GFE::ProductManifold<RealTuple<double,1>,Rotation<double,3>,UnitVector<double,2>> CrazyManifold;
test<CrazyManifold,1>(GeometryTypes::simplex(1));
//////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////
// Test functions on 2d simplex elements // Test functions on 2d simplex elements
...@@ -313,6 +316,8 @@ int main() ...@@ -313,6 +316,8 @@ int main()
test<UnitVector<double,3>,2>(GeometryTypes::simplex(2)); test<UnitVector<double,3>,2>(GeometryTypes::simplex(2));
test<Rotation<double,3>,2>(GeometryTypes::simplex(2)); test<Rotation<double,3>,2>(GeometryTypes::simplex(2));
test<RigidBodyMotion<double,3>,2>(GeometryTypes::simplex(2)); test<RigidBodyMotion<double,3>,2>(GeometryTypes::simplex(2));
typedef Dune::GFE::ProductManifold<RealTuple<double,1>,Rotation<double,3>,UnitVector<double,2>> CrazyManifold;
test<CrazyManifold,2>(GeometryTypes::simplex(2));
//////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////
// Test functions on 2d quadrilateral elements // Test functions on 2d quadrilateral elements
...@@ -323,5 +328,7 @@ int main() ...@@ -323,5 +328,7 @@ int main()
test<UnitVector<double,3>,2>(GeometryTypes::cube(2)); test<UnitVector<double,3>,2>(GeometryTypes::cube(2));
test<Rotation<double,3>,2>(GeometryTypes::cube(2)); test<Rotation<double,3>,2>(GeometryTypes::cube(2));
test<RigidBodyMotion<double,3>,2>(GeometryTypes::cube(2)); test<RigidBodyMotion<double,3>,2>(GeometryTypes::cube(2));
typedef Dune::GFE::ProductManifold<RealTuple<double,1>,Rotation<double,3>,UnitVector<double,2>> CrazyManifold;
test<CrazyManifold,2>(GeometryTypes::cube(2));
} }
test/localgfetestfunctiontest.cc
+ 3
− 0
View file @ 4a41903a
...@@ -10,6 +10,7 @@ ...@@ -10,6 +10,7 @@
#include <dune/gfe/rotation.hh> #include <dune/gfe/rotation.hh>
#include <dune/gfe/realtuple.hh> #include <dune/gfe/realtuple.hh>
#include <dune/gfe/unitvector.hh> #include <dune/gfe/unitvector.hh>
#include <dune/gfe/productmanifold.hh>
#include <dune/gfe/localgeodesicfefunction.hh> #include <dune/gfe/localgeodesicfefunction.hh>
#include <dune/gfe/localgfetestfunctionbasis.hh> #include <dune/gfe/localgfetestfunctionbasis.hh>
...@@ -79,6 +80,8 @@ int main() try ...@@ -79,6 +80,8 @@ int main() try
test<Rotation<double,3>, 1>(); test<Rotation<double,3>, 1>();
test<Rotation<double,3>, 2>(); test<Rotation<double,3>, 2>();
typedef Dune::GFE::ProductManifold<RealTuple<double,1>,Rotation<double,3>,UnitVector<double,2>> CrazyManifold;
test<CrazyManifold, 2>();
} catch (Exception& e) { } catch (Exception& e) {
std::cout << e << std::endl; std::cout << e << std::endl;
......
test/localprojectedfefunctiontest.cc
+ 21
− 0
View file @ 4a41903a
...@@ -15,6 +15,7 @@ ...@@ -15,6 +15,7 @@
#include <dune/gfe/rotation.hh> #include <dune/gfe/rotation.hh>
#include <dune/gfe/realtuple.hh> #include <dune/gfe/realtuple.hh>
#include <dune/gfe/unitvector.hh> #include <dune/gfe/unitvector.hh>
#include <dune/gfe/productmanifold.hh>
#include <dune/gfe/localprojectedfefunction.hh> #include <dune/gfe/localprojectedfefunction.hh>
#include "multiindex.hh" #include "multiindex.hh"
...@@ -105,6 +106,20 @@ void testDerivativeTangentiality(const RigidBodyMotion<double,3>& x, ...@@ -105,6 +106,20 @@ void testDerivativeTangentiality(const RigidBodyMotion<double,3>& x,
{ {
} }
// the columns of the derivative must be tangential to the manifold
template <int domainDim, int vectorDim,typename ...TargetSpaces>
void testDerivativeTangentiality(const Dune::GFE::ProductManifold<TargetSpaces...>& x,
const FieldMatrix<double,vectorDim,domainDim>& derivative)
{
size_t posHelper=0;
using namespace Dune::Hybrid;
forEach(integralRange(Dune::Hybrid::size(x)), [&](auto&& i) {
using Manifold = std::remove_reference_t<decltype(x[i])>;
testDerivativeTangentiality(x[i],Dune::GFE::blockAt<Manifold::embeddedDim,domainDim>( derivative,posHelper,0));
posHelper +=Manifold::embeddedDim;
});
}
/** \brief Test whether interpolation is invariant under permutation of the simplex vertices /** \brief Test whether interpolation is invariant under permutation of the simplex vertices
* \todo Implement this for all dimensions * \todo Implement this for all dimensions
*/ */
...@@ -255,6 +270,8 @@ int main() ...@@ -255,6 +270,8 @@ int main()
test<UnitVector<double,3>,1>(GeometryTypes::line); test<UnitVector<double,3>,1>(GeometryTypes::line);
test<Rotation<double,3>,1>(GeometryTypes::line); test<Rotation<double,3>,1>(GeometryTypes::line);
test<RigidBodyMotion<double,3>,1>(GeometryTypes::line); test<RigidBodyMotion<double,3>,1>(GeometryTypes::line);
typedef Dune::GFE::ProductManifold<RealTuple<double,1>,Rotation<double,3>,UnitVector<double,2>> CrazyManifold;
test<CrazyManifold, 1>(GeometryTypes::line);
//////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////
// Test functions on 2d simplex elements // Test functions on 2d simplex elements
...@@ -265,6 +282,8 @@ int main() ...@@ -265,6 +282,8 @@ int main()
test<UnitVector<double,3>,2>(GeometryTypes::triangle); test<UnitVector<double,3>,2>(GeometryTypes::triangle);
test<Rotation<double,3>,2>(GeometryTypes::triangle); test<Rotation<double,3>,2>(GeometryTypes::triangle);
test<RigidBodyMotion<double,3>,2>(GeometryTypes::triangle); test<RigidBodyMotion<double,3>,2>(GeometryTypes::triangle);
typedef Dune::GFE::ProductManifold<RealTuple<double,1>,Rotation<double,3>,UnitVector<double,2>> CrazyManifold;
test<CrazyManifold, 2>(GeometryTypes::triangle);
//////////////////////////////////////////////////////////////// ////////////////////////////////////////////////////////////////
// Test functions on 2d quadrilateral elements // Test functions on 2d quadrilateral elements
...@@ -275,5 +294,7 @@ int main() ...@@ -275,5 +294,7 @@ int main()
test<UnitVector<double,3>,2>(GeometryTypes::quadrilateral); test<UnitVector<double,3>,2>(GeometryTypes::quadrilateral);
test<Rotation<double,3>,2>(GeometryTypes::quadrilateral); test<Rotation<double,3>,2>(GeometryTypes::quadrilateral);
test<RigidBodyMotion<double,3>,2>(GeometryTypes::quadrilateral); test<RigidBodyMotion<double,3>,2>(GeometryTypes::quadrilateral);
typedef Dune::GFE::ProductManifold<RealTuple<double,1>,Rotation<double,3>,UnitVector<double,2>> CrazyManifold;
test<CrazyManifold, 2>(GeometryTypes::quadrilateral);
} }
test/targetspacetest.cc
+ 16
− 11
View file @ 4a41903a
...@@ -3,6 +3,7 @@ ...@@ -3,6 +3,7 @@
#include <dune/gfe/unitvector.hh> #include <dune/gfe/unitvector.hh>
#include <dune/gfe/realtuple.hh> #include <dune/gfe/realtuple.hh>
#include <dune/gfe/rotation.hh> #include <dune/gfe/rotation.hh>
#include <dune/gfe/productmanifold.hh>
#include <dune/gfe/hyperbolichalfspacepoint.hh> #include <dune/gfe/hyperbolichalfspacepoint.hh>
#include "valuefactory.hh" #include "valuefactory.hh"
...@@ -138,8 +139,8 @@ void testDerivativeOfSquaredDistance(const TargetSpace& a, const TargetSpace& b) ...@@ -138,8 +139,8 @@ void testDerivativeOfSquaredDistance(const TargetSpace& a, const TargetSpace& b)
// transform into embedded coordinates // transform into embedded coordinates
typename TargetSpace::EmbeddedTangentVector d2_fd_embedded; typename TargetSpace::EmbeddedTangentVector d2_fd_embedded;
B.mtv(d2_fd,d2_fd_embedded); B.mtv(d2_fd,d2_fd_embedded);
if ( (d2 - d2_fd_embedded).infinity_norm() > 100*eps ) { if ( (d2 - d2_fd_embedded).infinity_norm() > 200*eps ) {
std::cout << className(a) << ": Analytical gradient does not match fd approximation." << std::endl; std::cout << className(a) << ": Analytical gradient does not match fd approximation." << std::endl;
std::cout << "d2 Analytical: " << d2 << std::endl; std::cout << "d2 Analytical: " << d2 << std::endl;
std::cout << "d2 FD : " << d2_fd << std::endl; std::cout << "d2 FD : " << d2_fd << std::endl;
...@@ -156,14 +157,14 @@ void testHessianOfSquaredDistance(const TargetSpace& a, const TargetSpace& b) ...@@ -156,14 +157,14 @@ void testHessianOfSquaredDistance(const TargetSpace& a, const TargetSpace& b)
/////////////////////////////////////////////////////////////////// ///////////////////////////////////////////////////////////////////
// Test second derivative with respect to second argument // Test second derivative with respect to second argument
/////////////////////////////////////////////////////////////////// ///////////////////////////////////////////////////////////////////
FieldMatrix<double,embeddedDim,embeddedDim> d2d2 = TargetSpace::secondDerivativeOfDistanceSquaredWRTSecondArgument(a, b); FieldMatrix<double,embeddedDim,embeddedDim> d2d2 = (TargetSpace::secondDerivativeOfDistanceSquaredWRTSecondArgument(a, b)).matrix();
// finite-difference approximation // finite-difference approximation
FieldMatrix<double,embeddedDim,embeddedDim> d2d2_fd = getSecondDerivativeOfSecondArgumentFD<TargetSpace,embeddedDim>(a,b); FieldMatrix<double,embeddedDim,embeddedDim> d2d2_fd = getSecondDerivativeOfSecondArgumentFD<TargetSpace,embeddedDim>(a,b);
FieldMatrix<double,embeddedDim,embeddedDim> d2d2_diff = d2d2; FieldMatrix<double,embeddedDim,embeddedDim> d2d2_diff = d2d2;
d2d2_diff -= d2d2_fd; d2d2_diff -= d2d2_fd;
if ( (d2d2_diff).infinity_norm() > 100*eps) { if ( (d2d2_diff).infinity_norm() > 200*eps) {
std::cout << className(a) << ": Analytical second derivative does not match fd approximation." << std::endl; std::cout << className(a) << ": Analytical second derivative does not match fd approximation." << std::endl;
std::cout << "d2d2 Analytical:" << std::endl << d2d2 << std::endl; std::cout << "d2d2 Analytical:" << std::endl << d2d2 << std::endl;
std::cout << "d2d2 FD :" << std::endl << d2d2_fd << std::endl; std::cout << "d2d2 FD :" << std::endl << d2d2_fd << std::endl;
...@@ -207,7 +208,8 @@ void testMixedDerivativesOfSquaredDistance(const TargetSpace& a, const TargetSpa ...@@ -207,7 +208,8 @@ void testMixedDerivativesOfSquaredDistance(const TargetSpace& a, const TargetSpa
FieldMatrix<double,embeddedDim,embeddedDim> d1d2_diff = d1d2; FieldMatrix<double,embeddedDim,embeddedDim> d1d2_diff = d1d2;
d1d2_diff -= d1d2_fd; d1d2_diff -= d1d2_fd;
if ( (d1d2_diff).infinity_norm() > 100*eps ) {
if ( (d1d2_diff).infinity_norm() > 200*eps ) {
std::cout << className(a) << ": Analytical mixed second derivative does not match fd approximation." << std::endl; std::cout << className(a) << ": Analytical mixed second derivative does not match fd approximation." << std::endl;
std::cout << "d1d2 Analytical:" << std::endl << d1d2 << std::endl; std::cout << "d1d2 Analytical:" << std::endl << d1d2 << std::endl;
std::cout << "d1d2 FD :" << std::endl << d1d2_fd << std::endl; std::cout << "d1d2 FD :" << std::endl << d1d2_fd << std::endl;
...@@ -245,8 +247,8 @@ void testDerivativeOfHessianOfSquaredDistance(const TargetSpace& a, const Target ...@@ -245,8 +247,8 @@ void testDerivativeOfHessianOfSquaredDistance(const TargetSpace& a, const Target
d2d2d2_fd[i] /= 2*eps; d2d2d2_fd[i] /= 2*eps;
} }
if ( (d2d2d2 - d2d2d2_fd).infinity_norm() > 100*eps) { if ( (d2d2d2 - d2d2d2_fd).infinity_norm() > 200*eps) {
std::cout << className(a) << ": Analytical third derivative does not match fd approximation." << std::endl; std::cout << className(a) << ": Analytical third derivative does not match fd approximation." << std::endl;
std::cout << "d2d2d2 Analytical:" << std::endl << d2d2d2 << std::endl; std::cout << "d2d2d2 Analytical:" << std::endl << d2d2d2 << std::endl;
std::cout << "d2d2d2 FD :" << std::endl << d2d2d2_fd << std::endl; std::cout << "d2d2d2 FD :" << std::endl << d2d2d2_fd << std::endl;
...@@ -284,8 +286,8 @@ void testMixedDerivativeOfHessianOfSquaredDistance(const TargetSpace& a, const T ...@@ -284,8 +286,8 @@ void testMixedDerivativeOfHessianOfSquaredDistance(const TargetSpace& a, const T
d1d2d2_fd[i] /= 2*eps; d1d2d2_fd[i] /= 2*eps;
} }
if ( (d1d2d2 - d1d2d2_fd).infinity_norm() > 100*eps ) { if ( (d1d2d2 - d1d2d2_fd).infinity_norm() > 200*eps ) {
std::cout << className(a) << ": Analytical mixed third derivative does not match fd approximation." << std::endl; std::cout << className(a) << ": Analytical mixed third derivative does not match fd approximation." << std::endl;
std::cout << "d1d2d2 Analytical:" << std::endl << d1d2d2 << std::endl; std::cout << "d1d2d2 Analytical:" << std::endl << d1d2d2 << std::endl;
std::cout << "d1d2d2 FD :" << std::endl << d1d2d2_fd << std::endl; std::cout << "d1d2d2 FD :" << std::endl << d1d2d2_fd << std::endl;
...@@ -360,7 +362,7 @@ void test() ...@@ -360,7 +362,7 @@ void test()
std::cout << "p: " << testPointPair[0] << ", q: " << testPointPair[1] << std::endl; std::cout << "p: " << testPointPair[0] << ", q: " << testPointPair[1] << std::endl;
std::cout << TargetSpace::exp(p, TargetSpace::log(p,q)) << std::endl; std::cout << TargetSpace::exp(p, TargetSpace::log(p,q)) << std::endl;
//testDerivativesOfSquaredDistance<TargetSpace>(testPoints[i], testPoints[j]); testDerivativesOfSquaredDistance<TargetSpace>(testPoints[i], testPoints[j]);
} }
...@@ -378,6 +380,9 @@ int main() try ...@@ -378,6 +380,9 @@ int main() try
test<UnitVector<double,3> >(); test<UnitVector<double,3> >();
test<UnitVector<double,4> >(); test<UnitVector<double,4> >();
test<Dune::GFE::ProductManifold<RealTuple<double,1>,Rotation<double,3>,UnitVector<double,2>>>();
test<Dune::GFE::ProductManifold<Rotation<double,3>,UnitVector<double,5>>>();
// test<Rotation<double,3> >(); // test<Rotation<double,3> >();
// //
// test<RigidBodyMotion<double,3> >(); // test<RigidBodyMotion<double,3> >();
......
test/valuefactory.hh
+ 27
− 0
View file @ 4a41903a
...@@ -6,6 +6,7 @@ ...@@ -6,6 +6,7 @@
#include <dune/gfe/unitvector.hh> #include <dune/gfe/unitvector.hh>
#include <dune/gfe/rotation.hh> #include <dune/gfe/rotation.hh>
#include <dune/gfe/rigidbodymotion.hh> #include <dune/gfe/rigidbodymotion.hh>
#include <dune/gfe/productmanifold.hh>
#include <dune/gfe/orthogonalmatrix.hh> #include <dune/gfe/orthogonalmatrix.hh>
#include <dune/gfe/hyperbolichalfspacepoint.hh> #include <dune/gfe/hyperbolichalfspacepoint.hh>
...@@ -328,4 +329,30 @@ public: ...@@ -328,4 +329,30 @@ public:
/** \brief A class that creates sets of values of various types, to be used in unit tests
*
* This is the specialization for ProducManifold<...>
*/
template <typename ...TargetSpaces>
class ValueFactory<Dune::GFE::ProductManifold<TargetSpaces...> >
{
using TargetSpace = Dune::GFE::ProductManifold<TargetSpaces...>;
public:
static void get(std::vector<TargetSpace >& values) {
std::vector<typename TargetSpace::CoordinateType > testPoints(10);
std::generate(testPoints.begin(), testPoints.end(), [](){
return Dune::GFE::randomFieldVector<typename TargetSpace::field_type,TargetSpace::CoordinateType::dimension>(0.9,1.1) ; });
values.resize(testPoints.size());
std::transform(testPoints.cbegin(),testPoints.cend(),values.begin(),[](const auto& testPoint){return TargetSpace(testPoint);});
}
};
#endif #endif
\ No newline at end of file
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment