diff --git a/include/geos/Makefile.am b/include/geos/Makefile.am index 363494cf5e..e34f56d5ee 100644 --- a/include/geos/Makefile.am +++ b/include/geos/Makefile.am @@ -8,6 +8,7 @@ SUBDIRS = \ index \ io \ linearref \ + math \ noding \ operation \ planargraph \ diff --git a/include/geos/math/DD.h b/include/geos/math/DD.h new file mode 100644 index 0000000000..09ccd30396 --- /dev/null +++ b/include/geos/math/DD.h @@ -0,0 +1,94 @@ +/********************************************************************** + * + * GEOS - Geometry Engine Open Source + * http://geos.osgeo.org + * + * Copyright (C) 2001-2002 Vivid Solutions Inc. + * Copyright (C) 2020 Crunchy Data + * + * This is free software; you can redistribute and/or modify it under + * the terms of the GNU Lesser General Public Licence as published + * by the Free Software Foundation. + * See the COPYING file for more information. + * + **********************************************************************/ + +#ifndef GEOS_MATH_TOPOLOGYEXCEPTION_H +#define GEOS_MATH_TOPOLOGYEXCEPTION_H + +#include + +namespace geos { +namespace math { // geos.math + +/** + * \class DD + * + * \brief + * Wrapper for DoubleDouble higher precision mathematics + * operations. + */ +class GEOS_DLL DD { + private: + double hi; + double lo; + + DD parse(std::string &str); + int magnitude(double x); + + public: + DD(double p_hi, double p_lo) : hi(p_hi), lo(p_lo) {}; + DD(double x) : hi(x), lo(0.0) {}; + DD(DD &dd) : hi(dd.hi), lo(dd.lo) {}; + + bool operator==(DD const &rhs) const + { + return hi == rhs.hi && lo == rhs.lo; + } + + bool operator!=(DD const &rhs) const + { + return hi != rhs.hi || lo != rhs.lo; + } + + bool operator<(DD const &rhs) const + { + return (hi < rhs.hi) || (hi == rhs.hi && lo < rhs.lo); + } + + bool operator<=(DD const &rhs) const + { + return (hi < rhs.hi) || (hi == rhs.hi && lo <= rhs.lo); + } + + bool operator>(DD const &rhs) const + { + return (hi > rhs.hi) || (hi == rhs.hi && lo > rhs.lo); + } + + bool operator>=(DD const &rhs) const + { + return (hi > rhs.hi) || (hi == rhs.hi && lo >= rhs.lo); + } + + DD operator+(DD const &rhs) const; + DD operator+(double rhs) const; + + + bool isNaN(); + bool isNegative(); + bool isPositive(); + bool isZero(); + double doubleValue(); + int intValue(); + + void selfAdd(DD const &d); + void selfAdd(double p_hi, double p_lo); + void selfAdd(double y); +}; + +} // namespace geos::math +} // namespace geos + + +#endif // GEOS_MATH_TOPOLOGYEXCEPTION_H diff --git a/include/geos/math/Makefile.am b/include/geos/math/Makefile.am new file mode 100644 index 0000000000..0ebaad4a63 --- /dev/null +++ b/include/geos/math/Makefile.am @@ -0,0 +1,11 @@ +# +# This file is part of project GEOS (http://trac.osgeo.org/geos/) +# +#SUBDIRS = + +#EXTRA_DIST = + +geosdir = $(includedir)/geos/math + +geos_HEADERS = \ + DD.h diff --git a/src/Makefile.am b/src/Makefile.am index 1bba8e87ef..25d6edd0d2 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -7,6 +7,7 @@ SUBDIRS = \ index \ io \ linearref \ + math \ noding \ operation \ planargraph \ diff --git a/src/math/DD.cpp b/src/math/DD.cpp new file mode 100644 index 0000000000..1fe2a33051 --- /dev/null +++ b/src/math/DD.cpp @@ -0,0 +1,136 @@ +#include +#include +#include +#include + +using namespace std; + +#include +#include + +namespace geos { +namespace math { // geos.util + + +/* private */ +DD +DD::parse(std::string &str) +{ + DD d(0.0, 0.0); + return d; +} + +/* private */ +int +DD::magnitude(double x) +{ + double xAbs = std::fabs(x); + double xLog10 = std::log(xAbs) / std::log(10); + int xMag = (int) std::floor(xLog10); + /** + * Since log computation is inexact, there may be an off-by-one error + * in the computed magnitude. + * Following tests that magnitude is correct, and adjusts it if not + */ + double xApprox = std::pow(10, xMag); + if (xApprox * 10 <= xAbs) + xMag += 1; + + return xMag; +} + +/* public */ +bool DD::isNaN() +{ + return std::isnan(hi); +} +/* public */ +bool DD::isNegative() +{ + return hi < 0.0 || (hi == 0.0 && lo < 0.0); +} +/* public */ +bool DD::isPositive() +{ + return hi > 0.0 || (hi == 0.0 && lo > 0.0); +} +/* public */ +bool DD::isZero() +{ + return hi == 0.0 && lo == 0.0; +} + +/* public */ +double DD::doubleValue() +{ + return hi + lo; +} + +/* public */ +int DD::intValue() +{ + return (int) hi; +} + +/* public */ +void DD::selfAdd(DD const &y) +{ + return selfAdd(y.hi, y.lo); +} + +/* public */ +void DD::selfAdd(double yhi, double ylo) +{ + double H, h, T, t, S, s, e, f; + S = hi + yhi; + T = lo + ylo; + e = S - hi; + f = T - lo; + s = S-e; + t = T-f; + s = (yhi-e)+(hi-s); + t = (ylo-f)+(lo-t); + e = s+T; H = S+e; h = e+(S-H); e = t+h; + + double zhi = H + e; + double zlo = e + (H - zhi); + hi = zhi; + lo = zlo; + return; +} + +/* public */ +void DD::selfAdd(double y) +{ + double H, h, S, s, e, f; + S = hi + y; + e = S - hi; + s = S - e; + s = (y - e) + (hi - s); + f = s + lo; + H = S + f; + h = f + (S - H); + hi = H + h; + lo = h + (H - hi); + return; +} + +/* public */ +DD DD::operator+(DD const &rhs) const +{ + DD lhs(hi, lo); + lhs.selfAdd(rhs); + return lhs; +} + +/* public */ +DD DD::operator+(double rhs) const +{ + DD lhs(hi, lo); + lhs.selfAdd(rhs); + return lhs; +} + + +} +} diff --git a/src/math/Makefile.am b/src/math/Makefile.am new file mode 100644 index 0000000000..c30e97eac5 --- /dev/null +++ b/src/math/Makefile.am @@ -0,0 +1,13 @@ +# +# This file is part of project GEOS (http://trac.osgeo.org/geos/) +# +SUBDIRS = + +noinst_LTLIBRARIES = libmath.la + +AM_CPPFLAGS = -I$(top_srcdir)/include + +libmath_la_SOURCES = \ + DD.cpp + +libmath_la_LIBADD =