[Libreoffice-commits] core.git: sc/inc sc/source

[dante (via logerrit)]

 sc/inc/kahan.hxx                 |   40 ++++++++++++++++++++++++++++++++++++++-
 sc/source/core/tool/interpr3.cxx |   26 ++++++++++++++-----------
 2 files changed, 54 insertions(+), 12 deletions(-)

New commits:
commit 94977bb43d8dc91023bee7afa037f6319c36ccc3
Author:     dante <dante19031999 at gmail.com>
AuthorDate: Sat May 1 11:04:21 2021 +0200
Commit:     Mike Kaganski <mike.kaganski at collabora.com>
CommitDate: Tue May 4 20:34:11 2021 +0200

    tdf#137679 Use kahan summation for ScInterpreter::lcl_IterateInverse
    
    Change-Id: I02e108ac70ddd4ea8d8d97eb4f5fbc8996dd4bd1
    Reviewed-on: https://gerrit.libreoffice.org/c/core/+/114966
    Tested-by: Jenkins
    Reviewed-by: Mike Kaganski <mike.kaganski at collabora.com>

diff --git a/sc/inc/kahan.hxx b/sc/inc/kahan.hxx
index ffeeab0867df..23a29f6ee13f 100644
--- a/sc/inc/kahan.hxx
+++ b/sc/inc/kahan.hxx
@@ -53,6 +53,16 @@ public:
         add(fSum.m_fError);
     }
 
+    /**
+      * Substracts a value to the sum using Kahan summation.
+      * @param fSum
+      */
+    inline void subtract(const KahanSum& fSum)
+    {
+        add(-fSum.m_fSum);
+        add(-fSum.m_fError);
+    }
+
 public:
     constexpr KahanSum operator-() const
     {
@@ -75,7 +85,7 @@ public:
 
     inline void operator+=(double fSum) { add(fSum); }
 
-    inline void operator-=(const KahanSum& fSum) { add(-fSum); }
+    inline void operator-=(const KahanSum& fSum) { subtract(fSum); }
 
     inline void operator-=(double fSum) { add(-fSum); }
 
@@ -86,6 +96,13 @@ public:
         return fNSum;
     }
 
+    inline KahanSum operator+(const KahanSum& fSum) const
+    {
+        KahanSum fNSum(*this);
+        fNSum += fSum;
+        return fNSum;
+    }
+
     inline KahanSum operator-(double fSum) const
     {
         KahanSum fNSum(*this);
@@ -93,6 +110,13 @@ public:
         return fNSum;
     }
 
+    inline KahanSum operator-(const KahanSum& fSum) const
+    {
+        KahanSum fNSum(*this);
+        fNSum -= fSum;
+        return fNSum;
+    }
+
     /**
       * In some parts of the code of interpr_.cxx this may be used for
       * product instead of sum. This operator shall be used for that task.
@@ -109,6 +133,20 @@ public:
         m_fError /= fDivides;
     }
 
+    constexpr KahanSum operator*(double fTimes) const
+    {
+        KahanSum fSum(*this);
+        fSum *= fTimes;
+        return fSum;
+    }
+
+    constexpr KahanSum operator/(double fTimes) const
+    {
+        KahanSum fSum(*this);
+        fSum /= fTimes;
+        return fSum;
+    }
+
     constexpr bool operator<(const KahanSum& fSum) const { return get() < fSum.get(); }
 
     constexpr bool operator<(double fSum) const { return get() < fSum; }
diff --git a/sc/source/core/tool/interpr3.cxx b/sc/source/core/tool/interpr3.cxx
index 17fe63381362..c339a68dd80f 100644
--- a/sc/source/core/tool/interpr3.cxx
+++ b/sc/source/core/tool/interpr3.cxx
@@ -84,32 +84,36 @@ static double lcl_IterateInverse( const ScDistFunc& rFunction, double fAx, doubl
 
     //  find enclosing interval
 
+    KahanSum fkAx = fAx;
+    KahanSum fkBx = fBx;
     double fAy = rFunction.GetValue(fAx);
     double fBy = rFunction.GetValue(fBx);
-    double fTemp;
+    KahanSum fTemp;
     unsigned short nCount;
     for (nCount = 0; nCount < 1000 && !lcl_HasChangeOfSign(fAy,fBy); nCount++)
     {
         if (std::abs(fAy) <= std::abs(fBy))
         {
-            fTemp = fAx;
-            fAx += 2.0 * (fAx - fBx);
-            if (fAx < 0.0)
-                fAx = 0.0;
-            fBx = fTemp;
+            fTemp = fkAx;
+            fkAx += (fkAx - fkBx) * 2.0;
+            if (fkAx < 0.0)
+                fkAx = 0.0;
+            fkBx = fTemp;
             fBy = fAy;
-            fAy = rFunction.GetValue(fAx);
+            fAy = rFunction.GetValue(fkAx.get());
         }
         else
         {
-            fTemp = fBx;
-            fBx += 2.0 * (fBx - fAx);
-            fAx = fTemp;
+            fTemp = fkBx;
+            fkBx += (fkBx - fkAx) * 2.0;
+            fkAx = fTemp;
             fAy = fBy;
-            fBy = rFunction.GetValue(fBx);
+            fBy = rFunction.GetValue(fkBx.get());
         }
     }
 
+    fAx = fkAx.get();
+    fBx = fkBx.get();
     if (fAy == 0.0)
         return fAx;
     if (fBy == 0.0)