diff options
author | couriersud <couriersud@gmx.org> | 2019-11-03 01:33:48 +0100 |
---|---|---|
committer | couriersud <couriersud@gmx.org> | 2019-11-03 01:33:48 +0100 |
commit | cabfd7188ae3b4442f6d5b6346709ba3dc3be7f9 (patch) | |
tree | f937bc159c200133ff162bfa75cc6a846f652bdc /src | |
parent | 4c88923f463487236035bb16c6c32bdc4edafdb0 (diff) |
netlist: fix regression, preliminary __float128 support. (nw)
__float128 is a gnu extension delivering true 128 bit floating point
support. Currently not supported by clang. In addition, the quadmath
library needs to be linked. For the time being therefore disabled.
Diffstat (limited to 'src')
-rw-r--r-- | src/lib/netlist/nl_config.h | 47 | ||||
-rw-r--r-- | src/lib/netlist/plib/gmres.h | 12 | ||||
-rw-r--r-- | src/lib/netlist/plib/pconfig.h | 9 | ||||
-rw-r--r-- | src/lib/netlist/plib/pfmtlog.h | 23 | ||||
-rw-r--r-- | src/lib/netlist/plib/putil.h | 66 | ||||
-rw-r--r-- | src/lib/netlist/solver/nld_matrix_solver.h | 21 | ||||
-rw-r--r-- | src/lib/netlist/solver/nld_ms_direct.h | 2 | ||||
-rw-r--r-- | src/lib/netlist/solver/nld_ms_sor.h | 2 | ||||
-rw-r--r-- | src/lib/netlist/solver/nld_ms_sor_mat.h | 4 | ||||
-rw-r--r-- | src/lib/netlist/solver/nld_ms_w.h | 4 | ||||
-rw-r--r-- | src/lib/netlist/solver/nld_solver.cpp | 5 |
11 files changed, 172 insertions, 23 deletions
diff --git a/src/lib/netlist/nl_config.h b/src/lib/netlist/nl_config.h index 8d1cfa898f5..9baad43ac44 100644 --- a/src/lib/netlist/nl_config.h +++ b/src/lib/netlist/nl_config.h @@ -57,7 +57,8 @@ #define NL_USE_COPY_INSTEAD_OF_REFERENCE (0) #endif -/* +/*! Use the truthtable implementation of 7448 instead of the coded device + * * FIXME: Using truthtable is a lot slower than the explicit device * in breakout. Performance drops by 20%. This can be fixed by * setting param USE_DEACTIVATE for the device. @@ -67,7 +68,8 @@ #define NL_USE_TRUTHTABLE_7448 (0) #endif -/* +/*! Use the truthtable implementation of 74107 instead of the coded device + * * FIXME: The truthtable implementation of 74107 (JK-Flipflop) * is included for educational purposes to demonstrate how * to implement state holding devices as truthtables. @@ -78,6 +80,15 @@ #define NL_USE_TRUTHTABLE_74107 (0) #endif +/*! Use the __float128 type for matrix calculations. + * + * Defaults to PUSE_FLOAT128 + */ + +#ifndef NL_USE_FLOAT128 +#define NL_USE_FLOAT128 PUSE_FLOAT128 +#endif + //============================================================ // DEBUGGING //============================================================ @@ -104,7 +115,7 @@ static constexpr const auto NETLIST_CLOCK = NETLIST_INTERNAL_RES; // due to numeric issues //============================================================ -using nl_fptype = float; +using nl_fptype = double; namespace netlist { @@ -160,6 +171,36 @@ namespace netlist static inline constexpr const char * name() noexcept { return "float"; } static inline constexpr const char * suffix() noexcept { return "f"; } }; + +#if (NL_USE_FLOAT128) + /*! Specific constants for __float128 floating point type + */ + template <> + struct fp_constants<__float128> + { +#if 0 + // MAME compile doesn't support Q + static inline constexpr __float128 DIODE_MAXDIFF() noexcept { return 1e100Q; } + static inline constexpr __float128 DIODE_MAXVOLT() noexcept { return 300.0Q; } + + static inline constexpr __float128 TIMESTEP_MAXDIFF() noexcept { return 1e100Q; } + static inline constexpr __float128 TIMESTEP_MINDIV() noexcept { return 1e-60Q; } + + static inline constexpr const char * name() noexcept { return "__float128"; } + static inline constexpr const char * suffix() noexcept { return "Q"; } +#else + static inline constexpr __float128 DIODE_MAXDIFF() noexcept { return static_cast<__float128>(1e100L); } + static inline constexpr __float128 DIODE_MAXVOLT() noexcept { return static_cast<__float128>(300.0L); } + + static inline constexpr __float128 TIMESTEP_MAXDIFF() noexcept { return static_cast<__float128>(1e100L); } + static inline constexpr __float128 TIMESTEP_MINDIV() noexcept { return static_cast<__float128>(1e-60L); } + + static inline constexpr const char * name() noexcept { return "__float128"; } + static inline constexpr const char * suffix() noexcept { return "Q"; } + +#endif + }; +#endif } // namespace netlist #endif /* NLCONFIG_H_ */ diff --git a/src/lib/netlist/plib/gmres.h b/src/lib/netlist/plib/gmres.h index 840da0cfb45..686e1444d36 100644 --- a/src/lib/netlist/plib/gmres.h +++ b/src/lib/netlist/plib/gmres.h @@ -259,7 +259,7 @@ namespace plib m_ht[j][k] = vec_mult<FT>(n, m_v[kp1], m_v[j]); vec_add_mult_scalar(n, m_v[kp1], m_v[j], -m_ht[j][k]); } - m_ht[kp1][k] = std::sqrt(vec_mult2<FT>(n, m_v[kp1])); + m_ht[kp1][k] = plib::sqrt(vec_mult2<FT>(n, m_v[kp1])); // FIXME: comparison to zero if (m_ht[kp1][k] != plib::constants<FT>::zero()) @@ -268,7 +268,7 @@ namespace plib for (std::size_t j = 0; j < k; j++) givens_mult(m_c[j], m_s[j], m_ht[j][k], m_ht[j+1][k]); - const float_type mu = reciprocal(std::hypot(m_ht[k][k], m_ht[kp1][k])); + const float_type mu = reciprocal(plib::hypot(m_ht[k][k], m_ht[kp1][k])); m_c[k] = m_ht[k][k] * mu; m_s[k] = -m_ht[kp1][k] * mu; @@ -277,7 +277,7 @@ namespace plib givens_mult(m_c[k], m_s[k], m_g[k], m_g[kp1]); - FT rho = std::abs(m_g[kp1]); + FT rho = plib::abs(m_g[kp1]); // FIXME .. itr_used = itr_used + 1; @@ -359,12 +359,12 @@ namespace plib ops.solve_inplace(Ax); - const float_type rho_to_accuracy = std::sqrt(vec_mult2<FT>(n, Ax)) / accuracy; + const float_type rho_to_accuracy = plib::sqrt(vec_mult2<FT>(n, Ax)) / accuracy; rho_delta = accuracy * rho_to_accuracy; } else - rho_delta = accuracy * std::sqrt(static_cast<FT>(n)); + rho_delta = accuracy * plib::sqrt(static_cast<FT>(n)); /* * Using @@ -389,7 +389,7 @@ namespace plib ops.solve_inplace(residual); - rho = std::sqrt(vec_mult2<FT>(n, residual)); + rho = plib::sqrt(vec_mult2<FT>(n, residual)); if (rho < rho_delta) return itr_used + 1; diff --git a/src/lib/netlist/plib/pconfig.h b/src/lib/netlist/plib/pconfig.h index 102a991e613..58a9387c51c 100644 --- a/src/lib/netlist/plib/pconfig.h +++ b/src/lib/netlist/plib/pconfig.h @@ -33,6 +33,15 @@ #define PHAS_INT128 (0) #endif +/*! Add support for the __float128 floating point type. + * + * + */ + +#ifndef PUSE_FLOAT128 +#define PUSE_FLOAT128 (0) +#endif + /* * OpenMP adds about 10% to 20% performance for analog * netlists like kidniki. diff --git a/src/lib/netlist/plib/pfmtlog.h b/src/lib/netlist/plib/pfmtlog.h index 53ed539078e..fd1f073cb73 100644 --- a/src/lib/netlist/plib/pfmtlog.h +++ b/src/lib/netlist/plib/pfmtlog.h @@ -28,11 +28,22 @@ P_ENUM(plog_level, template <typename T> struct ptype_traits_base { - static const T cast(const T &x) { return x; } - static const bool is_signed = std::numeric_limits<T>::is_signed; + static constexpr const T cast(const T &x) { return x; } + static constexpr const bool is_signed = std::numeric_limits<T>::is_signed; static char32_t fmt_spec() { return 'u'; } }; +#if (PUSE_FLOAT128) +template <> +struct ptype_traits_base<__float128> +{ + // FIXME: need native support at some time + static constexpr const long double cast(const __float128 &x) { return static_cast<long double>(x); } + static constexpr const bool is_signed = true; + static char32_t fmt_spec() { return 'f'; } +}; +#endif + template <> struct ptype_traits_base<bool> { @@ -133,6 +144,14 @@ struct ptype_traits<long double> : ptype_traits_base<long double> static char32_t fmt_spec() { return 'f'; } }; +#if (PUSE_FLOAT128) +template<> +struct ptype_traits<__float128> : ptype_traits_base<__float128> +{ + static char32_t fmt_spec() { return 'f'; } +}; +#endif + template<> struct ptype_traits<char *> : ptype_traits_base<char *> { diff --git a/src/lib/netlist/plib/putil.h b/src/lib/netlist/plib/putil.h index e72eccfe01f..e575a9adc73 100644 --- a/src/lib/netlist/plib/putil.h +++ b/src/lib/netlist/plib/putil.h @@ -18,6 +18,11 @@ #include <locale> #include <sstream> #include <vector> +#include <cmath> + +#if (PUSE_FLOAT128) +#include <quadmath.h> +#endif #define PSTRINGIFY_HELP(y) # y #define PSTRINGIFY(x) PSTRINGIFY_HELP(x) @@ -281,6 +286,67 @@ namespace plib return constants<T>::one() / v; } + /*! abs function + * + * @tparam T type of the argument + * @param v argument + * @return absolute value of argument + */ + template <typename T> + static inline constexpr typename std::enable_if<std::is_floating_point<T>::value, T>::type + abs(T v) noexcept + { + return std::abs(v); + } + + /*! sqrt function + * + * @tparam T type of the argument + * @param v argument + * @return absolute value of argument + */ + template <typename T> + static inline constexpr typename std::enable_if<std::is_floating_point<T>::value, T>::type + sqrt(T v) noexcept + { + return std::sqrt(v); + } + + /*! hypot function + * + * @tparam T type of the argument + * @param v argument + * @return absolute value of argument + */ + template <typename T> + static inline constexpr typename std::enable_if<std::is_floating_point<T>::value, T>::type + hypot(T v1, T v2) noexcept + { + return std::hypot(v1, v2); + } + +#if (PUSE_FLOAT128) + static inline constexpr __float128 reciprocal(__float128 v) noexcept + { + return constants<__float128>::one() / v; + } + + static inline __float128 abs(__float128 v) noexcept + { + return fabsq(v); + } + + static inline __float128 sqrt(__float128 v) noexcept + { + return sqrtq(v); + } + + static inline __float128 hypot(__float128 v1, __float128 v2) noexcept + { + return hypotq(v1, v2); + } +#endif + static_assert(noexcept(constants<double>::one()) == true, "Not evaluated as constexpr"); diff --git a/src/lib/netlist/solver/nld_matrix_solver.h b/src/lib/netlist/solver/nld_matrix_solver.h index 09e650be1a3..26f6ea149a7 100644 --- a/src/lib/netlist/solver/nld_matrix_solver.h +++ b/src/lib/netlist/solver/nld_matrix_solver.h @@ -40,11 +40,20 @@ namespace solver GMRES ) +#if (NL_USE_FLOAT128) P_ENUM(matrix_fp_type_e, - FLOAT, - DOUBLE, - LONGDOUBLE + FLOAT + , DOUBLE + , LONGDOUBLE + , FLOAT128 ) +#else + P_ENUM(matrix_fp_type_e, + FLOAT + , DOUBLE + , LONGDOUBLE + ) +#endif struct solver_parameters_t { @@ -56,7 +65,7 @@ namespace solver , m_method(parent, "METHOD", matrix_type_e::MAT_CR) , m_fp_type(parent, "FPTYPE", matrix_fp_type_e::DOUBLE) , m_reltol(parent, "RELTOL", nlconst::magic(1e-3)) ///< SPICE RELTOL parameter - , m_vntol(parent, "VNTOL", nlconst::magic(1e-6)) ///< SPICE VNTOL parameter + , m_vntol(parent, "VNTOL", nlconst::magic(1e-7)) ///< SPICE VNTOL parameter , m_accuracy(parent, "ACCURACY", nlconst::magic(1e-7)) ///< Iterative solver accuracy , m_nr_loops(parent, "NR_LOOPS", 250) ///< Maximum number of Newton-Raphson loops , m_gs_loops(parent, "GS_LOOPS", 9) ///< Maximum number of Gauss-Seidel loops @@ -438,8 +447,8 @@ namespace solver { const auto vold(this->m_terms[i].template getV<FT>()); const auto vnew(m_new_V[i]); - const auto tol(vntol + reltol * std::max(std::abs(vnew),std::abs(vold))); - if (std::abs(vnew - vold) > tol) + const auto tol(vntol + reltol * std::max(plib::abs(vnew),plib::abs(vold))); + if (plib::abs(vnew - vold) > tol) return true; } return false; diff --git a/src/lib/netlist/solver/nld_ms_direct.h b/src/lib/netlist/solver/nld_ms_direct.h index 20f9fa23964..031cae34a27 100644 --- a/src/lib/netlist/solver/nld_ms_direct.h +++ b/src/lib/netlist/solver/nld_ms_direct.h @@ -89,7 +89,7 @@ namespace solver std::size_t maxrow = i; for (std::size_t j = i + 1; j < kN; j++) { - if (std::abs(m_A[j][i]) > std::abs(m_A[maxrow][i])) + if (plib::abs(m_A[j][i]) > plib::abs(m_A[maxrow][i])) //if (m_A[j][i] * m_A[j][i] > m_A[maxrow][i] * m_A[maxrow][i]) maxrow = j; } diff --git a/src/lib/netlist/solver/nld_ms_sor.h b/src/lib/netlist/solver/nld_ms_sor.h index ee79599a2ff..56075fa090e 100644 --- a/src/lib/netlist/solver/nld_ms_sor.h +++ b/src/lib/netlist/solver/nld_ms_sor.h @@ -134,7 +134,7 @@ namespace solver const float_type new_val = this->m_new_V[k] * one_m_w[k] + (Idrive + this->m_RHS[k]) * w[k]; - err = std::max(std::abs(new_val - this->m_new_V[k]), err); + err = std::max(plib::abs(new_val - this->m_new_V[k]), err); this->m_new_V[k] = new_val; } diff --git a/src/lib/netlist/solver/nld_ms_sor_mat.h b/src/lib/netlist/solver/nld_ms_sor_mat.h index e411b1f1e59..a5ded1a38ba 100644 --- a/src/lib/netlist/solver/nld_ms_sor_mat.h +++ b/src/lib/netlist/solver/nld_ms_sor_mat.h @@ -180,7 +180,7 @@ namespace solver FT gabs_t = plib::constants<FT>::zero(); for (std::size_t i = 0; i < e; i++) if (p[i] != k) - gabs_t = gabs_t + std::abs(this->m_A[k][p[i]]); + gabs_t = gabs_t + plib::abs(this->m_A[k][p[i]]); gabs_t *= plib::constants<FT>::one(); // derived by try and error if (gabs_t > this->m_A[k][k]) @@ -190,7 +190,7 @@ namespace solver } const float_type delta = w * (this->m_RHS[k] - Idrive) ; - cerr = std::max(cerr, std::abs(delta)); + cerr = std::max(cerr, plib::abs(delta)); this->m_new_V[k] += delta; } diff --git a/src/lib/netlist/solver/nld_ms_w.h b/src/lib/netlist/solver/nld_ms_w.h index b0c148bf252..5e47e93d171 100644 --- a/src/lib/netlist/solver/nld_ms_w.h +++ b/src/lib/netlist/solver/nld_ms_w.h @@ -340,8 +340,8 @@ namespace solver { tmp += A(i,j) * this->m_new_V[j]; } - if (std::abs(tmp-RHS(i)) > static_cast<float_type>(1e-6)) - plib::perrlogger("{} failed on row {}: {} RHS: {}\n", this->name(), i, std::abs(tmp-RHS(i)), RHS(i)); + if (plib::abs(tmp-RHS(i)) > static_cast<float_type>(1e-6)) + plib::perrlogger("{} failed on row {}: {} RHS: {}\n", this->name(), i, plib::abs(tmp-RHS(i)), RHS(i)); } bool err(false); diff --git a/src/lib/netlist/solver/nld_solver.cpp b/src/lib/netlist/solver/nld_solver.cpp index 86772fac6a0..ecbf48ab0c2 100644 --- a/src/lib/netlist/solver/nld_solver.cpp +++ b/src/lib/netlist/solver/nld_solver.cpp @@ -343,6 +343,11 @@ namespace devices case solver::matrix_fp_type_e::LONGDOUBLE: ms = create_solvers<long double>(sname, grp); break; +#if (NL_USE_FLOAT128) + case solver::matrix_fp_type_e::FLOAT128: + ms = create_solvers<__float128>(sname, grp); + break; +#endif } log().verbose("Solver {1}", ms->name()); |