Skip to content

Commit

Permalink
V3.78: adding a double-double arithmetic type (#426)
Browse files Browse the repository at this point in the history
* Setting SEMVER to v3.78

* Setting SEMVER to v3.78

* setting default devcontainer to clang16

* header updates for elastic integer

* header update for c_api

* header update for elastic float

* header update for elastic posit

* header updates for elastic unum

* header updates of mixed-precision applications

* header updates

* adding mathematics section to accuracy applications

* std::numbers inspired numerical constants, sans the C++20 concepts

* AppleClang compilation fix

* header changes for applications of accuracy

* header update and code hygiene for precision applications

* header updates for applications generating mathematical constants

* creating a mixed-precision midpoint/lerp skeleton test to lerp test complex ranges

* code hygiene on precision applications

* header updates and code hygiene for precision applications

* header updates and code hygiene for rest of application projects

* residual compensation algorithm using the quire

* code hygiene for dfloat skeleton

* code hygiene for dfloat skeleton

* code hygiene for bfloat16 skeleton

* initial skeleton for doubledouble arithmetic type

* gcc/clang compilation fix for dfloat template

* adding a doubledouble skeleton

* adding missing quad-double directory and build file

* adding maxpos calculation to doubledouble

* WIP: bootstrapping doubledouble algorithms

* WIP: experiments how to verify doubledouble arithmetic correctness

* compilation and test fixes of doubledouble regression tests

* disabling doubledouble regression tests until coverage

* WIP: bug fix in quad-double multiply operator

* WIP: implementing numerics and classification functions for doubledouble

* WIP: decimal string printing and parsing of doubledouble

* replacing MSVC specific _fpclass API with standard C++ for doubledouble fpclassify

* WIP: skeletons for doubledouble math function regression tests

* bug fix for compilation guard logic for checkNaN/checkInf constexpr functions

* removing BIT_CAST_SUPPORT and standardize on bit_cast.hpp auto configuration

* removing constexpr from bfloat[8|16].assign() as clang doesn't have constexpr std::string

* second attempt to get the BIT_CAST_CONSTEXPR compiler guards correct across the three compilers

* third attempt to get BIT_CAST and LONG_DOUBLE support working across all compilers

* defeaturing std::bit_cast<> from processing long doubles

* removing std::bit_cast<blob, long double>

* adding doubledouble constant generation example

* streamlining doubledouble api tests

* adding a skeleton for experimentation with doubledouble

* changing to_binary of a doubledouble to reflect the hidden bit of the second segment

* adding performance measurements for dd vs native and quad emulation

* marking all variables and intermediates in the Kahan sums and products volatile

* compilation fix for gcc

* experiment with infinite sequence constants

* updating build containers to contain gdb

* code hygiene on volatile qualifier

* implementing setbit on double-double

* WIP: std::stream functionality for double-double

* new utility to report on language features configured

* removing <format> as gcc doesn't support it yet

* quick readme update

* making <stdfloat> conditional as gcc and clang do not support it yet

* code hygiene for double-double

* completing arithmetic regression test skeletons for double-double

* removing unsafe sscanf from parse function of double-double

* WIP: code streamlining and documentation

* code hygiene

* adding pow, exp, and sqrt algorithms for double-double courtesy of Scibuilders

* adding double-double constants of note and a generator

* generating exp table

* code hygiene

* adding shims for the rest of the math library for double-double

* compilation fix for gcc and clang
  • Loading branch information
Ravenwater authored Aug 7, 2024
1 parent 25fa75d commit 15ce7cd
Show file tree
Hide file tree
Showing 264 changed files with 9,660 additions and 1,082 deletions.
2 changes: 1 addition & 1 deletion .devcontainer/devcontainer.json
Original file line number Diff line number Diff line change
@@ -1,3 +1,3 @@
{
"image": "stillwater/builders:clang14builder"
"image": "stillwater/builders:clang16builder"
}
2 changes: 1 addition & 1 deletion .github/workflows/cmake.yml
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@ name: CMake

on:
push:
branches: [ v3.77, dev, main ]
branches: [ v3.78, dev, main ]
pull_request:
branches: [ main ]

Expand Down
19 changes: 17 additions & 2 deletions CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,7 @@ if(NOT DEFINED UNIVERSAL_VERSION_MAJOR)
set(UNIVERSAL_VERSION_MAJOR 3)
endif()
if(NOT DEFINED UNIVERSAL_VERSION_MINOR)
set(UNIVERSAL_VERSION_MINOR 77)
set(UNIVERSAL_VERSION_MINOR 78)
endif()
if(NOT DEFINED UNIVERSAL_VERSION_PATCH)
set(UNIVERSAL_VERSION_PATCH 1)
Expand Down Expand Up @@ -130,6 +130,8 @@ option(BUILD_NUMBER_FIXPNTS "Set to ON to build static fixed-point
option(BUILD_NUMBER_BFLOATS "Set to ON to build static bfloat tests" OFF)
option(BUILD_NUMBER_CFLOATS "Set to ON to build static cfloat tests" OFF)
option(BUILD_NUMBER_DFLOATS "Set to ON to build static dfloat tests" OFF)
option(BUILD_NUMBER_DDS "Set to ON to build static double-double tests" OFF)
option(BUILD_NUMBER_QDS "Set to ON to build static quad-double tests" OFF)
option(BUILD_NUMBER_AREALS "Set to ON to build static areal tests" OFF)
option(BUILD_NUMBER_UNUM1S "Set to ON to build static unum type 1 tests" OFF)
option(BUILD_NUMBER_UNUM2S "Set to ON to build static unum type 2 tests" OFF)
Expand Down Expand Up @@ -449,7 +451,7 @@ macro (compile_all testing prefix folder)
get_filename_component (test ${source} NAME_WE)
string(REPLACE " " ";" new_source ${source})
set(test_name ${prefix}_${test})
message(STATUS "Add test ${test_name} from source ${new_source}.")
#message(STATUS "Add test ${test_name} from source ${new_source}.")
add_executable (${test_name} ${new_source})

#add_custom_target(valid SOURCES ${SOURCES})
Expand Down Expand Up @@ -660,6 +662,8 @@ if(BUILD_NUMBER_STATICS)
set(BUILD_NUMBER_BFLOATS ON)
set(BUILD_NUMBER_CFLOATS ON)
set(BUILD_NUMBER_DFLOATS ON)
set(BUILD_NUMBER_DDS ON)
set(BUILD_NUMBER_QDS ON)
set(BUILD_NUMBER_AREALS ON)
set(BUILD_NUMBER_UNUM1S ON)
set(BUILD_NUMBER_UNUM2S ON)
Expand Down Expand Up @@ -822,6 +826,16 @@ if(BUILD_NUMBER_DFLOATS)
add_subdirectory("static/dfloat")
endif(BUILD_NUMBER_DFLOATS)

# double-double floats
if(BUILD_NUMBER_DDS)
add_subdirectory("static/dd")
endif(BUILD_NUMBER_DDS)

# quad-double floats
if(BUILD_NUMBER_QDS)
add_subdirectory("static/qd")
endif(BUILD_NUMBER_QDS)

# conversion tests suites
if(BUILD_NUMBER_CONVERSIONS)
add_subdirectory("static/conversions")
Expand Down Expand Up @@ -925,6 +939,7 @@ add_subdirectory("applications/accuracy/optimization")
add_subdirectory("applications/accuracy/pde")
add_subdirectory("applications/accuracy/roots")
add_subdirectory("applications/accuracy/science")
add_subdirectory("applications/accuracy/mathematics")

add_subdirectory("applications/approximation/taylor")
add_subdirectory("applications/approximation/chebyshev")
Expand Down
1 change: 0 additions & 1 deletion Dockerfile
Original file line number Diff line number Diff line change
Expand Up @@ -27,7 +27,6 @@ RUN cmake -DBUILD_ALL=ON .. && make

# RELEASE stage
#FROM alpine:latest as release # hitting a segfault during startup of some playground programs
#FROM debian:buster-slim as release
FROM ubuntu:latest as release
LABEL Theodore Omtzigt

Expand Down
2 changes: 1 addition & 1 deletion README.md
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
# Universal: a header-only C++ template library for universal number arithmetic
# Universal: a header-only C++ template library of custom arithmetic plug-in types


| **System** | **Status** | **More information** |
Expand Down
3 changes: 2 additions & 1 deletion applications/accuracy/engineering/chem_equilibrium.cpp
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
// chem_equilibrium.cpp : example of calculating the chemical balance of a solution
//
// Copyright (C) 2017-2021 Stillwater Supercomputing, Inc.
// Copyright (C) 2017 Stillwater Supercomputing, Inc.
// SPDX-License-Identifier: MIT
//
// This file is part of the universal numbers project, which is released under an MIT Open Source license.
#include <iostream>
Expand Down
3 changes: 2 additions & 1 deletion applications/accuracy/engineering/water.cpp
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
// water.cpp: example program showing water chemical equilibrium calculation sensitivity
//
// Copyright (C) 2017-2021 Stillwater Supercomputing, Inc.
// Copyright (C) 2017 Stillwater Supercomputing, Inc.
// SPDX-License-Identifier: MIT
//
// This file is part of the universal numbers project, which is released under an MIT Open Source license.
//
Expand Down
3 changes: 3 additions & 0 deletions applications/accuracy/mathematics/CMakeLists.txt
Original file line number Diff line number Diff line change
@@ -0,0 +1,3 @@
file (GLOB SOURCES "./*.cpp")

compile_all("true" "math" "Applications/Accuracy/Mathematics" "${SOURCES}")
47 changes: 47 additions & 0 deletions applications/accuracy/mathematics/README.md
Original file line number Diff line number Diff line change
@@ -0,0 +1,47 @@
# Computational Mathematics Examples

This directory contains computational Mathematics examples that demonstrate and quantify the benefits
of using custom arithmetic to solve difficult numerical calculations using approximated real numbers.

## Accuracy

Posits control their precision by the difference between `nbits`and `es`. The number of fraction bits in the posit
is given by the equation:

`fraction_bits = nbits - sign_bit - regime_bits - exponent_bits`

Around `1` and `-1`, the regime field is encoded in its minimum form of 2 bits.
At the extreme regimes, the posit encoding will have no fraction bits, and effectively be an exponential.

## Dynamic Range

Posits control their dynamic range by the number of exponent bits. The number of exponent bits directly
controls the value of useed, which is the exponential shift, 2^2^es, of the regime.

```text
IEEE float configurations from numeric_limits<Ty>
float minexp scale -125 maxexp scale 128 minimum 1.17549e-38 maximum 3.40282e+38
double minexp scale -1021 maxexp scale 1024 minimum 2.22507e-308 maximum 1.79769e+308
long double minexp scale -1021 maxexp scale 1024 minimum 2.22507e-308 maximum 1.79769e+308
```

## Reproducibility

The crown jewel of posit arithmetic is the ability to restore associative and distributive laws of algebra.

With IEEE floating point, associativity, that is

```eq
y = (a + b) + c -> y = a + (b + c)
```

does not hold. The reason is that with IEEE floating point arithmetic, each operation immediately rounds
the result to the available precision. This caused two different rounding paths for the two equations.

Similarly, distribution, that is

```eq
y = (a + b) * c -> y = a*c + b*c
```

also does not hold, for the same reason.
223 changes: 223 additions & 0 deletions applications/accuracy/mathematics/numeric_constants.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,223 @@
// numeric_constants.cpp: experiments with mixed-precision representations of important numerical constants
//
// Copyright (C) 2017 Stillwater Supercomputing, Inc.
// SPDX-License-Identifier: MIT
//
// This file is part of the UNIVERSAL project, which is released under an MIT Open Source license.
#include <universal/utility/directives.hpp>
#include <universal/traits/arithmetic_traits.hpp>
#include <universal/number/cfloat/cfloat.hpp>
#include <universal/number/posit/posit.hpp>
#include <universal/number/integer/integer.hpp>
#include <numbers> // since C++20

namespace sw {
namespace universal {

/*
template <class _Ty>
struct _Invalid {
static_assert(_Always_false<_Ty>, "A program that instantiates a primary template of a mathematical constant "
"variable template is ill-formed. (N4950 [math.constants]/3)");
};
template <class _Ty>
inline constexpr _Ty e_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty log2e_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty log10e_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty pi_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty inv_pi_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty inv_sqrtpi_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty ln2_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty ln10_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty sqrt2_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty sqrt3_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty inv_sqrt3_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty egamma_v = _Invalid<_Ty>{};
template <class _Ty>
inline constexpr _Ty phi_v = _Invalid<_Ty>{};
*/

template <typename _Ty>
inline CONSTEXPRESSION _Ty e_v = static_cast<_Ty>(2.718281828459045);
template <typename _Ty>
inline CONSTEXPRESSION _Ty log2e_v = static_cast<_Ty>(1.4426950408889634);
template <typename _Ty>
inline CONSTEXPRESSION _Ty log10e_v = static_cast<_Ty>(0.4342944819032518);
template <typename _Ty>
inline CONSTEXPRESSION _Ty pi_v = static_cast<_Ty>(3.141592653589793);
template <typename _Ty>
inline CONSTEXPRESSION _Ty inv_pi_v = static_cast<_Ty>(0.3183098861837907);
template <typename _Ty>
inline CONSTEXPRESSION _Ty inv_sqrtpi_v = static_cast<_Ty>(0.5641895835477563);
template <typename _Ty>
inline CONSTEXPRESSION _Ty ln2_v = static_cast<_Ty>(0.6931471805599453);
template <typename _Ty>
inline CONSTEXPRESSION _Ty ln10_v = static_cast<_Ty>(2.302585092994046);
template <typename _Ty>
inline CONSTEXPRESSION _Ty sqrt2_v = static_cast<_Ty>(1.4142135623730951);
template <typename _Ty>
inline CONSTEXPRESSION _Ty sqrt3_v = static_cast<_Ty>(1.7320508075688772);
template <typename _Ty>
inline CONSTEXPRESSION _Ty inv_sqrt3_v = static_cast<_Ty>(0.5773502691896257);
template <typename _Ty>
inline CONSTEXPRESSION _Ty egamma_v = static_cast<_Ty>(0.5772156649015329);
template <typename _Ty>
inline CONSTEXPRESSION _Ty phi_v = static_cast<_Ty>(1.618033988749895);

inline CONSTEXPRESSION half e_h = e_v<half>;
inline CONSTEXPRESSION half log2e_h = log2e_v<half>;
inline CONSTEXPRESSION half log10e_h = log10e_v<half>;
inline CONSTEXPRESSION half pi_h = pi_v<half>;
inline CONSTEXPRESSION half inv_pi_h = inv_pi_v<half>;
inline CONSTEXPRESSION half inv_sqrtpi_h = inv_sqrtpi_v<half>;
inline CONSTEXPRESSION half ln2_h = ln2_v<half>;
inline CONSTEXPRESSION half ln10_h = ln10_v<half>;
inline CONSTEXPRESSION half sqrt2_h = sqrt2_v<half>;
inline CONSTEXPRESSION half sqrt3_h = sqrt3_v<half>;
inline CONSTEXPRESSION half inv_sqrt3_h = inv_sqrt3_v<half>;
inline CONSTEXPRESSION half egamma_h = egamma_v<half>;
inline CONSTEXPRESSION half phi_h = phi_v<half>;

inline CONSTEXPRESSION single e_f = e_v<single>;
inline CONSTEXPRESSION single log2e_f = log2e_v<single>;
inline CONSTEXPRESSION single log10e_f = log10e_v<single>;
inline CONSTEXPRESSION single pi_f = pi_v<single>;
inline CONSTEXPRESSION single inv_pi_f = inv_pi_v<single>;
inline CONSTEXPRESSION single inv_sqrtpi_f = inv_sqrtpi_v<single>;
inline CONSTEXPRESSION single ln2_f = ln2_v<single>;
inline CONSTEXPRESSION single ln10_f = ln10_v<single>;
inline CONSTEXPRESSION single sqrt2_f = sqrt2_v<single>;
inline CONSTEXPRESSION single sqrt3_f = sqrt3_v<single>;
inline CONSTEXPRESSION single inv_sqrt3_f = inv_sqrt3_v<single>;
inline CONSTEXPRESSION single egamma_f = egamma_v<single>;
inline CONSTEXPRESSION single phi_f = phi_v<single>;

inline CONSTEXPRESSION double e_d = std::numbers::e_v<double>;
inline CONSTEXPRESSION double log2e_d = std::numbers::log2e_v<double>;
inline CONSTEXPRESSION double log10e_d = std::numbers::log10e_v<double>;
inline CONSTEXPRESSION double pi_d = std::numbers::pi_v<double>;
inline CONSTEXPRESSION double inv_pi_d = std::numbers::inv_pi_v<double>;
inline CONSTEXPRESSION double inv_sqrtpi_d = std::numbers::inv_sqrtpi_v<double>;
inline CONSTEXPRESSION double ln2_d = std::numbers::ln2_v<double>;
inline CONSTEXPRESSION double ln10_d = std::numbers::ln10_v<double>;
inline CONSTEXPRESSION double sqrt2_d = std::numbers::sqrt2_v<double>;
inline CONSTEXPRESSION double sqrt3_d = std::numbers::sqrt3_v<double>;
inline CONSTEXPRESSION double inv_sqrt3_d = std::numbers::inv_sqrt3_v<double>;
inline CONSTEXPRESSION double egamma_d = std::numbers::egamma_v<double>;
inline CONSTEXPRESSION double phi_d = std::numbers::phi_v<double>;
}
}

template<typename _Ty>
void Compare(const std::string& tag, _Ty c, double ref) {
using namespace sw::universal;

std::cout << tag << " : ";
auto old = std::cout.precision();
std::cout << std::setprecision(std::numeric_limits<_Ty>::digits10);
std::cout << to_binary(c) << " : " << std::left << std::scientific << std::setw(std::numeric_limits<_Ty>::digits10+2) << std::setfill('0') << c << " : ";
double error = ref - double(c);
std::cout << std::left << std::scientific << std::setw(std::numeric_limits<_Ty>::digits10 + 3) << std::setfill('0') << std::abs(error) << '\n';
std::cout << std::setprecision(old);
}

void CompareHalfPrecisionConstants() {
using namespace sw::universal;
using namespace std::numbers;

std::cout << "constant : binary : value : error\n";
Compare("e_h ", e_h, e);
Compare("log2e_h ", log2e_h, log2e);
Compare("log10e_h ", log10e_h, log10e);
Compare("pi_h ", pi_h, pi);
Compare("inv_pi_h ", inv_pi_h, inv_pi);
Compare("inv_sqrtpi_h", inv_sqrtpi_h, inv_sqrtpi);
Compare("ln2_h ", ln2_h, ln2);
Compare("ln10_h ", ln10_h, ln10);
Compare("sqrt2_h ", sqrt2_h, sqrt2);
Compare("inv_sqrt3_h ", inv_sqrt3_h, inv_sqrt3);
Compare("egamma_h ", egamma_h, egamma);
Compare("phi_h ", phi_h, phi);
}

void CompareSinglePrecisionConstants() {
using namespace sw::universal;
using namespace std::numbers;

std::cout << "constant : binary : value : error\n";
Compare("e_f ", e_f, e);
Compare("log2e_f ", log2e_f, log2e);
Compare("log10e_f ", log10e_f, log10e);
Compare("pi_f ", pi_f, pi);
Compare("inv_pi_f ", inv_pi_f, inv_pi);
Compare("inv_sqrtpi_f", inv_sqrtpi_f, inv_sqrtpi);
Compare("ln2_f ", ln2_f, ln2);
Compare("ln10_f ", ln10_f, ln10);
Compare("sqrt2_f ", sqrt2_f, sqrt2);
Compare("inv_sqrt3_f ", inv_sqrt3_f, inv_sqrt3);
Compare("egamma_f ", egamma_f, egamma);
Compare("phi_f ", phi_f, phi);
}

void CompareDoublePrecisionConstants() {
using namespace sw::universal;
using namespace std::numbers;

std::cout << "constant : binary : value : error\n";
Compare("e_d ", e_d, e);
Compare("log2e_d ", log2e_d, log2e);
Compare("log10e_d ", log10e_d, log10e);
Compare("pi_d ", pi_d, pi);
Compare("inv_pi_d ", inv_pi_d, inv_pi);
Compare("inv_sqrtpi_d", inv_sqrtpi_d, inv_sqrtpi);
Compare("ln2_d ", ln2_d, ln2);
Compare("ln10_d ", ln10_d, ln10);
Compare("sqrt2_d ", sqrt2_d, sqrt2);
Compare("inv_sqrt3_d ", inv_sqrt3_d, inv_sqrt3);
Compare("egamma_d ", egamma_d, egamma);
Compare("phi_d ", phi_d, phi);
}

int main()
try {
using namespace sw::universal;

CompareHalfPrecisionConstants();
CompareSinglePrecisionConstants();
CompareDoublePrecisionConstants();

return EXIT_SUCCESS;
}
catch (char const* msg) {
std::cerr << "Caught ad-hoc exception: " << msg << std::endl;
return EXIT_FAILURE;
}
catch (const sw::universal::universal_arithmetic_exception& err) {
std::cerr << "Caught unexpected universal arithmetic exception: " << err.what() << std::endl;
return EXIT_FAILURE;
}
catch (const sw::universal::universal_internal_exception& err) {
std::cerr << "Caught unexpected universal internal exception: " << err.what() << std::endl;
return EXIT_FAILURE;
}
catch (std::runtime_error& err) {
std::cerr << "Caught unexpected runtime error: " << err.what() << std::endl;
return EXIT_FAILURE;
}
catch (...) {
std::cerr << "Caught unknown exception" << std::endl;
return EXIT_FAILURE;
}

3 changes: 2 additions & 1 deletion applications/accuracy/ode/convergence.cpp
Original file line number Diff line number Diff line change
@@ -1,6 +1,7 @@
// convergence.cpp: convergence analysis of ODE solvers
//
// Copyright (C) 2017-2021 Stillwater Supercomputing, Inc.
// Copyright (C) 2017 Stillwater Supercomputing, Inc.
// SPDX-License-Identifier: MIT
// Author: Jacob Todd [email protected]
//
// This file is part of the universal numbers project, which is released under an MIT Open Source license
Expand Down
Loading

0 comments on commit 15ce7cd

Please sign in to comment.