From e0f98c8386cde1ecaccb5a2c47c719ca25b97040 Mon Sep 17 00:00:00 2001
From: Oliver Sander <sander@igpm.rwth-aachen.de>
Date: Thu, 12 Jan 2012 17:37:29 +0000
Subject: [PATCH] implement fd approximation properly using the exponential map
[[Imported from SVN: r8378]]
---
test/targetspacetest.cc | 24 +++++++++++++++++-------
1 file changed, 17 insertions(+), 7 deletions(-)
diff --git a/test/targetspacetest.cc b/test/targetspacetest.cc
index a69d852f..cbea6cb5 100644
--- a/test/targetspacetest.cc
+++ b/test/targetspacetest.cc
@@ -121,16 +121,26 @@ void testDerivativeOfSquaredDistance(const TargetSpace& a, const TargetSpace& b)
typename TargetSpace::EmbeddedTangentVector d2 = TargetSpace::derivativeOfDistanceSquaredWRTSecondArgument(a, b);
// finite-difference approximation
- typename TargetSpace::EmbeddedTangentVector d2_fd;
- for (size_t i=0; i<embeddedDim; i++) {
- FieldVector<double,embeddedDim> bPlus = b.globalCoordinates();
- FieldVector<double,embeddedDim> bMinus = b.globalCoordinates();
- bPlus[i] += eps;
- bMinus[i] -= eps;
+ Dune::FieldMatrix<double,TargetSpace::TangentVector::dimension,embeddedDim> B = b.orthonormalFrame();
+
+ typename TargetSpace::TangentVector d2_fd;
+ for (size_t i=0; i<TargetSpace::TangentVector::dimension; i++) {
+
+ typename TargetSpace::EmbeddedTangentVector fwVariation = B[i];
+ typename TargetSpace::EmbeddedTangentVector bwVariation = B[i];
+ fwVariation *= eps;
+ bwVariation *= -eps;
+ TargetSpace bPlus = TargetSpace::exp(b,fwVariation);
+ TargetSpace bMinus = TargetSpace::exp(b,bwVariation);
+
d2_fd[i] = (energy(a,bPlus) - energy(a,bMinus)) / (2*eps);
}
+
+ // transform into embedded coordinates
+ typename TargetSpace::EmbeddedTangentVector d2_fd_embedded;
+ B.mtv(d2_fd,d2_fd_embedded);
- if ( (d2 - d2_fd).infinity_norm() > 100*eps ) {
+ if ( (d2 - d2_fd_embedded).infinity_norm() > 100*eps ) {
std::cout << className(a) << ": Analytical gradient does not match fd approximation." << std::endl;
std::cout << "d2 Analytical: " << d2 << std::endl;
std::cout << "d2 FD : " << d2_fd << std::endl;
--
GitLab