From 80aafc23e235b9e9864abe8e108bc673f3d5fa5c Mon Sep 17 00:00:00 2001 From: Oliver Sander <sander@igpm.rwth-aachen.de> Date: Wed, 6 Apr 2011 07:21:48 +0000 Subject: [PATCH] more test points and some cleanup [[Imported from SVN: r7099]] --- test/localgeodesicfefunctiontest.cc | 121 ++++++++++++++++++---------- 1 file changed, 78 insertions(+), 43 deletions(-) diff --git a/test/localgeodesicfefunctiontest.cc b/test/localgeodesicfefunctiontest.cc index b2333bf5..08e8babd 100644 --- a/test/localgeodesicfefunctiontest.cc +++ b/test/localgeodesicfefunctiontest.cc @@ -15,6 +15,8 @@ // Domain dimension const int dim = 2; +const double eps = 1e-6; + using namespace Dune; void testDerivativeTangentiality(const RealTuple<1>& x, @@ -36,8 +38,8 @@ void testDerivativeTangentiality(const UnitVector<vectorDim>& x, for (int j=0; j<vectorDim; j++) sp += x.globalCoordinates()[j] * derivative[j][i]; - - std::cout << "Column: " << i << ", product: " << sp << std::endl; + if (std::fabs(sp) > 1e-8) + DUNE_THROW(Dune::Exception, "Derivative is not tangential: Column: " << i << ", product: " << sp); } @@ -84,8 +86,8 @@ void testPermutationInvariance(const std::vector<TargetSpace>& corners) TargetSpace v2 = f2.evaluate(l2); // Check that they are all equal - assert(TargetSpace::distance(v0,v1) < 1e-5); - assert(TargetSpace::distance(v0,v2) < 1e-5); + assert(TargetSpace::distance(v0,v1) < eps); + assert(TargetSpace::distance(v0,v2) < eps); } @@ -113,8 +115,14 @@ void testDerivative(const std::vector<TargetSpace>& corners) // evaluate fd approximation of derivative Dune::FieldMatrix<double, TargetSpace::EmbeddedTangentVector::size, dim> fdDerivative = f.evaluateDerivativeFD(quadPos); - std::cout << "Analytical: " << std::endl << derivative << std::endl; - std::cout << "FD: " << std::endl << fdDerivative << std::endl; + Dune::FieldMatrix<double, TargetSpace::EmbeddedTangentVector::size, dim> diff = derivative; + diff -= fdDerivative; + + if ( diff.infinity_norm() > 100*eps ) { + std::cout << className(corners[0]) << ": Analytical gradient does not match fd approximation." << std::endl; + std::cout << "Analytical: " << derivative << std::endl; + std::cout << "FD : " << fdDerivative << std::endl; + } testDerivativeTangentiality(f.evaluate(quadPos), derivative); @@ -135,54 +143,81 @@ void testRealTuples() testDerivative(corners); } -void testUnitVectors() +void testUnitVector2d() { - std::cout << " --- Testing UnitVector<3> ---" << std::endl; + std::cout << " --- Testing UnitVector<2> ---" << std::endl; - typedef UnitVector<3> TargetSpace; + typedef UnitVector<2> TargetSpace; + int nTestPoints = 10; + double testPoints[10][2] = {{1,0}, {0.5,0.5}, {0,1}, {-0.5,0.5}, {-1,0}, {-0.5,-0.5}, {0,-1}, {0.5,-0.5}, {0.1,1}, {1,.1}}; + assert(dim==2); + + // Set up elements of S^1 std::vector<TargetSpace> corners(dim+1); - // test some simplex - FieldVector<double,3> input; - input[0] = 1; input[1] = 0; input[2] = 0; - corners[0] = input; - input[0] = 0; input[1] = 1; input[2] = 0; - corners[1] = input; - input[0] = 0; input[1] = 0; input[2] = 1; - corners[2] = input; - - testPermutationInvariance(corners); - testDerivative(corners); - - // test the constant function, i.e., everything is mapped onto a single point - input[0] = 1; input[1] = 0; input[2] = 0; - corners[0] = input; - corners[1] = input; - corners[2] = input; - - testPermutationInvariance(corners); - testDerivative(corners); + for (int i=0; i<nTestPoints; i++) { + + for (int j=0; j<nTestPoints; j++) { + + for (int k=0; k<nTestPoints; k++) { + + Dune::array<double,2> w0 = {testPoints[i][0], testPoints[i][1]}; + corners[0] = UnitVector<2>(w0); + Dune::array<double,2> w1 = {testPoints[j][0], testPoints[j][1]}; + corners[1] = UnitVector<2>(w1); + Dune::array<double,2> w2 = {testPoints[k][0], testPoints[k][1]}; + corners[2] = UnitVector<2>(w2); + + if (UnitVector<2>::distance(corners[0],corners[1]) > M_PI*0.98 + or UnitVector<2>::distance(corners[1],corners[2]) > M_PI*0.98 + or UnitVector<2>::distance(corners[2],corners[0]) > M_PI*0.98) + continue; + + testPermutationInvariance(corners); + testDerivative(corners); + + } + } + } } -void testUnitVectors2() +void testUnitVector3d() { - std::cout << " --- Testing UnitVector<2> ---" << std::endl; + std::cout << " --- Testing UnitVector<3> ---" << std::endl; - typedef UnitVector<2> TargetSpace; + typedef UnitVector<3> TargetSpace; + int nTestPoints = 10; + double testPoints[10][3] = {{1,0,0}, {0,1,0}, {-0.838114,0.356751,-0.412667}, + {-0.490946,-0.306456,0.81551},{-0.944506,0.123687,-0.304319}, + {-0.6,0.1,-0.2},{0.45,0.12,0.517}, + {-0.1,0.3,-0.1},{-0.444506,0.123687,0.104319},{-0.7,-0.123687,-0.304319}}; + assert(dim==2); + + // Set up elements of S^1 std::vector<TargetSpace> corners(dim+1); - FieldVector<double,2> input; - input[0] = 1; input[1] = 0; - corners[0] = input; - input[0] = 1; input[1] = 0; - corners[1] = input; - input[0] = 0; input[1] = 1; - corners[2] = input; + for (int i=0; i<nTestPoints; i++) { + + for (int j=0; j<nTestPoints; j++) { + + for (int k=0; k<nTestPoints; k++) { + + Dune::array<double,3> w0 = {testPoints[i][0], testPoints[i][1], testPoints[i][2]}; + corners[0] = UnitVector<3>(w0); + Dune::array<double,3> w1 = {testPoints[j][0], testPoints[j][1], testPoints[j][2]}; + corners[1] = UnitVector<3>(w1); + Dune::array<double,3> w2 = {testPoints[k][0], testPoints[k][1], testPoints[k][2]}; + corners[2] = UnitVector<3>(w2); + + testPermutationInvariance(corners); + testDerivative(corners); + + } + } + } - testPermutationInvariance(corners); - testDerivative(corners); } void testRotations() @@ -215,7 +250,7 @@ int main() feenableexcept(FE_INVALID); //testRealTuples(); - testUnitVectors(); - testUnitVectors2(); + testUnitVector2d(); + testUnitVector3d(); //testRotations(); } -- GitLab