summaryrefslogtreecommitdiffstatshomepage
diff options
context:
space:
mode:
author couriersud <couriersud@gmx.org>2019-11-02 12:25:50 +0100
committer couriersud <couriersud@gmx.org>2019-11-02 12:25:50 +0100
commit6573037934e93af6d103fd532049b28af0d01c8f (patch)
tree8b749a840e9cc34b4fa4e48a7841c6679ba9596c
parent248d3629d086f821bd01c6439fb6fb2ba53a68bd (diff)
netlist: Support for float, double and long double solvers. [Couriersud]
- Added new solver parameter FPTYPE. This determines in which floating point domain the linear system is solved. May be one of "FLOAT", "DOUBLE" or "LONGDOUBLE" - Added option "--fperr" to nltool. This enables floating point exceptions. This helps debugging the code under gdb. The purpose of this going forward is to have more choice in optimization. Non-dynamic systems should be just fine in the float domain. Dynamic systems (i.e. diodes, bjts, mosfets) should in general work with double. Certain edge cases may require long double resolution.
-rw-r--r--src/lib/netlist/analog/nld_bjt.cpp8
-rw-r--r--src/lib/netlist/analog/nld_generic_models.h4
-rw-r--r--src/lib/netlist/analog/nld_mosfet.cpp4
-rw-r--r--src/lib/netlist/analog/nld_switches.cpp2
-rw-r--r--src/lib/netlist/analog/nlid_twoterm.cpp2
-rw-r--r--src/lib/netlist/analog/nlid_twoterm.h4
-rw-r--r--src/lib/netlist/devices/nld_74123.cpp2
-rw-r--r--src/lib/netlist/devices/nld_74ls629.cpp2
-rw-r--r--src/lib/netlist/devices/nld_log.cpp6
-rw-r--r--src/lib/netlist/devices/nld_mm5837.cpp4
-rw-r--r--src/lib/netlist/devices/nld_r2r_dac.cpp2
-rw-r--r--src/lib/netlist/devices/nld_system.cpp2
-rw-r--r--src/lib/netlist/devices/nlid_system.h22
-rw-r--r--src/lib/netlist/nl_config.h26
-rw-r--r--src/lib/netlist/plib/gmres.h25
-rw-r--r--src/lib/netlist/plib/mat_cr.h8
-rw-r--r--src/lib/netlist/plib/pfmtlog.h6
-rw-r--r--src/lib/netlist/plib/ptime.h28
-rw-r--r--src/lib/netlist/plib/putil.h13
-rw-r--r--src/lib/netlist/prg/nltool.cpp15
-rw-r--r--src/lib/netlist/solver/nld_matrix_solver.cpp10
-rw-r--r--src/lib/netlist/solver/nld_matrix_solver.h55
-rw-r--r--src/lib/netlist/solver/nld_ms_direct.h12
-rw-r--r--src/lib/netlist/solver/nld_ms_direct1.h4
-rw-r--r--src/lib/netlist/solver/nld_ms_direct2.h4
-rw-r--r--src/lib/netlist/solver/nld_ms_gcr.h6
-rw-r--r--src/lib/netlist/solver/nld_ms_gmres.h10
-rw-r--r--src/lib/netlist/solver/nld_ms_sm.h30
-rw-r--r--src/lib/netlist/solver/nld_ms_sor.h36
-rw-r--r--src/lib/netlist/solver/nld_ms_sor_mat.h16
-rw-r--r--src/lib/netlist/solver/nld_ms_w.h43
-rw-r--r--src/lib/netlist/solver/nld_solver.cpp212
-rw-r--r--src/lib/netlist/solver/nld_solver.h4
33 files changed, 364 insertions, 263 deletions
diff --git a/src/lib/netlist/analog/nld_bjt.cpp b/src/lib/netlist/analog/nld_bjt.cpp
index 0402f94e8b9..f1d66c2344f 100644
--- a/src/lib/netlist/analog/nld_bjt.cpp
+++ b/src/lib/netlist/analog/nld_bjt.cpp
@@ -21,18 +21,18 @@ namespace analog
class diode
{
public:
- diode() : m_Is(1e-15), m_VT(0.0258), m_VT_inv(1.0 / m_VT) {}
+ diode() : m_Is(1e-15), m_VT(0.0258), m_VT_inv(plib::reciprocal(m_VT)) {}
diode(const nl_fptype Is, const nl_fptype n)
{
m_Is = Is;
m_VT = 0.0258 * n;
- m_VT_inv = 1.0 / m_VT;
+ m_VT_inv = plib::reciprocal(m_VT);
}
void set(const nl_fptype Is, const nl_fptype n)
{
m_Is = Is;
m_VT = 0.0258 * n;
- m_VT_inv = 1.0 / m_VT;
+ m_VT_inv = plib::reciprocal(m_VT);
}
nl_fptype I(const nl_fptype V) const { return m_Is * std::exp(V * m_VT_inv) - m_Is; }
nl_fptype g(const nl_fptype V) const { return m_Is * m_VT_inv * std::exp(V * m_VT_inv); }
@@ -334,7 +334,7 @@ namespace analog
/* Base current is 0.005 / beta
* as a rough estimate, we just scale the conductance down */
- m_gB = 1.0 / (m_V/(0.005 / BF));
+ m_gB = plib::reciprocal((m_V/(0.005 / BF)));
//m_gB = d.gI(0.005 / alpha);
diff --git a/src/lib/netlist/analog/nld_generic_models.h b/src/lib/netlist/analog/nld_generic_models.h
index f33168594fb..fe56e25e674 100644
--- a/src/lib/netlist/analog/nld_generic_models.h
+++ b/src/lib/netlist/analog/nld_generic_models.h
@@ -71,7 +71,7 @@ namespace analog
void timestep(nl_fptype cap, nl_fptype v, nl_fptype step)
{
- m_h = 1.0 / step;
+ m_h = plib::reciprocal(step);
m_c = cap;
m_v = v;
}
@@ -108,7 +108,7 @@ namespace analog
void timestep(nl_fptype cap, nl_fptype v, nl_fptype step)
{
plib::unused_var(cap);
- m_h = 1.0 / step;
+ m_h = plib::reciprocal(step);
m_v = v;
}
void setparams(nl_fptype gmin) { m_gmin = gmin; }
diff --git a/src/lib/netlist/analog/nld_mosfet.cpp b/src/lib/netlist/analog/nld_mosfet.cpp
index 3bd4d4fa2c4..e78de46252f 100644
--- a/src/lib/netlist/analog/nld_mosfet.cpp
+++ b/src/lib/netlist/analog/nld_mosfet.cpp
@@ -449,9 +449,9 @@ namespace analog
const nl_fptype k = 3.5; // see "Circuit Simulation", page 185
nl_fptype d = (Vgs - m_Vgs);
- Vgs = m_Vgs + 1.0/k * (d < 0 ? -1.0 : 1.0) * std::log1p(k * std::abs(d));
+ Vgs = m_Vgs + plib::reciprocal(k) * (d < 0 ? -1.0 : 1.0) * std::log1p(k * std::abs(d));
d = (Vgd - m_Vgd);
- Vgd = m_Vgd + 1.0/k * (d < 0 ? -1.0 : 1.0) * std::log1p(k * std::abs(d));
+ Vgd = m_Vgd + plib::reciprocal(k) * (d < 0 ? -1.0 : 1.0) * std::log1p(k * std::abs(d));
m_Vgs = Vgs;
m_Vgd = Vgd;
diff --git a/src/lib/netlist/analog/nld_switches.cpp b/src/lib/netlist/analog/nld_switches.cpp
index 0e62fe13952..9b73eb9650f 100644
--- a/src/lib/netlist/analog/nld_switches.cpp
+++ b/src/lib/netlist/analog/nld_switches.cpp
@@ -12,7 +12,7 @@
/* FIXME : convert to parameters */
-#define R_OFF (1.0 / exec().gmin())
+#define R_OFF (plib::reciprocal(exec().gmin()))
#define R_ON 0.01
namespace netlist
diff --git a/src/lib/netlist/analog/nlid_twoterm.cpp b/src/lib/netlist/analog/nlid_twoterm.cpp
index bc63e625c1a..fe928782810 100644
--- a/src/lib/netlist/analog/nlid_twoterm.cpp
+++ b/src/lib/netlist/analog/nlid_twoterm.cpp
@@ -53,7 +53,7 @@ NETLIB_UPDATE(twoterm)
NETLIB_RESET(R_base)
{
NETLIB_NAME(twoterm)::reset();
- set_R(1.0 / exec().gmin());
+ set_R(plib::reciprocal(exec().gmin()));
}
// ----------------------------------------------------------------------------------------
diff --git a/src/lib/netlist/analog/nlid_twoterm.h b/src/lib/netlist/analog/nlid_twoterm.h
index 2618c3eb20e..cb80e15a34a 100644
--- a/src/lib/netlist/analog/nlid_twoterm.h
+++ b/src/lib/netlist/analog/nlid_twoterm.h
@@ -411,7 +411,7 @@ namespace analog
{
m_t += step;
m_funcparam[0] = m_t;
- this->set_G_V_I(1.0 / m_R(),
+ this->set_G_V_I(plib::reciprocal(m_R()),
m_compiled.evaluate(m_funcparam),
0.0);
}
@@ -422,7 +422,7 @@ namespace analog
NETLIB_RESETI()
{
NETLIB_NAME(twoterm)::reset();
- this->set_G_V_I(1.0 / m_R(), m_V(), 0.0);
+ this->set_G_V_I(plib::reciprocal(m_R()), m_V(), 0.0);
}
private:
diff --git a/src/lib/netlist/devices/nld_74123.cpp b/src/lib/netlist/devices/nld_74123.cpp
index 213fa9cc396..067e6ca1918 100644
--- a/src/lib/netlist/devices/nld_74123.cpp
+++ b/src/lib/netlist/devices/nld_74123.cpp
@@ -255,7 +255,7 @@ namespace netlist
NETLIB_RESET(74123)
{
- m_KP = 1.0 / (1.0 + exp(m_K()));
+ m_KP = plib::reciprocal(1.0 + exp(m_K()));
m_RP.reset();
m_RN.reset();
diff --git a/src/lib/netlist/devices/nld_74ls629.cpp b/src/lib/netlist/devices/nld_74ls629.cpp
index 0f1c64fc74a..42c9292315c 100644
--- a/src/lib/netlist/devices/nld_74ls629.cpp
+++ b/src/lib/netlist/devices/nld_74ls629.cpp
@@ -212,7 +212,7 @@ namespace netlist
// FIXME: we need a possibility to remove entries from queue ...
// or an exact model ...
- m_clock.m_inc = netlist_time::from_double(0.5 / freq);
+ m_clock.m_inc = netlist_time::from_fp<nl_fptype>(0.5 / freq);
//m_clock.update();
//NL_VERBOSE_OUT(("{1} {2} {3} {4}\n", name(), v_freq, v_rng, freq));
diff --git a/src/lib/netlist/devices/nld_log.cpp b/src/lib/netlist/devices/nld_log.cpp
index 3b173060986..03165a573ef 100644
--- a/src/lib/netlist/devices/nld_log.cpp
+++ b/src/lib/netlist/devices/nld_log.cpp
@@ -31,7 +31,7 @@ namespace netlist
NETLIB_UPDATEI()
{
/* use pstring::sprintf, it is a LOT faster */
- m_writer.writeline(plib::pfmt("{1:.9} {2}").e(exec().time().as_double()).e(static_cast<nl_fptype>(m_I())));
+ m_writer.writeline(plib::pfmt("{1:.9} {2}").e(exec().time().as_fp<nl_fptype>()).e(static_cast<nl_fptype>(m_I())));
}
NETLIB_RESETI() { }
@@ -50,7 +50,7 @@ namespace netlist
NETLIB_UPDATEI()
{
- m_writer.writeline(plib::pfmt("{1:.9} {2}").e(exec().time().as_double()).e(static_cast<nl_fptype>(m_I() - m_I2())));
+ m_writer.writeline(plib::pfmt("{1:.9} {2}").e(exec().time().as_fp<nl_fptype>()).e(static_cast<nl_fptype>(m_I() - m_I2())));
}
NETLIB_RESETI() { }
@@ -83,7 +83,7 @@ namespace netlist
NETLIB_UPDATE(wav)
{
- fprintf(m_file, "%e %e\n", netlist().time().as_double(), m_I());
+ fprintf(m_file, "%e %e\n", netlist().time().as_fp<nl_fptype>(), m_I());
}
NETLIB_NAME(log)::~NETLIB_NAME(wav)()
diff --git a/src/lib/netlist/devices/nld_mm5837.cpp b/src/lib/netlist/devices/nld_mm5837.cpp
index 2c58354fd3f..5d32c8ff8cd 100644
--- a/src/lib/netlist/devices/nld_mm5837.cpp
+++ b/src/lib/netlist/devices/nld_mm5837.cpp
@@ -70,7 +70,7 @@ namespace netlist
//m_V0.initial(0.0);
//m_RV.do_reset();
m_RV.set_G_V_I(plib::constants<nl_fptype>::one() / R_LOW, 0.0, 0.0);
- m_inc = netlist_time::from_double(1.0 / m_FREQ());
+ m_inc = netlist_time::from_fp(plib::reciprocal(m_FREQ()));
if (m_FREQ() < 24000 || m_FREQ() > 56000)
log().warning(MW_FREQUENCY_OUTSIDE_OF_SPECS_1(m_FREQ()));
@@ -80,7 +80,7 @@ namespace netlist
NETLIB_UPDATE_PARAM(MM5837_dip)
{
- m_inc = netlist_time::from_double(1.0 / m_FREQ());
+ m_inc = netlist_time::from_fp(plib::reciprocal(m_FREQ()));
if (m_FREQ() < 24000 || m_FREQ() > 56000)
log().warning(MW_FREQUENCY_OUTSIDE_OF_SPECS_1(m_FREQ()));
}
diff --git a/src/lib/netlist/devices/nld_r2r_dac.cpp b/src/lib/netlist/devices/nld_r2r_dac.cpp
index e2695766332..c1cc4462b7b 100644
--- a/src/lib/netlist/devices/nld_r2r_dac.cpp
+++ b/src/lib/netlist/devices/nld_r2r_dac.cpp
@@ -44,7 +44,7 @@ namespace netlist
nl_fptype V = m_VIN() / static_cast<nl_fptype>(1 << m_num())
* static_cast<nl_fptype>(m_val());
- this->set_G_V_I(1.0 / m_R(), V, 0.0);
+ this->set_G_V_I(plib::reciprocal(m_R()), V, 0.0);
}
} //namespace analog
diff --git a/src/lib/netlist/devices/nld_system.cpp b/src/lib/netlist/devices/nld_system.cpp
index 427c2f3dcae..54ef14fa1d0 100644
--- a/src/lib/netlist/devices/nld_system.cpp
+++ b/src/lib/netlist/devices/nld_system.cpp
@@ -25,7 +25,7 @@ namespace netlist
NETLIB_RESET(extclock)
{
m_cnt = 0;
- m_off = netlist_time::from_double(m_offset());
+ m_off = netlist_time::from_fp<decltype(m_offset())>(m_offset());
m_feedback.m_delegate = NETLIB_DELEGATE(extclock, update);
//m_feedback.m_delegate .set(&NETLIB_NAME(extclock)::update, this);
diff --git a/src/lib/netlist/devices/nlid_system.h b/src/lib/netlist/devices/nlid_system.h
index c51f5c36c56..d57de8f4750 100644
--- a/src/lib/netlist/devices/nlid_system.h
+++ b/src/lib/netlist/devices/nlid_system.h
@@ -54,7 +54,7 @@ namespace devices
, m_Q(*this, "Q")
, m_freq(*this, "FREQ", 7159000.0 * 5)
{
- m_inc = netlist_time::from_double(1.0 / (m_freq()*2.0));
+ m_inc = netlist_time::from_fp<decltype(m_freq())>(1.0 / (m_freq()*2.0));
}
NETLIB_RESETI()
@@ -64,7 +64,7 @@ namespace devices
NETLIB_UPDATE_PARAMI()
{
- m_inc = netlist_time::from_double(1.0 / (m_freq()*2.0));
+ m_inc = netlist_time::from_fp<decltype(m_freq())>(1.0 / (m_freq()*2.0));
}
NETLIB_UPDATEI()
@@ -93,7 +93,7 @@ namespace devices
, m_Q(*this, "Q")
, m_freq(*this, "FREQ", 7159000.0 * 5.0)
{
- m_inc = netlist_time::from_double(1.0 / (m_freq()*2.0));
+ m_inc = netlist_time::from_fp<decltype(m_freq())>(1.0 / (m_freq()*2.0));
connect(m_feedback, m_Q);
}
@@ -101,7 +101,7 @@ namespace devices
NETLIB_UPDATE_PARAMI()
{
- m_inc = netlist_time::from_double(1.0 / (m_freq() * 2.0));
+ m_inc = netlist_time::from_fp<decltype(m_freq())>(1.0 / (m_freq() * 2.0));
}
NETLIB_UPDATEI()
@@ -139,8 +139,8 @@ namespace devices
NETLIB_UPDATEI()
{
- m_funcparam[0] = exec().time().as_double();
- const netlist_time m_inc = netlist_time::from_double(m_compiled.evaluate(m_funcparam));
+ m_funcparam[0] = exec().time().as_fp<nl_fptype>();
+ const netlist_time m_inc = netlist_time::from_fp(m_compiled.evaluate(m_funcparam));
m_Q.push(!m_feedback(), m_inc);
}
@@ -168,13 +168,13 @@ namespace devices
, m_cnt(*this, "m_cnt", 0)
, m_off(*this, "m_off", netlist_time::zero())
{
- m_inc[0] = netlist_time::from_double(1.0 / (m_freq() * 2.0));
+ m_inc[0] = netlist_time::from_fp<decltype(m_freq())>(1.0 / (m_freq() * 2.0));
connect(m_feedback, m_Q);
- netlist_time base = netlist_time::from_double(1.0 / (m_freq()*2.0));
+ netlist_time base = netlist_time::from_fp<decltype(m_freq())>(1.0 / (m_freq()*2.0));
std::vector<pstring> pat(plib::psplit(m_pattern(),","));
- m_off = netlist_time::from_double(m_offset());
+ m_off = netlist_time::from_fp(m_offset());
std::array<std::int64_t, 32> pati = { 0 };
@@ -331,8 +331,8 @@ namespace devices
NETLIB_RESETI()
{
- m_RIN.set_G_V_I(1.0 / m_p_RIN(),0,0);
- m_ROUT.set_G_V_I(1.0 / m_p_ROUT(),0,0);
+ m_RIN.set_G_V_I(plib::reciprocal(m_p_RIN()),0,0);
+ m_ROUT.set_G_V_I(plib::reciprocal(m_p_ROUT()),0,0);
}
NETLIB_UPDATEI()
diff --git a/src/lib/netlist/nl_config.h b/src/lib/netlist/nl_config.h
index 80f74cb1582..1633efc523b 100644
--- a/src/lib/netlist/nl_config.h
+++ b/src/lib/netlist/nl_config.h
@@ -105,10 +105,6 @@ static constexpr const auto NETLIST_CLOCK = NETLIST_INTERNAL_RES;
//============================================================
using nl_fptype = double;
-//using nl_fptype = float;
-
-using nl_mat_fptype = nl_fptype;
-//using nl_mat_fptype = float;
namespace netlist
{
@@ -120,6 +116,20 @@ namespace netlist
struct fp_constants
{ };
+ /*! Specific constants for long double floating point type
+ */
+ template <>
+ struct fp_constants<long double>
+ {
+ static inline constexpr long double DIODE_MAXDIFF() noexcept { return 1e100L; }
+ static inline constexpr long double DIODE_MAXVOLT() noexcept { return 300.0L; }
+
+ static inline constexpr long double TIMESTEP_MAXDIFF() noexcept { return 1e100L; }
+ static inline constexpr long double TIMESTEP_MINDIV() noexcept { return 1e-60L; }
+
+ static inline constexpr const char * name() noexcept { return "long double"; }
+ };
+
/*! Specific constants for double floating point type
*/
template <>
@@ -130,6 +140,8 @@ namespace netlist
static inline constexpr double TIMESTEP_MAXDIFF() noexcept { return 1e100; }
static inline constexpr double TIMESTEP_MINDIV() noexcept { return 1e-60; }
+
+ static inline constexpr const char * name() noexcept { return "double"; }
};
/*! Specific constants for float floating point type
@@ -137,11 +149,13 @@ namespace netlist
template <>
struct fp_constants<float>
{
- static inline constexpr float DIODE_MAXDIFF() noexcept { return 1e5; }
- static inline constexpr float DIODE_MAXVOLT() noexcept { return 30.0; }
+ static inline constexpr float DIODE_MAXDIFF() noexcept { return 1e2; }
+ static inline constexpr float DIODE_MAXVOLT() noexcept { return 20.0; }
static inline constexpr float TIMESTEP_MAXDIFF() noexcept { return 1e30f; }
static inline constexpr float TIMESTEP_MINDIV() noexcept { return 1e-8f; }
+
+ static inline constexpr const char * name() noexcept { return "float"; }
};
} // namespace netlist
diff --git a/src/lib/netlist/plib/gmres.h b/src/lib/netlist/plib/gmres.h
index 8c5c724c581..840da0cfb45 100644
--- a/src/lib/netlist/plib/gmres.h
+++ b/src/lib/netlist/plib/gmres.h
@@ -126,7 +126,7 @@ namespace plib
// 136%
for (std::size_t j = m_mat.row_idx[i]; j< m_mat.row_idx[i+1]; j++)
v += m_mat.A[j] * m_mat.A[j];
- m_diag[i] = 1.0 / std::sqrt(v);
+ m_diag[i] = reciprocal(std::sqrt(v));
#elif 0
// works halfway, i.e. Mame perforamnce 50%
// 147% - lowest average solution time with 7.094
@@ -151,14 +151,14 @@ namespace plib
// 151%
for (std::size_t j = m_mat.row_idx[i]; j< m_mat.row_idx[i+1]; j++)
v += std::abs(m_mat.A[j]);
- m_diag[i] = 1.0 / v;
+ m_diag[i] = reciprocal(v);
#else
// 124%
for (std::size_t j = m_mat.row_idx[i]; j< m_mat.row_idx[i+1]; j++)
v = std::max(v, std::abs(m_mat.A[j]));
- m_diag[i] = 1.0 / v;
+ m_diag[i] = reciprocal(v);
#endif
- //m_diag[i] = 1.0 / m_mat.A[m_mat.diag[i]];
+ //m_diag[i] = reciprocal(m_mat.A[m_mat.diag[i]]);
}
}
@@ -261,18 +261,19 @@ namespace plib
}
m_ht[kp1][k] = std::sqrt(vec_mult2<FT>(n, m_v[kp1]));
- if (m_ht[kp1][k] != 0.0)
- vec_scale(n, m_v[kp1], constants<FT>::one() / m_ht[kp1][k]);
+ // FIXME: comparison to zero
+ if (m_ht[kp1][k] != plib::constants<FT>::zero())
+ vec_scale(n, m_v[kp1], reciprocal(m_ht[kp1][k]));
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 = 1.0 / std::hypot(m_ht[k][k], m_ht[kp1][k]);
+ const float_type mu = reciprocal(std::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;
m_ht[k][k] = m_c[k] * m_ht[k][k] - m_s[k] * m_ht[kp1][k];
- m_ht[kp1][k] = 0.0;
+ m_ht[kp1][k] = plib::constants<FT>::zero();
givens_mult(m_c[k], m_s[k], m_g[k], m_g[kp1]);
@@ -287,7 +288,7 @@ namespace plib
/* x += m_v[j] * m_y[j] */
for (std::size_t i = k + 1; i-- > 0;)
{
- double tmp = m_g[i];
+ auto tmp(m_g[i]);
for (std::size_t j = i + 1; j <= k; j++)
tmp -= m_ht[i][j] * m_y[j];
m_y[i] = tmp / m_ht[i][i];
@@ -334,7 +335,7 @@ namespace plib
*------------------------------------------------------------------------*/
std::size_t itr_used = 0;
- double rho_delta = 0.0;
+ auto rho_delta(plib::constants<float_type>::zero());
const std::size_t n = size();
@@ -400,7 +401,7 @@ namespace plib
vec_set_scalar(RESTART+1, m_g, +constants<FT>::zero());
m_g[0] = rho;
- vec_mult_scalar(n, m_v[0], residual, constants<FT>::one() / rho);
+ vec_mult_scalar(n, m_v[0], residual, reciprocal(rho));
if (do_k<RESTART-1>(ops, x, itr_used, rho_delta, true))
// converged
@@ -504,7 +505,7 @@ namespace plib
else
{
beta = alpha * ( c / 2.0)*( c / 2.0);
- alpha = 1.0 / (d - beta);
+ alpha = reciprocal(d - beta);
for (std::size_t k = 0; k < size(); k++)
p[k] = residual[k] + beta * p[k];
}
diff --git a/src/lib/netlist/plib/mat_cr.h b/src/lib/netlist/plib/mat_cr.h
index 3a6c59f24ab..132cab5d96f 100644
--- a/src/lib/netlist/plib/mat_cr.h
+++ b/src/lib/netlist/plib/mat_cr.h
@@ -167,7 +167,7 @@ namespace plib
const std::size_t oe = nz_num;
while (k < oe)
{
- T tmp = 0.0;
+ T tmp = plib::constants<T>::zero();
const std::size_t e = row_idx[row+1];
for (; k < e; k++)
tmp += A[k] * x[col_idx[k]];
@@ -294,7 +294,7 @@ namespace plib
{
std::size_t nzbdp = 0;
std::size_t pi = base::diag[i];
- const typename base::value_type f = 1.0 / base::A[pi++];
+ auto f = reciprocal(base::A[pi++]);
const std::size_t piie = base::row_idx[i+1];
const auto &nz = base::nzbd[i];
@@ -345,7 +345,7 @@ namespace plib
{
std::size_t nzbdp = 0;
std::size_t pi = base::diag[i];
- const typename base::value_type f = 1.0 / base::A[pi++];
+ const auto f = reciprocal(base::A[pi++]);
const std::size_t piie = base::row_idx[i+1];
const auto &nz = base::nzbd[i];
@@ -358,7 +358,7 @@ namespace plib
while (base::col_idx[pj] < i)
pj++;
- const typename base::value_type f1 = - base::A[pj++] * f;
+ auto f1 = - base::A[pj++] * f;
// subtract row i from j
// fill-in available assumed, i.e. matrix was prepared
diff --git a/src/lib/netlist/plib/pfmtlog.h b/src/lib/netlist/plib/pfmtlog.h
index 13ec2396a8a..17ca3f1d7ad 100644
--- a/src/lib/netlist/plib/pfmtlog.h
+++ b/src/lib/netlist/plib/pfmtlog.h
@@ -128,6 +128,12 @@ struct ptype_traits<double> : ptype_traits_base<double>
};
template<>
+struct ptype_traits<long double> : ptype_traits_base<long double>
+{
+ static char32_t fmt_spec() { return 'f'; }
+};
+
+template<>
struct ptype_traits<char *> : ptype_traits_base<char *>
{
static const char *cast(const char *x) { return x; }
diff --git a/src/lib/netlist/plib/ptime.h b/src/lib/netlist/plib/ptime.h
index 4044f518e24..c0a6cd5a34b 100644
--- a/src/lib/netlist/plib/ptime.h
+++ b/src/lib/netlist/plib/ptime.h
@@ -99,11 +99,18 @@ namespace plib
}
constexpr internal_type as_raw() const noexcept { return m_time; }
- constexpr double as_double() const noexcept
+
+ template <typename FT>
+ constexpr typename std::enable_if<std::is_floating_point<FT>::value, FT>::type
+ as_fp() const noexcept
{
- return static_cast<double>(m_time) * inv_res;
+ return static_cast<FT>(m_time) * inv_res<FT>();
}
+ constexpr double as_double() const noexcept { return as_fp<double>(); }
+ constexpr double as_float() const noexcept { return as_fp<float>(); }
+ constexpr double as_long_double() const noexcept { return as_fp<long double>(); }
+
// for save states ....
C14CONSTEXPR internal_type *get_internaltype_ptr() noexcept { return &m_time; }
@@ -113,7 +120,19 @@ namespace plib
static constexpr const ptime from_sec(const internal_type s) noexcept { return ptime(s, UINT64_C( 1)); }
static constexpr const ptime from_hz(const internal_type hz) noexcept { return ptime(1 , hz); }
static constexpr const ptime from_raw(const internal_type raw) noexcept { return ptime(raw); }
- static constexpr const ptime from_double(const double t) noexcept { return ptime(static_cast<internal_type>(std::floor(t * static_cast<double>(RES) + 0.5)), RES); }
+
+ template <typename FT>
+ static constexpr const typename std::enable_if<std::is_floating_point<FT>::value, ptime>::type
+ from_fp(const FT t) noexcept { return ptime(static_cast<internal_type>(std::floor(t * static_cast<FT>(RES) + 0.5)), RES); }
+
+ static constexpr const ptime from_double(const double t) noexcept
+ { return from_fp<double>(t); }
+
+ static constexpr const ptime from_float(const float t) noexcept
+ { return from_fp<float>(t); }
+
+ static constexpr const ptime from_long_double(const long double t) noexcept
+ { return from_fp<long double>(t); }
static constexpr const ptime zero() noexcept { return ptime(0, RES); }
static constexpr const ptime quantum() noexcept { return ptime(1, RES); }
@@ -126,7 +145,8 @@ namespace plib
constexpr internal_type in_sec() const noexcept { return m_time / (RES / UINT64_C( 1)); }
private:
- static constexpr const double inv_res = 1.0 / static_cast<double>(RES);
+ template <typename FT>
+ static constexpr FT inv_res() noexcept { return static_cast<FT>(1.0) / static_cast<FT>(RES); }
internal_type m_time;
};
diff --git a/src/lib/netlist/plib/putil.h b/src/lib/netlist/plib/putil.h
index 9a21a223ba9..ae358b960df 100644
--- a/src/lib/netlist/plib/putil.h
+++ b/src/lib/netlist/plib/putil.h
@@ -259,6 +259,19 @@ namespace plib
static inline constexpr const T cast(V &&v) noexcept { return static_cast<T>(v); }
};
+ /*! typesafe reciprocal function
+ *
+ * @tparam T type of the argument
+ * @param v argument
+ * @return reciprocal of argument
+ */
+ template <typename T>
+ static inline constexpr typename std::enable_if<std::is_floating_point<T>::value, T>::type
+ reciprocal(T v) noexcept
+ {
+ return constants<T>::one() / v;
+ }
+
static_assert(noexcept(constants<double>::one()) == true, "Not evaluated as constexpr");
diff --git a/src/lib/netlist/prg/nltool.cpp b/src/lib/netlist/prg/nltool.cpp
index 20a37d24ff1..9ca4579ae4d 100644
--- a/src/lib/netlist/prg/nltool.cpp
+++ b/src/lib/netlist/prg/nltool.cpp
@@ -53,6 +53,8 @@ public:
opt_inp(*this, "i", "input", "", "input file to process (default is none)"),
opt_loadstate(*this,"", "loadstate", "", "load state from file and continue from there"),
opt_savestate(*this,"", "savestate", "", "save state to file at end of run"),
+ opt_fperr(*this, "", "fperr",
+ "raise exception on floating point errors. This is intended to be used during debugging."),
opt_grp5(*this, "Options for convert command", "These options are only used by the convert command."),
opt_type(*this, "y", "type", 0, std::vector<pstring>({"spice","eagle","rinf"}), "type of file to be converted: spice,eagle,rinf"),
@@ -97,6 +99,7 @@ public:
plib::option_str opt_inp;
plib::option_str opt_loadstate;
plib::option_str opt_savestate;
+ plib::option_bool opt_fperr;
plib::option_group opt_grp5;
plib::option_str_limit<unsigned> opt_type;
plib::option_group opt_grp6;
@@ -319,7 +322,7 @@ struct input_t
if (e != 3)
throw netlist::nl_exception(plib::pfmt("error {1} scanning line {2}\n")(e)(line));
m_value = static_cast<nl_fptype>(val);
- m_time = netlist::netlist_time::from_double(t);
+ m_time = netlist::netlist_time::from_fp(t);
m_param = setup.find_param(pstring(buf.data()), true);
}
@@ -395,7 +398,7 @@ void tool_app_t::run()
nt.reset();
inps = read_input(nt.setup(), opt_inp());
- ttr = netlist::netlist_time::from_double(opt_ttr());
+ ttr = netlist::netlist_time::from_fp(opt_ttr());
}
@@ -831,10 +834,6 @@ int tool_app_t::execute()
{
tool_app_t opts;
- /* make SIGFPE actually deliver signals on supoorted platforms */
- plib::fpsignalenabler::global_enable(true);
- plib::fpsignalenabler sigen(plib::FP_ALL & ~plib::FP_INEXACT & ~plib::FP_UNDERFLOW);
-
if (opt_help())
{
pout(usage());
@@ -869,6 +868,9 @@ int tool_app_t::execute()
try
{
+ plib::fpsignalenabler::global_enable(opt_fperr());
+ plib::fpsignalenabler fpprotect(plib::FP_DIVBYZERO | plib::FP_UNDERFLOW | plib::FP_OVERFLOW | plib::FP_INVALID);
+
pstring cmd = opt_cmd.as_string();
if (cmd == "listdevices")
listdevices();
@@ -902,6 +904,7 @@ int tool_app_t::execute()
perr("plib exception caught: {}\n", e.text());
return 2;
}
+
#if 0
std::cout.imbue(std::locale("de_DE.utf8"));
std::cout.imbue(std::locale("C.UTF-8"));
diff --git a/src/lib/netlist/solver/nld_matrix_solver.cpp b/src/lib/netlist/solver/nld_matrix_solver.cpp
index afcdec647cd..a1bb4f2f88e 100644
--- a/src/lib/netlist/solver/nld_matrix_solver.cpp
+++ b/src/lib/netlist/solver/nld_matrix_solver.cpp
@@ -397,13 +397,13 @@ namespace solver
if (m_params.m_dynamic_ts && has_timestep_devices())
{
- m_Q_sync.net().toggle_and_push_to_queue(netlist_time::from_double(m_params.m_min_timestep));
+ m_Q_sync.net().toggle_and_push_to_queue(netlist_time::from_fp(m_params.m_min_timestep));
}
}
void matrix_solver_t::step(const netlist_time &delta)
{
- const nl_fptype dd = delta.as_double();
+ const nl_fptype dd = delta.as_fp<nl_fptype>();
for (auto &d : m_step_devices)
d->timestep(dd);
}
@@ -439,7 +439,7 @@ namespace solver
if (this_resched > 1 && !m_Q_sync.net().is_queued())
{
log().warning(MW_NEWTON_LOOPS_EXCEEDED_ON_NET_1(this->name()));
- m_Q_sync.net().toggle_and_push_to_queue(netlist_time::from_double(m_params.m_nr_recalc_delay));
+ m_Q_sync.net().toggle_and_push_to_queue(netlist_time::from_fp(m_params.m_nr_recalc_delay()));
}
}
else
@@ -447,7 +447,7 @@ namespace solver
this->vsolve_non_dynamic(false);
}
- const netlist_time next_time_step = compute_next_timestep(delta.as_double());
+ const netlist_time next_time_step = compute_next_timestep(delta.as_fp<nl_fptype>());
return next_time_step;
}
@@ -558,7 +558,7 @@ namespace solver
static_cast<nl_fptype>(this->m_stat_newton_raphson) / static_cast<nl_fptype>(this->m_stat_vsolver_calls));
log().verbose(" {1:10} invocations ({2:6.0} Hz) {3:10} gs fails ({4:6.2} %) {5:6.3} average",
this->m_stat_calculations,
- static_cast<nl_fptype>(this->m_stat_calculations) / this->exec().time().as_double(),
+ static_cast<nl_fptype>(this->m_stat_calculations) / this->exec().time().as_fp<nl_fptype>(),
this->m_iterative_fail,
100.0 * static_cast<nl_fptype>(this->m_iterative_fail)
/ static_cast<nl_fptype>(this->m_stat_calculations),
diff --git a/src/lib/netlist/solver/nld_matrix_solver.h b/src/lib/netlist/solver/nld_matrix_solver.h
index 632d038594f..341374b1df0 100644
--- a/src/lib/netlist/solver/nld_matrix_solver.h
+++ b/src/lib/netlist/solver/nld_matrix_solver.h
@@ -40,6 +40,12 @@ namespace solver
GMRES
)
+ P_ENUM(matrix_fp_type_e,
+ FLOAT,
+ DOUBLE,
+ LONGDOUBLE
+ )
+
struct solver_parameters_t
{
solver_parameters_t(device_t &parent)
@@ -48,6 +54,7 @@ namespace solver
/* iteration parameters */
, m_gs_sor(parent, "SOR_FACTOR", 1.059)
, m_method(parent, "METHOD", matrix_type_e::MAT_CR)
+ , m_fp_type(parent, "FPTYPE", matrix_fp_type_e::DOUBLE)
, m_accuracy(parent, "ACCURACY", 1e-7)
, m_gs_loops(parent, "GS_LOOPS", 9) // Gauss-Seidel loops
@@ -55,7 +62,7 @@ namespace solver
, m_gmin(parent, "GMIN", 1e-9)
, m_pivot(parent, "PIVOT", false) // use pivoting - on supported solvers
, m_nr_loops(parent, "NR_LOOPS", 250) // Newton-Raphson loops
- , m_nr_recalc_delay(parent, "NR_RECALC_DELAY", netlist_time::quantum().as_double()) // Delay to next solve attempt if nr loops exceeded
+ , m_nr_recalc_delay(parent, "NR_RECALC_DELAY", netlist_time::quantum().as_fp<nl_fptype>()) // Delay to next solve attempt if nr loops exceeded
, m_parallel(parent, "PARALLEL", 0)
/* automatic time step */
@@ -72,7 +79,7 @@ namespace solver
{
m_min_timestep = m_dynamic_min_ts();
- m_max_timestep = netlist_time::from_double(1.0 / m_freq()).as_double();
+ m_max_timestep = netlist_time::from_fp(plib::reciprocal(m_freq())).as_fp<decltype(m_max_timestep)>();
if (m_dynamic_ts)
{
@@ -87,6 +94,7 @@ namespace solver
param_fp_t m_freq;
param_fp_t m_gs_sor;
param_enum_t<matrix_type_e> m_method;
+ param_enum_t<matrix_fp_type_e> m_fp_type;
param_fp_t m_accuracy;
param_num_t<std::size_t> m_gs_loops;
param_fp_t m_gmin;
@@ -122,9 +130,11 @@ namespace solver
terminal_t **terms() noexcept { return m_terms.data(); }
- nl_fptype getV() const noexcept { return m_net->Q_Analog(); }
+ template <typename FT, typename = std::enable_if<std::is_floating_point<FT>::value, void>>
+ FT getV() const noexcept { return static_cast<FT>(m_net->Q_Analog()); }
- void setV(nl_fptype v) noexcept { m_net->set_Q_Analog(v); }
+ template <typename FT, typename = std::enable_if<std::is_floating_point<FT>::value, void>>
+ void setV(FT v) noexcept { m_net->set_Q_Analog(static_cast<nl_fptype>(v)); }
bool isNet(const analog_net_t * net) const noexcept { return net == m_net; }
@@ -240,7 +250,7 @@ namespace solver
using vtype = typename std::decay<decltype(V[0])>::type;
vtype cerr = 0;
for (std::size_t i = 0; i < iN; i++)
- cerr = std::max(cerr, std::abs(V[i] - static_cast<vtype>(this->m_terms[i].getV())));
+ cerr = std::max(cerr, std::abs(V[i] - this->m_terms[i].getV<vtype>()));
return cerr;
}
@@ -383,9 +393,9 @@ namespace solver
: matrix_solver_t(anetlist, name, nets, params)
, m_dim(size)
, m_mat_ptr(size, this->max_railstart() + 1)
- , m_last_V(size, plib::constants<float_type>::zero())
- , m_DD_n_m_1(size, plib::constants<float_type>::zero())
- , m_h_n_m_1(size, plib::constants<float_type>::zero())
+ , m_last_V(size, plib::constants<nl_fptype>::zero())
+ , m_DD_n_m_1(size, plib::constants<nl_fptype>::zero())
+ , m_h_n_m_1(size, plib::constants<nl_fptype>::zero())
{
/*
* save states
@@ -407,11 +417,11 @@ namespace solver
plib::parray2D<float_type *, SIZE, 0> m_mat_ptr;
/* state - variable time_stepping */
PALIGNAS_VECTOROPT()
- plib::parray<float_type, SIZE> m_last_V;
+ plib::parray<nl_fptype, SIZE> m_last_V;
PALIGNAS_VECTOROPT()
- plib::parray<float_type, SIZE> m_DD_n_m_1;
+ plib::parray<nl_fptype, SIZE> m_DD_n_m_1;
PALIGNAS_VECTOROPT()
- plib::parray<float_type, SIZE> m_h_n_m_1;
+ plib::parray<nl_fptype, SIZE> m_h_n_m_1;
constexpr std::size_t size() const noexcept { return (SIZE > 0) ? static_cast<std::size_t>(SIZE) : m_dim; }
@@ -428,7 +438,7 @@ namespace solver
// avoid floating point exceptions
const nl_fptype DD_n = std::max(-fp_constants<nl_fptype>::TIMESTEP_MAXDIFF(),
- std::min(+fp_constants<nl_fptype>::TIMESTEP_MAXDIFF(),(t.getV() - m_last_V[k])));
+ std::min(+fp_constants<nl_fptype>::TIMESTEP_MAXDIFF(),(t.template getV<nl_fptype>() - m_last_V[k])));
const nl_fptype hn = cur_ts;
//printf("%g %g %g %g\n", DD_n, hn, t.m_DD_n_m_1, t.m_h_n_m_1);
@@ -445,7 +455,7 @@ namespace solver
if (new_net_timestep < new_solver_timestep)
new_solver_timestep = new_net_timestep;
- m_last_V[k] = t.getV();
+ m_last_V[k] = t.template getV<nl_fptype>();
}
if (new_solver_timestep < m_params.m_min_timestep)
{
@@ -457,7 +467,7 @@ namespace solver
/*
* FIXME: Factor 2 below is important. Without, we get timing issues. This must be a bug elsewhere.
*/
- return std::max(netlist_time::from_double(new_solver_timestep), netlist_time::quantum() * 2);
+ return std::max(netlist_time::from_fp(new_solver_timestep), netlist_time::quantum() * 2);
}
@@ -488,12 +498,12 @@ namespace solver
void clear_square_mat(M &m)
{
const std::size_t n = size();
-
for (std::size_t k=0; k < n; k++)
{
auto *p = &(m[k][0]);
+ using mat_elem_type = typename std::decay<decltype(*p)>::type;
for (std::size_t i=0; i < n; i++)
- p[i] = 0.0;
+ p[i] = plib::constants<mat_elem_type>::zero();
}
}
@@ -514,11 +524,14 @@ namespace solver
const auto &Idr = m_Idrn[k];
const auto &cnV = m_connected_net_Vn[k];
+ // FIXME: gonn, gtn and Idr - which float types should they have?
+
for (std::size_t i = 0; i < railstart; i++)
- *tcr_r[i] += go[i];
+ *tcr_r[i] += static_cast<FT>(go[i]);
- typename RT::value_type gtot_t = 0.0;
- typename RT::value_type RHS_t = 0.0;
+ /* use native floattype for now */
+ auto gtot_t(plib::constants<nl_fptype>::zero());
+ auto RHS_t (plib::constants<nl_fptype>::zero());
for (std::size_t i = 0; i < term_count; i++)
{
@@ -536,9 +549,9 @@ namespace solver
RHS_t += (/*m_Idr[i]*/ (- go[i]) * *cnV[i]);
}
- RHS[k] = RHS_t;
+ RHS[k] = static_cast<FT>(RHS_t);
// update diagonal element ...
- *tcr_r[railstart] += gtot_t; //mat.A[mat.diag[k]] += gtot_t;
+ *tcr_r[railstart] += static_cast<FT>(gtot_t); //mat.A[mat.diag[k]] += gtot_t;
}
}
diff --git a/src/lib/netlist/solver/nld_ms_direct.h b/src/lib/netlist/solver/nld_ms_direct.h
index 62332dbe6b5..dc7aea3e743 100644
--- a/src/lib/netlist/solver/nld_ms_direct.h
+++ b/src/lib/netlist/solver/nld_ms_direct.h
@@ -72,7 +72,7 @@ namespace solver
for (std::size_t i = 0; i < kN; i++)
{
/* FIXME: Singular matrix? */
- const FT f = 1.0 / m_A[i][i];
+ const FT f = plib::reciprocal(m_A[i][i]);
const auto &nzrd = this->m_terms[i].m_nzrd;
const auto &nzbd = this->m_terms[i].m_nzbd;
@@ -101,13 +101,13 @@ namespace solver
if (maxrow != i)
{
/* Swap the maxrow and ith row */
- for (std::size_t k = 0; k < kN + 1; k++) {
+ for (std::size_t k = 0; k < kN; k++) {
std::swap(m_A[i][k], m_A[maxrow][k]);
}
- //std::swap(RHS(i), RHS(maxrow));
+ std::swap(m_RHS[i], m_RHS[maxrow]);
}
/* FIXME: Singular matrix? */
- const FT f = 1.0 / m_A[i][i];
+ const FT f = plib::reciprocal(m_A[i][i]);
/* Eliminate column i from row j */
@@ -168,9 +168,9 @@ namespace solver
this->LE_solve();
this->LE_back_subst(m_new_V);
- const FT err = (newton_raphson ? this->delta(m_new_V) : 0.0);
+ const FT err = (newton_raphson ? this->delta(m_new_V) : plib::constants<FT>::zero());
this->store(m_new_V);
- return (err > this->m_params.m_accuracy) ? 2 : 1;
+ return (err > static_cast<FT>(this->m_params.m_accuracy)) ? 2 : 1;
}
template <typename FT, int SIZE>
diff --git a/src/lib/netlist/solver/nld_ms_direct1.h b/src/lib/netlist/solver/nld_ms_direct1.h
index 2931aaec29c..fa4d50677fb 100644
--- a/src/lib/netlist/solver/nld_ms_direct1.h
+++ b/src/lib/netlist/solver/nld_ms_direct1.h
@@ -39,9 +39,9 @@ namespace solver
std::array<FT, 1> new_V = { this->m_RHS[0] / this->m_A[0][0] };
- const FT err = (newton_raphson ? this->delta(new_V) : 0.0);
+ const FT err = (newton_raphson ? this->delta(new_V) : plib::constants<FT>::zero());
this->store(new_V);
- return (err > this->m_params.m_accuracy) ? 2 : 1;
+ return (err > static_cast<FT>(this->m_params.m_accuracy)) ? 2 : 1;
}
};
diff --git a/src/lib/netlist/solver/nld_ms_direct2.h b/src/lib/netlist/solver/nld_ms_direct2.h
index 03050083158..a24ceb9ef71 100644
--- a/src/lib/netlist/solver/nld_ms_direct2.h
+++ b/src/lib/netlist/solver/nld_ms_direct2.h
@@ -47,9 +47,9 @@ namespace solver
std::array<float_type, 2> new_V = {v0, v1};
this->m_stat_calculations++;
- const float_type err = (newton_raphson ? this->delta(new_V) : 0.0);
+ const float_type err = (newton_raphson ? this->delta(new_V) : plib::constants<FT>::zero());
this->store(new_V);
- return (err > this->m_params.m_accuracy) ? 2 : 1;
+ return (err > static_cast<float_type>(this->m_params.m_accuracy)) ? 2 : 1;
}
};
diff --git a/src/lib/netlist/solver/nld_ms_gcr.h b/src/lib/netlist/solver/nld_ms_gcr.h
index 66b56397a3d..d7f9716a5d7 100644
--- a/src/lib/netlist/solver/nld_ms_gcr.h
+++ b/src/lib/netlist/solver/nld_ms_gcr.h
@@ -227,7 +227,7 @@ namespace solver
unsigned matrix_solver_GCR_t<FT, SIZE>::vsolve_non_dynamic(const bool newton_raphson)
{
/* populate matrix */
- mat.set_scalar(0.0);
+ mat.set_scalar(plib::constants<FT>::zero());
this->fill_matrix(RHS);
/* now solve it */
@@ -247,9 +247,9 @@ namespace solver
this->m_stat_calculations++;
- const FT err = (newton_raphson ? this->delta(new_V) : 0.0);
+ const FT err = (newton_raphson ? this->delta(new_V) : plib::constants<FT>::zero());
this->store(new_V);
- return (err > this->m_params.m_accuracy) ? 2 : 1;
+ return (err > static_cast<FT>(this->m_params.m_accuracy)) ? 2 : 1;
}
} // namespace solver
diff --git a/src/lib/netlist/solver/nld_ms_gmres.h b/src/lib/netlist/solver/nld_ms_gmres.h
index 9a3c4e2f1bc..6d9e7929fe3 100644
--- a/src/lib/netlist/solver/nld_ms_gmres.h
+++ b/src/lib/netlist/solver/nld_ms_gmres.h
@@ -103,17 +103,17 @@ namespace solver
plib::parray<FT, SIZE> RHS(iN);
- m_ops.m_mat.set_scalar(0.0);
+ m_ops.m_mat.set_scalar(plib::constants<FT>::zero());
/* populate matrix and V for first estimate */
this->fill_matrix(RHS);
for (std::size_t k = 0; k < iN; k++)
{
- this->m_new_V[k] = this->m_terms[k].getV();
+ this->m_new_V[k] = this->m_terms[k].template getV<float_type>();
}
- const float_type accuracy = this->m_params.m_accuracy;
+ const float_type accuracy = static_cast<float_type>(this->m_params.m_accuracy);
auto iter = std::max(plib::constants<std::size_t>::one(), this->m_params.m_gs_loops());
auto gsl = m_gmres.solve(m_ops, this->m_new_V, RHS, iter, accuracy);
@@ -127,9 +127,9 @@ namespace solver
return matrix_solver_direct_t<FT, SIZE>::vsolve_non_dynamic(newton_raphson);
}
- const float_type err = (newton_raphson ? this->delta(this->m_new_V) : 0.0);
+ const float_type err = (newton_raphson ? this->delta(this->m_new_V) : plib::constants<FT>::zero());
this->store(this->m_new_V);
- return (err > this->m_params.m_accuracy) ? 2 : 1;
+ return (err > static_cast<float_type>(this->m_params.m_accuracy)) ? 2 : 1;
}
diff --git a/src/lib/netlist/solver/nld_ms_sm.h b/src/lib/netlist/solver/nld_ms_sm.h
index 4828dbfe72a..679819a2627 100644
--- a/src/lib/netlist/solver/nld_ms_sm.h
+++ b/src/lib/netlist/solver/nld_ms_sm.h
@@ -124,15 +124,15 @@ namespace solver
for (std::size_t j = 0; j < kN; j++)
{
W(i,j) = lA(i,j) = A(i,j);
- Ainv(i,j) = 0.0;
+ Ainv(i,j) = plib::constants<FT>::zero();
}
- Ainv(i,i) = 1.0;
+ Ainv(i,i) = plib::constants<FT>::one();
}
/* down */
for (std::size_t i = 0; i < kN; i++)
{
/* FIXME: Singular matrix? */
- const float_type f = 1.0 / W(i,i);
+ const float_type f = plib::reciprocal(W(i,i));
const auto * const p = this->m_terms[i].m_nzrd.data();
const std::size_t e = this->m_terms[i].m_nzrd.size();
@@ -144,7 +144,8 @@ namespace solver
{
const unsigned j = pb[jb];
const float_type f1 = - W(j,i) * f;
- if (f1 != 0.0)
+ // FIXME: comparison to zero
+ if (f1 != plib::constants<float_type>::zero())
{
for (std::size_t k = 0; k < e; k++)
W(j,p[k]) += W(i,p[k]) * f1;
@@ -157,11 +158,12 @@ namespace solver
for (std::size_t i = kN; i-- > 0; )
{
/* FIXME: Singular matrix? */
- const float_type f = 1.0 / W(i,i);
+ const float_type f = plib::reciprocal(W(i,i));
for (std::size_t j = i; j-- > 0; )
{
const float_type f1 = - W(j,i) * f;
- if (f1 != 0.0)
+ // FIXME: comparison to zero
+ if (f1 != plib::constants<float_type>::zero())
{
for (std::size_t k = i; k < kN; k++)
W(j,k) += W(i,k) * f1;
@@ -185,7 +187,7 @@ namespace solver
const std::size_t kN = this->size();
for (std::size_t i=0; i<kN; i++)
- x[i] = 0.0;
+ x[i] = plib::constants<FT>::zero();
for (std::size_t k=0; k<kN; k++)
{
@@ -234,13 +236,14 @@ namespace solver
v[col] = A(row,col) - lA(row,col);
if (incremental)
lA(row,col) = A(row,col);
- if (v[col] != 0.0)
+ // FIXME: comparison to zero
+ if (v[col] != plib::constants<float_type>::zero())
cols[colcount++] = col;
}
if (colcount > 0)
{
- float_type lamba = 0.0;
+ auto lamba(plib::constants<FT>::zero());
std::array<float_type, m_pitch> w = {0};
/* compute w and lamba */
@@ -258,11 +261,12 @@ namespace solver
w[k] += Ainv(col,k) * f; /* Transpose(Ainv) * v */
}
- lamba = -1.0 / (1.0 + lamba);
+ lamba = -plib::reciprocal(plib::constants<float_type>::one() + lamba);
for (std::size_t i=0; i<iN; i++)
{
const float_type f = lamba * z[i];
- if (f != 0.0)
+ // FIXME: comparison to zero
+ if (f != plib::constants<float_type>::zero())
for (std::size_t k = 0; k < iN; k++)
Ainv(i,k) += f * w[k];
}
@@ -275,9 +279,9 @@ namespace solver
this->LE_compute_x(new_V);
- const float_type err = (newton_raphson ? this->delta(new_V) : 0.0);
+ const float_type err = (newton_raphson ? this->delta(new_V) : plib::constants<FT>::zero());
this->store(new_V);
- return (err > this->m_params.m_accuracy) ? 2 : 1;
+ return (err > static_cast<FT>(this->m_params.m_accuracy)) ? 2 : 1;
}
template <typename FT, int SIZE>
diff --git a/src/lib/netlist/solver/nld_ms_sor.h b/src/lib/netlist/solver/nld_ms_sor.h
index 8bf4bf48721..08d8a358a45 100644
--- a/src/lib/netlist/solver/nld_ms_sor.h
+++ b/src/lib/netlist/solver/nld_ms_sor.h
@@ -34,9 +34,9 @@ namespace solver
const solver_parameters_t *params, const std::size_t size)
: matrix_solver_direct_t<FT, SIZE>(anetlist, name, nets, params, size)
, m_lp_fact(*this, "m_lp_fact", 0)
- , w(size, 0.0)
- , one_m_w(size, 0.0)
- , RHS(size, 0.0)
+ , w(size, plib::constants<FT>::zero())
+ , one_m_w(size, plib::constants<FT>::zero())
+ , RHS(size, plib::constants<FT>::zero())
//, new_V(size, 0.0)
{
}
@@ -70,13 +70,13 @@ namespace solver
* omega = 2.0 / (1.0 + std::sqrt(1-rho))
*/
- const float_type ws = this->m_params.m_gs_sor;
+ const auto ws(static_cast<float_type>(this->m_params.m_gs_sor));
for (std::size_t k = 0; k < iN; k++)
{
- float_type gtot_t = 0.0;
- float_type gabs_t = 0.0;
- float_type RHS_t = 0.0;
+ nl_fptype gtot_t = 0.0;
+ nl_fptype gabs_t = 0.0;
+ nl_fptype RHS_t = 0.0;
const std::size_t term_count = this->m_terms[k].count();
const nl_fptype * const gt = this->m_gtn[k];
@@ -84,7 +84,7 @@ namespace solver
const nl_fptype * const Idr = this->m_Idrn[k];
auto other_cur_analog = this->m_connected_net_Vn[k];
- this->m_new_V[k] = this->m_terms[k].getV();
+ this->m_new_V[k] = this->m_terms[k].template getV<float_type>();
for (std::size_t i = 0; i < term_count; i++)
{
@@ -95,7 +95,7 @@ namespace solver
for (std::size_t i = this->m_terms[k].railstart(); i < term_count; i++)
RHS_t = RHS_t - go[i] * *other_cur_analog[i];
- RHS[k] = RHS_t;
+ RHS[k] = static_cast<float_type>(RHS_t);
if (this->m_params.m_use_gabs)
{
@@ -105,23 +105,23 @@ namespace solver
gabs_t *= plib::constants<nl_fptype>::cast(0.5); // derived by try and error
if (gabs_t <= gtot_t)
{
- w[k] = ws / gtot_t;
+ w[k] = ws / static_cast<float_type>(gtot_t);
one_m_w[k] = plib::constants<FT>::one() - ws;
}
else
{
- w[k] = plib::constants<FT>::one() / (gtot_t + gabs_t);
- one_m_w[k] = plib::constants<FT>::one() - plib::constants<FT>::one() * gtot_t / (gtot_t + gabs_t);
+ w[k] = plib::reciprocal(static_cast<float_type>(gtot_t + gabs_t));
+ one_m_w[k] = plib::constants<FT>::one() - plib::constants<FT>::one() * static_cast<FT>(gtot_t / (gtot_t + gabs_t));
}
}
else
{
- w[k] = ws / gtot_t;
+ w[k] = ws / static_cast<float_type>(gtot_t);
one_m_w[k] = plib::constants<FT>::one() - ws;
}
}
- const float_type accuracy = this->m_params.m_accuracy;
+ const auto accuracy(static_cast<float_type>(this->m_params.m_accuracy));
do {
resched = false;
@@ -132,9 +132,9 @@ namespace solver
const std::size_t railstart = this->m_terms[k].railstart();
const nl_fptype * go = this->m_gonn[k];
- float_type Idrive = 0.0;
+ float_type Idrive = plib::constants<float_type>::zero();
for (std::size_t i = 0; i < railstart; i++)
- Idrive = Idrive - go[i] * this->m_new_V[static_cast<std::size_t>(net_other[i])];
+ Idrive = Idrive - static_cast<float_type>(go[i]) * this->m_new_V[static_cast<std::size_t>(net_other[i])];
const float_type new_val = this->m_new_V[k] * one_m_w[k] + (Idrive + RHS[k]) * w[k];
@@ -158,9 +158,9 @@ namespace solver
return matrix_solver_direct_t<FT, SIZE>::vsolve_non_dynamic(newton_raphson);
}
- const float_type err = (newton_raphson ? this->delta(this->m_new_V) : 0.0);
+ const float_type err = (newton_raphson ? this->delta(this->m_new_V) : plib::constants<FT>::zero());
this->store(this->m_new_V);
- return (err > this->m_params.m_accuracy) ? 2 : 1;
+ return (err > static_cast<float_type>(this->m_params.m_accuracy)) ? 2 : 1;
}
} // namespace solver
diff --git a/src/lib/netlist/solver/nld_ms_sor_mat.h b/src/lib/netlist/solver/nld_ms_sor_mat.h
index 8a1575ebfc7..da862e182ec 100644
--- a/src/lib/netlist/solver/nld_ms_sor_mat.h
+++ b/src/lib/netlist/solver/nld_ms_sor_mat.h
@@ -37,7 +37,7 @@ namespace solver
const solver_parameters_t *params, std::size_t size)
: matrix_solver_direct_t<FT, SIZE>(anetlist, name, nets, params, size)
, m_Vdelta(*this, "m_Vdelta", std::vector<float_type>(size))
- , m_omega(*this, "m_omega", params->m_gs_sor)
+ , m_omega(*this, "m_omega", static_cast<float_type>(params->m_gs_sor))
, m_lp_fact(*this, "m_lp_fact", 0)
{
}
@@ -84,7 +84,7 @@ namespace solver
{
float_type sq = 0;
float_type sqo = 0;
- const float_type rez_cts = 1.0 / this->current_timestep();
+ const float_type rez_cts = plib::reciprocal(this->current_timestep());
for (unsigned k = 0; k < this->size(); k++)
{
const analog_net_t *n = this->m_nets[k];
@@ -158,11 +158,11 @@ namespace solver
#endif
for (std::size_t k = 0; k < iN; k++)
- this->m_new_V[k] = this->m_terms[k].getV();
+ this->m_new_V[k] = this->m_terms[k].template getV<FT>();
do {
resched = false;
- float_type cerr = 0.0;
+ FT cerr = plib::constants<FT>::zero();
for (std::size_t k = 0; k < iN; k++)
{
@@ -177,7 +177,7 @@ namespace solver
FT w = m_omega / this->m_A[k][k];
if (this->m_params.m_use_gabs)
{
- FT gabs_t = 0.0;
+ 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]]);
@@ -194,7 +194,7 @@ namespace solver
this->m_new_V[k] += delta;
}
- if (cerr > this->m_params.m_accuracy)
+ if (cerr > static_cast<float_type>(this->m_params.m_accuracy))
{
resched = true;
}
@@ -210,9 +210,9 @@ namespace solver
return matrix_solver_direct_t<FT, SIZE>::solve_non_dynamic(newton_raphson);
}
- const float_type err = (newton_raphson ? this->delta(this->m_new_V) : 0.0);
+ const float_type err = (newton_raphson ? this->delta(this->m_new_V) : plib::constants<FT>::zero());
this->store(this->m_new_V);
- return (err > this->m_params.m_accuracy) ? 2 : 1;
+ return (err > static_cast<float_type>(this->m_params.m_accuracy)) ? 2 : 1;
}
diff --git a/src/lib/netlist/solver/nld_ms_w.h b/src/lib/netlist/solver/nld_ms_w.h
index 3470cd86f22..b89edd70b74 100644
--- a/src/lib/netlist/solver/nld_ms_w.h
+++ b/src/lib/netlist/solver/nld_ms_w.h
@@ -134,15 +134,15 @@ namespace solver
for (std::size_t j = 0; j < kN; j++)
{
W(i,j) = lA(i,j) = A(i,j);
- Ainv(i,j) = 0.0;
+ Ainv(i,j) = plib::constants<FT>::zero();
}
- Ainv(i,i) = 1.0;
+ Ainv(i,i) = plib::constants<FT>::one();
}
/* down */
for (std::size_t i = 0; i < kN; i++)
{
/* FIXME: Singular matrix? */
- const float_type f = 1.0 / W(i,i);
+ const float_type f = plib::reciprocal(W(i,i));
const auto * const p = this->m_terms[i].m_nzrd.data();
const size_t e = this->m_terms[i].m_nzrd.size();
@@ -154,7 +154,8 @@ namespace solver
{
const auto j = pb[jb];
const float_type f1 = - W(j,i) * f;
- if (f1 != 0.0)
+ // FIXME: comparison to zero
+ if (f1 != plib::constants<float_type>::zero())
{
for (std::size_t k = 0; k < e; k++)
W(j,p[k]) += W(i,p[k]) * f1;
@@ -167,11 +168,12 @@ namespace solver
for (std::size_t i = kN; i-- > 0; )
{
/* FIXME: Singular matrix? */
- const float_type f = 1.0 / W(i,i);
+ const float_type f = plib::reciprocal(W(i,i));
for (std::size_t j = i; j-- > 0; )
{
const float_type f1 = - W(j,i) * f;
- if (f1 != 0.0)
+ // FIXME: comparison to zero
+ if (f1 != plib::constants<float_type>::zero())
{
for (std::size_t k = i; k < kN; k++)
W(j,k) += W(i,k) * f1;
@@ -194,7 +196,7 @@ namespace solver
const std::size_t kN = this->size();
for (std::size_t i=0; i<kN; i++)
- x[i] = 0.0;
+ x[i] = plib::constants<FT>::zero();
for (std::size_t k=0; k<kN; k++)
{
@@ -257,7 +259,7 @@ namespace solver
for (unsigned i = 0; i < rowcount; i++)
{
const unsigned r = rows[i];
- nl_fptype tmp = 0.0;
+ FT tmp = plib::constants<FT>::zero();
for (unsigned k = 0; k < iN; k++)
tmp += VT(r,k) * new_V[k];
w[i] = tmp;
@@ -265,17 +267,18 @@ namespace solver
for (unsigned i = 0; i < rowcount; i++)
for (unsigned k=0; k< rowcount; k++)
- H[i][k] = 0.0;
+ H[i][k] = plib::constants<FT>::zero();
for (unsigned i = 0; i < rowcount; i++)
- H[i][i] = 1.0;
+ H[i][i] = plib::constants<FT>::one();
/* Construct H = (I + VT*Z) */
for (unsigned i = 0; i < rowcount; i++)
for (unsigned k=0; k< colcount[i]; k++)
{
const unsigned col = cols[i][k];
float_type f = VT(rows[i],col);
- if (f!=0.0)
+ // FIXME: comparison to zero
+ if (f != plib::constants<float_type>::zero())
for (unsigned j= 0; j < rowcount; j++)
H[i][j] += f * Ainv(col,rows[j]);
}
@@ -283,14 +286,16 @@ namespace solver
/* Gaussian elimination of H */
for (unsigned i = 0; i < rowcount; i++)
{
- if (H[i][i] == 0.0)
+ // FIXME: comparison to zero
+ if (H[i][i] == plib::constants<float_type>::zero())
plib::perrlogger("{} H singular\n", this->name());
- const float_type f = 1.0 / H[i][i];
+ const float_type f = plib::reciprocal(H[i][i]);
for (unsigned j = i+1; j < rowcount; j++)
{
const float_type f1 = - f * H[j][i];
- if (f1!=0.0)
+ // FIXME: comparison to zero
+ if (f1 != plib::constants<float_type>::zero())
{
float_type *pj = &H[j][i+1];
const float_type *pi = &H[i][i+1];
@@ -317,7 +322,7 @@ namespace solver
/* x = y - Zt */
for (unsigned i=0; i<iN; i++)
{
- float_type tmp = 0.0;
+ float_type tmp = plib::constants<FT>::zero();
for (unsigned j=0; j<rowcount;j++)
{
const unsigned row = rows[j];
@@ -332,18 +337,18 @@ namespace solver
if (false)
for (unsigned i=0; i<iN; i++)
{
- float_type tmp = 0.0;
+ float_type tmp = plib::constants<FT>::zero();
for (unsigned j=0; j<iN; j++)
{
tmp += A(i,j) * new_V[j];
}
- if (std::abs(tmp-RHS(i)) > 1e-6)
+ 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));
}
- const float_type err = (newton_raphson ? this->delta(new_V) : 0.0);
+ const float_type err = (newton_raphson ? this->delta(new_V) : plib::constants<FT>::zero());
this->store(new_V);
- return (err > this->m_params.m_accuracy) ? 2 : 1;
+ return (err > static_cast<float_type>(this->m_params.m_accuracy)) ? 2 : 1;
}
template <typename FT, int SIZE>
diff --git a/src/lib/netlist/solver/nld_solver.cpp b/src/lib/netlist/solver/nld_solver.cpp
index 506bff88dab..86772fac6a0 100644
--- a/src/lib/netlist/solver/nld_solver.cpp
+++ b/src/lib/netlist/solver/nld_solver.cpp
@@ -77,7 +77,7 @@ namespace devices
netlist_time now(exec().time());
/* force solving during start up if there are no time-step devices */
/* FIXME: Needs a more elegant solution */
- bool force_solve = (now < netlist_time::from_double(2 * m_params.m_max_timestep));
+ bool force_solve = (now < netlist_time::from_fp<decltype(m_params.m_max_timestep)>(2 * m_params.m_max_timestep));
std::size_t nthreads = std::min(static_cast<std::size_t>(m_params.m_parallel()), plib::omp::get_max_threads());
@@ -105,7 +105,7 @@ namespace devices
/* step circuit */
if (!m_Q_step.net().is_queued())
{
- m_Q_step.net().toggle_and_push_to_queue(netlist_time::from_double(m_params.m_max_timestep));
+ m_Q_step.net().toggle_and_push_to_queue(netlist_time::from_fp(m_params.m_max_timestep));
}
}
@@ -151,6 +151,112 @@ namespace devices
return plib::unique_ptr<solver::matrix_solver_t>();
}
+ template <typename FT>
+ plib::unique_ptr<solver::matrix_solver_t> NETLIB_NAME(solver)::create_solvers(
+ const pstring &sname,
+ analog_net_t::list_t &nets)
+ {
+ std::size_t net_count = nets.size();
+ switch (net_count)
+ {
+ case 1:
+ return plib::make_unique<solver::matrix_solver_direct1_t<FT>>(state(), sname, nets, &m_params);
+ break;
+ case 2:
+ return plib::make_unique<solver::matrix_solver_direct2_t<FT>>(state(), sname, nets, &m_params);
+ break;
+#if 0
+ case 3:
+ return create_solver<FT, 3>(3, sname, nets);
+ break;
+ case 4:
+ return create_solver<FT, 4>(4, sname, nets);
+ break;
+ case 5:
+ return create_solver<FT, 5>(5, sname, nets);
+ break;
+ case 6:
+ return create_solver<FT, 6>(6, sname, nets);
+ break;
+ case 7:
+ return create_solver<FT, 7>(7, sname, nets);
+ break;
+ case 8:
+ return create_solver<FT, 8>(8, sname, nets);
+ break;
+ case 9:
+ return create_solver<FT, 9>(9, sname, nets);
+ break;
+ case 10:
+ return create_solver<FT, 10>(10, sname, nets);
+ break;
+#if 0
+ case 11:
+ return create_solver<FT, 11>(11, sname);
+ break;
+ case 12:
+ return create_solver<FT, 12>(12, sname);
+ break;
+ case 15:
+ return create_solver<FT, 15>(15, sname);
+ break;
+ case 31:
+ return create_solver<FT, 31>(31, sname);
+ break;
+ case 35:
+ return create_solver<FT, 35>(35, sname);
+ break;
+ case 43:
+ return create_solver<FT, 43>(43, sname);
+ break;
+ case 49:
+ return create_solver<FT, 49>(49, sname);
+ break;
+#endif
+#if 1
+ case 87:
+ return create_solver<FT,86>(86, sname, nets);
+ break;
+#endif
+#endif
+ default:
+ log().info(MI_NO_SPECIFIC_SOLVER(net_count));
+ if (net_count <= 8)
+ {
+ return create_solver<FT, -8>(net_count, sname, nets);
+ }
+ else if (net_count <= 16)
+ {
+ return create_solver<FT, -16>(net_count, sname, nets);
+ }
+ else if (net_count <= 32)
+ {
+ return create_solver<FT, -32>(net_count, sname, nets);
+ }
+ else if (net_count <= 64)
+ {
+ return create_solver<FT, -64>(net_count, sname, nets);
+ }
+ else if (net_count <= 128)
+ {
+ return create_solver<FT, -128>(net_count, sname, nets);
+ }
+ else if (net_count <= 256)
+ {
+ return create_solver<FT, -256>(net_count, sname, nets);
+ }
+ else if (net_count <= 512)
+ {
+ return create_solver<FT, -512>(net_count, sname, nets);
+ }
+ else
+ {
+ return create_solver<FT, 0>(net_count, sname, nets);
+ }
+ break;
+ }
+ }
+
struct net_splitter
{
@@ -212,7 +318,6 @@ namespace devices
void NETLIB_NAME(solver)::post_start()
{
-
log().verbose("Scanning net groups ...");
// determine net groups
@@ -225,105 +330,18 @@ namespace devices
for (auto & grp : splitter.groups)
{
plib::unique_ptr<solver::matrix_solver_t> ms;
- std::size_t net_count = grp.size();
pstring sname = plib::pfmt("Solver_{1}")(m_mat_solvers.size());
- switch (net_count)
+ switch (m_params.m_fp_type())
{
- case 1:
- ms = plib::make_unique<solver::matrix_solver_direct1_t<nl_mat_fptype>>(state(), sname, grp, &m_params);
- break;
- case 2:
- ms = plib::make_unique<solver::matrix_solver_direct2_t<nl_mat_fptype>>(state(), sname, grp, &m_params);
- break;
-#if 0
- case 3:
- ms = create_solver<nl_mat_fptype, 3>(3, sname, grp);
- break;
- case 4:
- ms = create_solver<nl_mat_fptype, 4>(4, sname, grp);
- break;
- case 5:
- ms = create_solver<nl_mat_fptype, 5>(5, sname, grp);
- break;
- case 6:
- ms = create_solver<nl_mat_fptype, 6>(6, sname, grp);
- break;
- case 7:
- ms = create_solver<nl_mat_fptype, 7>(7, sname, grp);
- break;
- case 8:
- ms = create_solver<nl_mat_fptype, 8>(8, sname, grp);
- break;
- case 9:
- ms = create_solver<nl_mat_fptype, 9>(9, sname, grp);
- break;
- case 10:
- ms = create_solver<nl_mat_fptype, 10>(10, sname, grp);
+ case solver::matrix_fp_type_e::FLOAT:
+ ms = create_solvers<float>(sname, grp);
break;
- #if 0
- case 11:
- ms = create_solver<nl_mat_fptype, 11>(11, sname);
+ case solver::matrix_fp_type_e::DOUBLE:
+ ms = create_solvers<double>(sname, grp);
break;
- case 12:
- ms = create_solver<nl_mat_fptype, 12>(12, sname);
- break;
- case 15:
- ms = create_solver<nl_mat_fptype, 15>(15, sname);
- break;
- case 31:
- ms = create_solver<nl_mat_fptype, 31>(31, sname);
- break;
- case 35:
- ms = create_solver<nl_mat_fptype, 35>(35, sname);
- break;
- case 43:
- ms = create_solver<nl_mat_fptype, 43>(43, sname);
- break;
- case 49:
- ms = create_solver<nl_mat_fptype, 49>(49, sname);
- break;
- #endif
- #if 1
- case 87:
- ms = create_solver<nl_mat_fptype,86>(86, sname, grp);
- break;
- #endif
- #endif
- default:
- log().info(MI_NO_SPECIFIC_SOLVER(net_count));
- if (net_count <= 8)
- {
- ms = create_solver<nl_mat_fptype, -8>(net_count, sname, grp);
- }
- else if (net_count <= 16)
- {
- ms = create_solver<nl_mat_fptype, -16>(net_count, sname, grp);
- }
- else if (net_count <= 32)
- {
- ms = create_solver<nl_mat_fptype, -32>(net_count, sname, grp);
- }
- else if (net_count <= 64)
- {
- ms = create_solver<nl_mat_fptype, -64>(net_count, sname, grp);
- }
- else if (net_count <= 128)
- {
- ms = create_solver<nl_mat_fptype, -128>(net_count, sname, grp);
- }
- else if (net_count <= 256)
- {
- ms = create_solver<nl_mat_fptype, -256>(net_count, sname, grp);
- }
- else if (net_count <= 512)
- {
- ms = create_solver<nl_mat_fptype, -512>(net_count, sname, grp);
- }
- else
- {
- ms = create_solver<nl_mat_fptype, 0>(net_count, sname, grp);
- }
+ case solver::matrix_fp_type_e::LONGDOUBLE:
+ ms = create_solvers<long double>(sname, grp);
break;
}
diff --git a/src/lib/netlist/solver/nld_solver.h b/src/lib/netlist/solver/nld_solver.h
index dba6b69f358..02b231655ee 100644
--- a/src/lib/netlist/solver/nld_solver.h
+++ b/src/lib/netlist/solver/nld_solver.h
@@ -60,6 +60,10 @@ namespace devices
template <typename FT, int SIZE>
plib::unique_ptr<solver::matrix_solver_t> create_solver(std::size_t size,
const pstring &solvername, analog_net_t::list_t &nets);
+
+ template <typename FT>
+ plib::unique_ptr<solver::matrix_solver_t> create_solvers(
+ const pstring &sname, analog_net_t::list_t &nets);
};
} // namespace devices