From 0e34d2640ba99a97d680f9449904394e67d8e877 Mon Sep 17 00:00:00 2001
From: Oliver Sander <sander@igpm.rwth-aachen.de>
Date: Wed, 6 Apr 2011 09:19:15 +0000
Subject: [PATCH] Implement test for the derivative of function gradients wrt
coefficients.
This test itself is not tested yet and may still contain bugs.
Also, I out-commented the permutation-invariance test. I have to make
sure somehow that test points are close enough together.
[[Imported from SVN: r7100]]
---
test/localgeodesicfefunctiontest.cc | 68 ++++++++++++++++++++++++++++-
1 file changed, 66 insertions(+), 2 deletions(-)
diff --git a/test/localgeodesicfefunctiontest.cc b/test/localgeodesicfefunctiontest.cc
index 08e8babd..1470414e 100644
--- a/test/localgeodesicfefunctiontest.cc
+++ b/test/localgeodesicfefunctiontest.cc
@@ -129,6 +129,68 @@ void testDerivative(const std::vector<TargetSpace>& corners)
}
}
+template <class TargetSpace>
+void testDerivativeOfGradientWRTCoefficients(const std::vector<TargetSpace>& corners)
+{
+ // Make local fe function to be tested
+ LocalGeodesicFEFunction<2,double,TargetSpace> f(corners);
+
+ // A quadrature rule as a set of test points
+ int quadOrder = 3;
+
+ const Dune::QuadratureRule<double, dim>& quad
+ = Dune::QuadratureRules<double, dim>::rule(GeometryType(GeometryType::simplex,dim), quadOrder);
+
+ for (size_t pt=0; pt<quad.size(); pt++) {
+
+ const Dune::FieldVector<double,dim>& quadPos = quad[pt].position();
+
+ // loop over the coefficients
+ for (size_t i=0; i<corners.size(); i++) {
+
+ // evaluate actual derivative
+ Tensor3<double, TargetSpace::EmbeddedTangentVector::size, TargetSpace::EmbeddedTangentVector::size, dim> derivative;
+ f.evaluateDerivativeOfGradientWRTCoefficient(quadPos, i, derivative);
+
+ // evaluate fd approximation of derivative
+ Tensor3<double, TargetSpace::EmbeddedTangentVector::size, TargetSpace::EmbeddedTangentVector::size, dim> fdDerivative;
+
+ for (int j=0; j<TargetSpace::EmbeddedTangentVector::size; j++) {
+
+ std::vector<TargetSpace> cornersPlus = corners;
+ std::vector<TargetSpace> cornersMinus = corners;
+ FieldVector<double,dim> aPlus = corners[i].globalCoordinates();
+ FieldVector<double,dim> aMinus = corners[i].globalCoordinates();
+ aPlus[j] += eps;
+ aMinus[j] -= eps;
+ cornersPlus[i] = TargetSpace(aPlus);
+ cornersMinus[i] = TargetSpace(aMinus);
+ LocalGeodesicFEFunction<2,double,TargetSpace> fPlus(cornersPlus);
+ LocalGeodesicFEFunction<2,double,TargetSpace> fMinus(cornersMinus);
+
+ FieldMatrix<double,TargetSpace::EmbeddedTangentVector::size,dim> hPlus = fPlus.evaluateDerivative(quadPos);
+ FieldMatrix<double,TargetSpace::EmbeddedTangentVector::size,dim> hMinus = fMinus.evaluateDerivative(quadPos);
+
+ fdDerivative[j] = hPlus;
+ fdDerivative[j] -= hMinus;
+ fdDerivative[j] /= 2*eps;
+
+ }
+
+ if ( (derivative - fdDerivative).infinity_norm() > eps ) {
+ std::cout << className(corners[0]) << ": Analytical derivative of 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);
+
+ }
+
+ }
+}
+
+
void testRealTuples()
{
std::cout << " --- Testing RealTuple<1> ---" << std::endl;
@@ -174,8 +236,9 @@ void testUnitVector2d()
or UnitVector<2>::distance(corners[2],corners[0]) > M_PI*0.98)
continue;
- testPermutationInvariance(corners);
+ //testPermutationInvariance(corners);
testDerivative(corners);
+ testDerivativeOfGradientWRTCoefficients(corners);
}
}
@@ -211,8 +274,9 @@ void testUnitVector3d()
Dune::array<double,3> w2 = {testPoints[k][0], testPoints[k][1], testPoints[k][2]};
corners[2] = UnitVector<3>(w2);
- testPermutationInvariance(corners);
+ //testPermutationInvariance(corners);
testDerivative(corners);
+ testDerivativeOfGradientWRTCoefficients(corners);
}
}
--
GitLab