diff --git a/dune/tectonic/circularconvexfunction.hh b/dune/tectonic/circularconvexfunction.hh index 8e64cc6d2c737a69c3e9f2aa139fa9a2d4955cd6..c7fdc54736498f71e40e0882fecf29e3fc33c202 100644 --- a/dune/tectonic/circularconvexfunction.hh +++ b/dune/tectonic/circularconvexfunction.hh @@ -13,7 +13,13 @@ template <class NonlinearityType> class CircularConvexFunction { CircularConvexFunction(MatrixType const &A, VectorType const &b, NonlinearityType const &phi, VectorType const &x, VectorType const &dir) - : A(A), b(b), phi(phi), x(x), dir(dir) {} + : A(A), + b(b), + phi(phi), + x(x), + dir(dir), + xnorm(x.two_norm()), + dnorm(dir.two_norm()) {} double quadraticPart() const { return 0; } @@ -43,7 +49,7 @@ template <class NonlinearityType> class CircularConvexFunction { void cartesian(double m, VectorType &y) const { y = 0; y.axpy(std::cos(m), x); - y.axpy(std::sin(m), dir); + y.axpy(std::sin(m) * xnorm / dnorm, dir); } private: @@ -53,6 +59,9 @@ template <class NonlinearityType> class CircularConvexFunction { VectorType const &x; VectorType const &dir; + double const dnorm; + double const xnorm; + /* If x and d were normalised, we would have cartesian(a) = cos(a) * x + sin(a) * d and @@ -63,8 +72,8 @@ template <class NonlinearityType> class CircularConvexFunction { */ void tangent(double m, VectorType &y) const { y = 0; - y.axpy(-std::sin(m) * dir.two_norm2(), x); - y.axpy(std::cos(m) * x.two_norm2(), dir); + y.axpy(-std::sin(m) * dnorm, x); + y.axpy(std::cos(m) * xnorm, dir); } }; } diff --git a/src/test-gradient-sample-steep.cc b/src/test-gradient-sample-steep.cc index 42b8a54ba3c1bf3bc410b18e72b8b8b67eabdd0d..7cca14ab09719a71111499cd9a16ed091efca5ef 100644 --- a/src/test-gradient-sample-steep.cc +++ b/src/test-gradient-sample-steep.cc @@ -36,18 +36,18 @@ int main() { start[0] = 0; start[1] = 1; - double const ret1 = functionTester(J, start, 5); + double const ret1 = functionTester(J, start, 3); // Something random start[0] = 279; start[1] = -96; - double const ret2 = functionTester(J, start, 5); - assert(std::abs(ret1 - ret2) < 1e-7); + double const ret2 = functionTester(J, start, 3); + assert(std::abs(ret1 - ret2) < 1e-5); start[0] = 0; start[1] = 0; double const ret3 = functionTester(J, start, 1); - assert(std::abs(ret1 - ret3) < 1e-7); + assert(std::abs(ret1 - ret3) < 1e-5); } diff --git a/src/test-gradient-sample-steep2.cc b/src/test-gradient-sample-steep2.cc index 894cb3abb5da3dac5864b3c9ecf35789ec1d63ee..d5a40ad4f01b00914b9e035d5f69812458f1d6f2 100644 --- a/src/test-gradient-sample-steep2.cc +++ b/src/test-gradient-sample-steep2.cc @@ -42,12 +42,12 @@ int main() { start[0] = 279; start[1] = -96; - double const ret2 = functionTester(J, start, 5); - assert(std::abs(ret1 - ret2) < 1e-7); + double const ret2 = functionTester(J, start, 3); + assert(std::abs(ret1 - ret2) < 1e-5); start[0] = 0; start[1] = 0; - double const ret3 = functionTester(J, start, 2); - assert(std::abs(ret1 - ret3) < 1e-7); + double const ret3 = functionTester(J, start, 1); + assert(std::abs(ret1 - ret3) < 1e-5); } diff --git a/src/test-gradient-sample-verysteep.cc b/src/test-gradient-sample-verysteep.cc index 2552f6ecbd75c213c13dda4d70c422ebc71ef32f..0da5ad9f37165dd432669d0895458f9801e2fe6d 100644 --- a/src/test-gradient-sample-verysteep.cc +++ b/src/test-gradient-sample-verysteep.cc @@ -42,12 +42,12 @@ int main() { start[0] = 279; start[1] = -96; - double const ret2 = functionTester(J, start, 5); - assert(std::abs(ret1 - ret2) < 1e-7); + double const ret2 = functionTester(J, start, 4); + assert(std::abs(ret1 - ret2) < 1e-8); start[0] = 0; start[1] = 0; - double const ret3 = functionTester(J, start, 2); - assert(std::abs(ret1 - ret3) < 1e-7); + double const ret3 = functionTester(J, start, 1); + assert(std::abs(ret1 - ret3) < 1e-5); }