53 lines
1.7 KiB
Diff
53 lines
1.7 KiB
Diff
From 5d37981978494681a6e55224ac6e6f765eb59edb Mon Sep 17 00:00:00 2001
|
|
From: Roberto Di Remigio <roberto.diremigio@gmail.com>
|
|
Date: Mon, 23 Aug 2021 08:23:23 +0200
|
|
Subject: [PATCH] Always use a temporary of large enough size
|
|
|
|
Should fix #151
|
|
---
|
|
src/functionals/brx.cpp | 15 +++++++++++----
|
|
1 file changed, 11 insertions(+), 4 deletions(-)
|
|
|
|
diff --git a/src/functionals/brx.cpp b/src/functionals/brx.cpp
|
|
index 260a0ed..acaaafc 100644
|
|
--- a/src/functionals/brx.cpp
|
|
+++ b/src/functionals/brx.cpp
|
|
@@ -50,11 +50,15 @@ static double BR(double z) {
|
|
// Obtain the Taylor expansion of x(y), which is the
|
|
// inverse of BR_y. Use linear method for simplicity.
|
|
template <typename T, int Ndeg>
|
|
-void BR_taylor(const T & z0, taylor<T, 1, Ndeg> & t) {
|
|
- taylor<T, 1, Ndeg> f, d;
|
|
+taylor <T, 1, Ndeg> BR_taylor(const T & z0) {
|
|
+ static_assert(Ndeg >= 3);
|
|
+
|
|
+ taylor<T, 1, Ndeg> t;
|
|
t = 0;
|
|
t[0] = BR(z0);
|
|
t[1] = 1;
|
|
+
|
|
+ taylor<T, 1, Ndeg> f;
|
|
f = BR_z(t);
|
|
t[1] = 1 / f[1];
|
|
// Linear method, for quadratic see i.e. Brent & Kung ~197x
|
|
@@ -62,6 +66,8 @@ void BR_taylor(const T & z0, taylor<T, 1, Ndeg> & t) {
|
|
f = BR_z(t);
|
|
t[i] = -f[i] * t[1];
|
|
}
|
|
+
|
|
+ return t;
|
|
}
|
|
|
|
/* This is a fully differentiable solver for Eq.(21) in
|
|
@@ -70,8 +76,9 @@ void BR_taylor(const T & z0, taylor<T, 1, Ndeg> & t) {
|
|
*/
|
|
template <typename T, int Nvar>
|
|
static ctaylor<T, Nvar> BR(const ctaylor<T, Nvar> & t) {
|
|
- taylor<T, 1, Nvar> tmp;
|
|
- BR_taylor(t.c[0], tmp);
|
|
+ // temporary has dimension 3 at least. See: https://github.com/dftlibs/xcfun/issues/151
|
|
+ // in C++14 and later can use std::max(Nvar, 3) for the second template argument
|
|
+ auto tmp = BR_taylor<T, (Nvar >= 3) ? Nvar : 3>(t.c[0]);
|
|
|
|
ctaylor<T, Nvar> res = tmp[0];
|
|
for (int i = 1; i <= Nvar; i++)
|