netlist: code maintenance (nw)

- clang lint and pedantic fixes
- mat_cr.h: separate solving linear systems from underlying matrix
This commit is contained in:
couriersud 2019-10-12 19:36:50 +02:00
parent df53d0ef17
commit 8f83e4392f
26 changed files with 605 additions and 499 deletions

View File

@ -73,6 +73,7 @@ project "netlist"
MAME_DIR .. "src/lib/netlist/plib/pstate.h",
MAME_DIR .. "src/lib/netlist/plib/pstring.cpp",
MAME_DIR .. "src/lib/netlist/plib/pstring.h",
MAME_DIR .. "src/lib/netlist/plib/pstrutil.h",
MAME_DIR .. "src/lib/netlist/plib/pstream.h",
MAME_DIR .. "src/lib/netlist/plib/ptime.h",
MAME_DIR .. "src/lib/netlist/plib/ptypes.h",

View File

@ -476,7 +476,7 @@ namespace analog
const nl_double Vctrl = (is_forward ? Vgs : Vgd) - Vth;
nl_double Ids, gm, gds, gmb;
nl_double Ids(0), gm(0), gds(0), gmb(0);
const nl_double absVds = std::abs(Vds);
if (Vctrl <= 0.0)

View File

@ -24,6 +24,8 @@ TIDY_FLAGSX += -cppcoreguidelines-macro-usage,
TIDY_FLAGSX += -cppcoreguidelines-non-private-member-variables-in-classes,-misc-non-private-member-variables-in-classes,
TIDY_FLAGSX += -bugprone-macro-parentheses,-misc-macro-parentheses,
TIDY_FLAGSX += -modernize-use-trailing-return-type
TIDY_FLAGSX += -bugprone-too-small-loop-variable
TIDY_FLAGSX += -cppcoreguidelines-pro-bounds-array-to-pointer-decay
space :=
space +=

View File

@ -182,7 +182,7 @@ namespace netlist
NETLIB_UPDATE(74123)
{
netlist_sig_t m_trig;
netlist_sig_t m_trig(0);
netlist_sig_t res = !m_CLRQ();
netlist_time t_AB_to_Q = NLTIME_FROM_NS(10);
netlist_time t_C_to_Q = NLTIME_FROM_NS(10);

View File

@ -175,8 +175,6 @@ namespace netlist
{
{
// recompute
nl_double freq;
nl_double v_freq_2, v_freq_3, v_freq_4;
nl_double v_freq = m_FC();
nl_double v_rng = m_RNG();
@ -196,10 +194,10 @@ namespace netlist
/* Polyfunctional3D_model created by zunzun.com using sum of squared absolute error */
v_freq_2 = v_freq * v_freq;
v_freq_3 = v_freq_2 * v_freq;
v_freq_4 = v_freq_3 * v_freq;
freq = k1;
nl_double v_freq_2 = v_freq * v_freq;
nl_double v_freq_3 = v_freq_2 * v_freq;
nl_double v_freq_4 = v_freq_3 * v_freq;
nl_double freq = k1;
freq += k2 * v_freq;
freq += k3 * v_freq_2;
freq += k4 * v_freq_3;

View File

@ -416,7 +416,7 @@ nl_double parser_t::eval_param(const token_t &tok)
static std::array<pstring, 6> macs = {"", "RES_K", "RES_M", "CAP_U", "CAP_N", "CAP_P"};
static std::array<nl_double, 6> facs = {1, 1e3, 1e6, 1e-6, 1e-9, 1e-12};
std::size_t f=0;
nl_double ret;
nl_double ret(0);
for (std::size_t i=1; i<macs.size();i++)
if (tok.str() == macs[i])
@ -429,7 +429,7 @@ nl_double parser_t::eval_param(const token_t &tok)
}
else
{
bool err;
bool err(false);
ret = plib::pstonum_ne<nl_double, true>(tok.str(), err);
if (err)
error(plib::pfmt("Parameter value <{1}> not double \n")(tok.str()));

View File

@ -35,10 +35,11 @@ namespace plib
template <typename FT, int SIZE>
struct mat_precondition_ILU
{
using mat_type = plib::matrix_compressed_rows_t<FT, SIZE>;
using mat_type = plib::pmatrix_cr_t<FT, SIZE>;
using matLU_type = plib::pLUmatrix_cr_t<mat_type>;
mat_precondition_ILU(std::size_t size, std::size_t ilu_scale = 4
, std::size_t bw = plib::matrix_compressed_rows_t<FT, SIZE>::FILL_INFINITY)
, std::size_t bw = plib::pmatrix_cr_t<FT, SIZE>::FILL_INFINITY)
: m_mat(static_cast<typename mat_type::index_type>(size))
, m_LU(static_cast<typename mat_type::index_type>(size))
, m_ILU_scale(static_cast<std::size_t>(ilu_scale))
@ -50,9 +51,7 @@ namespace plib
void build(M &fill)
{
m_mat.build_from_fill_mat(fill, 0);
m_LU.gaussian_extend_fill_mat(fill);
m_LU.build_from_fill_mat(fill, m_ILU_scale, m_band_width); // ILU(2)
//m_LU.build_from_fill_mat(fill, 9999, 20); // Band matrix width 20
m_LU.build(fill, m_ILU_scale);
}
@ -64,11 +63,7 @@ namespace plib
void precondition()
{
if (m_ILU_scale < 1)
m_LU.raw_copy_from(m_mat);
else
m_LU.reduction_copy_from(m_mat);
m_LU.incomplete_LU_factorization();
m_LU.incomplete_LU_factorization(m_mat);
}
template<typename V>
@ -80,7 +75,7 @@ namespace plib
PALIGNAS_VECTOROPT()
mat_type m_mat;
PALIGNAS_VECTOROPT()
mat_type m_LU;
matLU_type m_LU;
std::size_t m_ILU_scale;
std::size_t m_band_width;
};
@ -174,7 +169,7 @@ namespace plib
v[i] = v[i] * m_diag[i];
}
plib::matrix_compressed_rows_t<FT, SIZE> m_mat;
plib::pmatrix_cr_t<FT, SIZE> m_mat;
plib::parray<FT, SIZE> m_diag;
plib::parray<std::vector<std::size_t>, SIZE > nzcol;
};
@ -210,7 +205,7 @@ namespace plib
plib::unused_var(v);
}
plib::matrix_compressed_rows_t<FT, SIZE> m_mat;
plib::pmatrix_cr_t<FT, SIZE> m_mat;
};
/* FIXME: hardcoding RESTART to 20 becomes an issue on very large

View File

@ -31,12 +31,15 @@ namespace plib
// template<typename T, int N, typename C = std::size_t>
template<typename T, int N, typename C = uint16_t>
struct matrix_compressed_rows_t
struct pmatrix_cr_t
{
using index_type = C;
using value_type = T;
COPYASSIGNMOVE(matrix_compressed_rows_t, default)
static constexpr const int NSQ = (N > 0 ? -N * N : N * N);
static constexpr const int Np1 = (N == 0) ? 0 : (N < 0 ? N - 1 : N + 1);
COPYASSIGNMOVE(pmatrix_cr_t, default)
enum constants_e
{
@ -44,34 +47,27 @@ namespace plib
};
parray<index_type, N> diag; // diagonal index pointer n
parray<index_type, (N == 0) ? 0 : (N < 0 ? N - 1 : N + 1)> row_idx; // row index pointer n + 1
parray<index_type, N < 0 ? -N * N : N *N> col_idx; // column index array nz_num, initially (n * n)
parray<value_type, N < 0 ? -N * N : N *N> A; // Matrix elements nz_num, initially (n * n)
//parray<C, N < 0 ? -N * (N-1) / 2 : N * (N+1) / 2 > nzbd; // Support for gaussian elimination
parray<std::vector<index_type>, N > nzbd; // Support for gaussian elimination
// contains elimination rows below the diagonal
// FIXME: convert to pvector
parray<index_type, (N == 0) ? 0 : (N < 0 ? N - 1 : N + 1)> ilu_rows;
std::vector<std::vector<index_type>> m_ge_par; // parallel execution support for Gauss
parray<index_type, Np1> row_idx; // row index pointer n + 1
parray<index_type, NSQ> col_idx; // column index array nz_num, initially (n * n)
parray<value_type, NSQ> A; // Matrix elements nz_num, initially (n * n)
index_type nz_num;
explicit matrix_compressed_rows_t(const index_type n)
explicit pmatrix_cr_t(const index_type n)
: diag(n)
, row_idx(n+1)
, col_idx(n*n)
, A(n*n)
, nz_num(0)
//, nzbd(n * (n+1) / 2)
, nzbd(n)
, ilu_rows(n+1)
, nz_num(0)
, m_size(n)
{
for (index_type i=0; i<n+1; i++)
row_idx[i] = 0;
}
~matrix_compressed_rows_t() = default;
~pmatrix_cr_t() = default;
constexpr index_type size() const { return static_cast<index_type>((N>0) ? N : m_size); }
@ -113,39 +109,6 @@ namespace plib
}
}
template <typename M>
std::pair<std::size_t, std::size_t> gaussian_extend_fill_mat(M &fill)
{
std::size_t ops = 0;
std::size_t fill_max = 0;
for (std::size_t k = 0; k < fill.size(); k++)
{
ops++; // 1/A(k,k)
for (std::size_t row = k + 1; row < fill.size(); row++)
{
if (fill[row][k] < FILL_INFINITY)
{
ops++;
for (std::size_t col = k + 1; col < fill[row].size(); col++)
//if (fill[k][col] < FILL_INFINITY)
{
auto f = std::min(fill[row][col], 1 + fill[row][k] + fill[k][col]);
if (f < FILL_INFINITY)
{
if (f > fill_max)
fill_max = f;
ops += 2;
}
fill[row][col] = f;
}
}
}
}
build_parallel_gaussian_execution_scheme(fill);
return { fill_max, ops };
}
template <typename M>
void build_from_fill_mat(const M &f, std::size_t max_fill = FILL_INFINITY - 1,
std::size_t band_width = FILL_INFINITY)
@ -180,142 +143,8 @@ namespace plib
nzbd[k].push_back(0); // end of sequence
}
/* build ilu_rows */
index_type p(0);
for (index_type k=1; k < size(); k++)
if (row_idx[k] < diag[k])
ilu_rows[p++] = k;
ilu_rows[p] = 0; // end of array
}
template <typename V>
void gaussian_elimination(V & RHS)
{
const std::size_t iN = size();
for (std::size_t i = 0; i < iN - 1; i++)
{
std::size_t nzbdp = 0;
std::size_t pi = diag[i];
const value_type f = 1.0 / A[pi++];
const std::size_t piie = row_idx[i+1];
const auto &nz = nzbd[i];
while (auto j = nz[nzbdp++]) // NOLINT(bugprone-infinite-loop)
{
// proceed to column i
std::size_t pj = row_idx[j];
std::size_t pje = row_idx[j+1];
while (col_idx[pj] < i)
pj++;
const value_type f1 = - A[pj++] * f;
// subtract row i from j
// fill-in available assumed, i.e. matrix was prepared
for (std::size_t pii = pi; pii<piie && pj < pje; pii++)
{
while (col_idx[pj] < col_idx[pii])
pj++;
if (col_idx[pj] == col_idx[pii])
A[pj++] += A[pii] * f1;
}
RHS[j] += f1 * RHS[i];
}
}
}
int get_parallel_level(std::size_t k) const
{
for (std::size_t i = 0; i < m_ge_par.size(); i++)
if (plib::container::contains( m_ge_par[i], k))
return static_cast<int>(i);
return -1;
}
template <typename V>
void gaussian_elimination_parallel(V & RHS)
{
//printf("omp: %ld %d %d\n", m_ge_par.size(), nz_num, (int)m_ge_par[m_ge_par.size()-2].size());
for (auto l = 0ul; l < m_ge_par.size(); l++)
plib::omp::for_static(nz_num, 0ul, m_ge_par[l].size(), [this, &RHS, &l] (unsigned ll)
{
auto &i = m_ge_par[l][ll];
{
std::size_t nzbdp = 0;
std::size_t pi = diag[i];
const value_type f = 1.0 / A[pi++];
const std::size_t piie = row_idx[i+1];
const auto &nz = nzbd[i];
while (auto j = nz[nzbdp++])
{
// proceed to column i
std::size_t pj = row_idx[j];
while (col_idx[pj] < i)
pj++;
const value_type f1 = - A[pj++] * f;
// subtract row i from j
// fill-in available assumed, i.e. matrix was prepared
for (std::size_t pii = pi; pii<piie; pii++)
{
while (col_idx[pj] < col_idx[pii])
pj++;
if (col_idx[pj] == col_idx[pii])
A[pj++] += A[pii] * f1;
}
RHS[j] += f1 * RHS[i];
}
}
});
}
template <typename V1, typename V2>
void gaussian_back_substitution(V1 &V, const V2 &RHS)
{
const std::size_t iN = size();
/* row n-1 */
V[iN - 1] = RHS[iN - 1] / A[diag[iN - 1]];
for (std::size_t j = iN - 1; j-- > 0;)
{
value_type tmp = 0;
const auto jdiag = diag[j];
const std::size_t e = row_idx[j+1];
for (std::size_t pk = jdiag + 1; pk < e; pk++)
tmp += A[pk] * V[col_idx[pk]];
V[j] = (RHS[j] - tmp) / A[jdiag];
}
}
template <typename V1>
void gaussian_back_substitution(V1 &V)
{
const std::size_t iN = size();
/* row n-1 */
V[iN - 1] = V[iN - 1] / A[diag[iN - 1]];
for (std::size_t j = iN - 1; j-- > 0;)
{
value_type tmp = 0;
const auto jdiag = diag[j];
const std::size_t e = row_idx[j+1];
for (std::size_t pk = jdiag + 1; pk < e; pk++)
tmp += A[pk] * V[col_idx[pk]];
V[j] = (V[j] - tmp) / A[jdiag];
}
}
template <typename VTV, typename VTR>
void mult_vec(VTR & res, const VTV & x)
{
@ -394,7 +223,7 @@ namespace plib
}
}
/* checks at all - may crash */
/* no checks at all - may crash */
template <typename LUMAT>
void raw_copy_from(LUMAT & src)
{
@ -402,7 +231,305 @@ namespace plib
A[k] = src.A[k];
}
void incomplete_LU_factorization()
protected:
parray<std::vector<index_type>, N > nzbd; // Support for gaussian elimination
private:
//parray<C, N < 0 ? -N * (N-1) / 2 : N * (N+1) / 2 > nzbd; // Support for gaussian elimination
index_type m_size;
};
template<typename B>
struct pGEmatrix_cr_t : public B
{
using base = B;
using index_type = typename base::index_type;
COPYASSIGNMOVE(pGEmatrix_cr_t, default)
explicit pGEmatrix_cr_t(const index_type n)
: B(n)
{
}
~pGEmatrix_cr_t() = default;
template <typename M>
std::pair<std::size_t, std::size_t> gaussian_extend_fill_mat(M &fill)
{
std::size_t ops = 0;
std::size_t fill_max = 0;
for (std::size_t k = 0; k < fill.size(); k++)
{
ops++; // 1/A(k,k)
for (std::size_t row = k + 1; row < fill.size(); row++)
{
if (fill[row][k] < base::FILL_INFINITY)
{
ops++;
for (std::size_t col = k + 1; col < fill[row].size(); col++)
//if (fill[k][col] < FILL_INFINITY)
{
auto f = std::min(fill[row][col], 1 + fill[row][k] + fill[k][col]);
if (f < base::FILL_INFINITY)
{
if (f > fill_max)
fill_max = f;
ops += 2;
}
fill[row][col] = f;
}
}
}
}
build_parallel_gaussian_execution_scheme(fill);
return { fill_max, ops };
}
template <typename V>
void gaussian_elimination(V & RHS)
{
const std::size_t iN = base::size();
for (std::size_t i = 0; i < iN - 1; i++)
{
std::size_t nzbdp = 0;
std::size_t pi = base::diag[i];
const typename base::value_type f = 1.0 / base::A[pi++];
const std::size_t piie = base::row_idx[i+1];
const auto &nz = base::nzbd[i];
while (auto j = nz[nzbdp++]) // NOLINT(bugprone-infinite-loop)
{
// proceed to column i
std::size_t pj = base::row_idx[j];
std::size_t pje = base::row_idx[j+1];
while (base::col_idx[pj] < i)
pj++;
const typename base::value_type f1 = - base::A[pj++] * f;
// subtract row i from j
// fill-in available assumed, i.e. matrix was prepared
for (std::size_t pii = pi; pii<piie && pj < pje; pii++)
{
while (base::col_idx[pj] < base::col_idx[pii])
pj++;
if (base::col_idx[pj] == base::col_idx[pii])
base::A[pj++] += base::A[pii] * f1;
}
RHS[j] += f1 * RHS[i];
}
}
}
int get_parallel_level(std::size_t k) const
{
for (std::size_t i = 0; i < m_ge_par.size(); i++)
if (plib::container::contains( m_ge_par[i], k))
return static_cast<int>(i);
return -1;
}
template <typename V>
void gaussian_elimination_parallel(V & RHS)
{
//printf("omp: %ld %d %d\n", m_ge_par.size(), nz_num, (int)m_ge_par[m_ge_par.size()-2].size());
for (auto l = 0ul; l < m_ge_par.size(); l++)
plib::omp::for_static(base::nz_num, 0ul, m_ge_par[l].size(), [this, &RHS, &l] (unsigned ll)
{
auto &i = m_ge_par[l][ll];
{
std::size_t nzbdp = 0;
std::size_t pi = base::diag[i];
const typename base::value_type f = 1.0 / base::A[pi++];
const std::size_t piie = base::row_idx[i+1];
const auto &nz = base::nzbd[i];
while (auto j = nz[nzbdp++])
{
// proceed to column i
std::size_t pj = base::row_idx[j];
while (base::col_idx[pj] < i)
pj++;
const typename base::value_type f1 = - base::A[pj++] * f;
// subtract row i from j
// fill-in available assumed, i.e. matrix was prepared
for (std::size_t pii = pi; pii<piie; pii++)
{
while (base::col_idx[pj] < base::col_idx[pii])
pj++;
if (base::col_idx[pj] == base::col_idx[pii])
base::A[pj++] += base::A[pii] * f1;
}
RHS[j] += f1 * RHS[i];
}
}
});
}
template <typename V1, typename V2>
void gaussian_back_substitution(V1 &V, const V2 &RHS)
{
const std::size_t iN = base::size();
/* row n-1 */
V[iN - 1] = RHS[iN - 1] / base::A[base::diag[iN - 1]];
for (std::size_t j = iN - 1; j-- > 0;)
{
typename base::value_type tmp = 0;
const auto jdiag = base::diag[j];
const std::size_t e = base::row_idx[j+1];
for (std::size_t pk = jdiag + 1; pk < e; pk++)
tmp += base::A[pk] * V[base::col_idx[pk]];
V[j] = (RHS[j] - tmp) / base::A[jdiag];
}
}
template <typename V1>
void gaussian_back_substitution(V1 &V)
{
const std::size_t iN = base::size();
/* row n-1 */
V[iN - 1] = V[iN - 1] / base::A[base::diag[iN - 1]];
for (std::size_t j = iN - 1; j-- > 0;)
{
typename base::value_type tmp = 0;
const auto jdiag = base::diag[j];
const std::size_t e = base::row_idx[j+1];
for (std::size_t pk = jdiag + 1; pk < e; pk++)
tmp += base::A[pk] * V[base::col_idx[pk]];
V[j] = (V[j] - tmp) / base::A[jdiag];
}
}
private:
template <typename M>
void build_parallel_gaussian_execution_scheme(const M &fill)
{
// calculate parallel scheme for gaussian elimination
std::vector<std::vector<index_type>> rt(base::size());
for (index_type k = 0; k < base::size(); k++)
{
for (index_type j = k+1; j < base::size(); j++)
{
if (fill[j][k] < base::FILL_INFINITY)
{
rt[k].push_back(j);
}
}
}
std::vector<index_type> levGE(base::size(), 0);
index_type cl = 0;
for (index_type k = 0; k < base::size(); k++ )
{
if (levGE[k] >= cl)
{
std::vector<index_type> t = rt[k];
for (index_type j = k+1; j < base::size(); j++ )
{
bool overlap = false;
// is there overlap
if (plib::container::contains(t, j))
overlap = true;
for (auto &x : rt[j])
if (plib::container::contains(t, x))
{
overlap = true;
break;
}
if (overlap)
levGE[j] = cl + 1;
else
{
t.push_back(j);
for (auto &x : rt[j])
t.push_back(x);
}
}
cl++;
}
}
m_ge_par.clear();
m_ge_par.resize(cl+1);
for (index_type k = 0; k < base::size(); k++)
m_ge_par[levGE[k]].push_back(k);
//for (std::size_t k = 0; k < m_ge_par.size(); k++)
// printf("%d %d\n", (int) k, (int) m_ge_par[k].size());
}
// contains elimination rows below the diagonal
std::vector<std::vector<index_type>> m_ge_par; // parallel execution support for Gauss
};
template<typename B>
struct pLUmatrix_cr_t : public B
{
using base = B;
using index_type = typename base::index_type;
COPYASSIGNMOVE(pLUmatrix_cr_t, default)
explicit pLUmatrix_cr_t(const index_type n)
: B(n)
, ilu_rows(n+1)
, m_ILUp(0)
{
}
~pLUmatrix_cr_t() = default;
template <typename M>
void build(M &fill, std::size_t ilup)
{
index_type p(0);
/* build ilu_rows */
for (index_type i=1; i < fill.size(); i++)
{
bool found(false);
for (index_type k = 0; k < i; k++)
{
// if (fill[i][k] < base::FILL_INFINITY)
if (fill[i][k] <= ilup)
{
// assume A[k][k]!=0
for (index_type j=k+1; j < fill.size(); j++)
{
auto f = std::min(fill[i][j], 1 + fill[i][k] + fill[k][j]);
//if (f < base::FILL_INFINITY)
if (f <= ilup)
{
#if 0
if (f > fill_max)
fill_max = f;
ops += 2;
#endif
fill[i][j] = f;
}
}
found = true;
}
}
if (found)
ilu_rows[p++] = i;
}
ilu_rows[p] = 0; // end of array
this->build_from_fill_mat(fill, ilup); //, m_band_width); // ILU(2)
m_ILUp = ilup;
}
void incomplete_LU_factorization(const base &mat)
{
/*
* incomplete LU Factorization according to http://de.wikipedia.org/wiki/ILU-Zerlegung
@ -421,31 +548,35 @@ namespace plib
* i=i+1
*
*/
if (m_ILUp < 1)
this->raw_copy_from(mat);
else
this->reduction_copy_from(mat);
index_type p(0);
while (auto i = ilu_rows[p++]) // NOLINT(bugprone-infinite-loop)
{
const auto p_i_end = row_idx[i + 1];
const auto p_i_end = base::row_idx[i + 1];
// loop over all columns k left of diag in row i
//if (row_idx[i] < diag[i])
// printf("occ %d\n", (int)i);
for (auto i_k = row_idx[i]; i_k < diag[i]; i_k++)
for (auto i_k = base::row_idx[i]; i_k < base::diag[i]; i_k++)
{
const auto k(col_idx[i_k]);
const auto p_k_end(row_idx[k + 1]);
const T LUp_i_k = A[i_k] = A[i_k] / A[diag[k]];
const auto k(base::col_idx[i_k]);
const auto p_k_end(base::row_idx[k + 1]);
const typename base::value_type LUp_i_k = base::A[i_k] = base::A[i_k] / base::A[base::diag[k]];
auto k_j(diag[k] + 1);
auto i_j(i_k + 1);
index_type k_j(base::diag[k] + 1);
index_type i_j(i_k + 1);
while (i_j < p_i_end && k_j < p_k_end ) // pj = (i, j)
{
// we can assume that within a row ja increases continuously */
const auto c_i_j(col_idx[i_j]); // row i, column j
const auto c_k_j(col_idx[k_j]); // row k, column j
const index_type c_i_j(base::col_idx[i_j]); // row i, column j
const auto c_k_j(base::col_idx[k_j]); // row k, column j
if (c_k_j == c_i_j)
A[i_j] -= LUp_i_k * A[k_j];
base::A[i_j] -= LUp_i_k * base::A[k_j];
k_j += (c_k_j <= c_i_j ? 1 : 0);
i_j += (c_k_j >= c_i_j ? 1 : 0);
@ -478,86 +609,30 @@ namespace plib
* This can be solved for x using backwards elimination in U.
*
*/
for (index_type i = 1; i < size(); ++i )
for (index_type i = 1; i < base::size(); ++i )
{
T tmp = 0.0;
const auto j1(row_idx[i]);
const auto j2(diag[i]);
typename base::value_type tmp(0);
const auto j1(base::row_idx[i]);
const auto j2(base::diag[i]);
for (auto j = j1; j < j2; ++j )
tmp += A[j] * r[col_idx[j]];
tmp += base::A[j] * r[base::col_idx[j]];
r[i] -= tmp;
}
// i now is equal to n;
for (auto i = size(); i-- > 0; )
for (index_type i = base::size(); i-- > 0; )
{
T tmp = 0.0;
const auto di(diag[i]);
const auto j2(row_idx[i+1]);
for (auto j = di + 1; j < j2; j++ )
tmp += A[j] * r[col_idx[j]];
r[i] = (r[i] - tmp) / A[di];
typename base::value_type tmp(0);
const auto di(base::diag[i]);
const auto j2(base::row_idx[i+1]);
for (index_type j = di + 1; j < j2; j++ )
tmp += base::A[j] * r[base::col_idx[j]];
r[i] = (r[i] - tmp) / base::A[di];
}
}
private:
template <typename M>
void build_parallel_gaussian_execution_scheme(const M &fill)
{
// calculate parallel scheme for gaussian elimination
std::vector<std::vector<index_type>> rt(size());
for (index_type k = 0; k < size(); k++)
{
for (index_type j = k+1; j < size(); j++)
{
if (fill[j][k] < FILL_INFINITY)
{
rt[k].push_back(j);
}
}
}
std::vector<index_type> levGE(size(), 0);
index_type cl = 0;
for (index_type k = 0; k < size(); k++ )
{
if (levGE[k] >= cl)
{
std::vector<index_type> t = rt[k];
for (index_type j = k+1; j < size(); j++ )
{
bool overlap = false;
// is there overlap
if (plib::container::contains(t, j))
overlap = true;
for (auto &x : rt[j])
if (plib::container::contains(t, x))
{
overlap = true;
break;
}
if (overlap)
levGE[j] = cl + 1;
else
{
t.push_back(j);
for (auto &x : rt[j])
t.push_back(x);
}
}
cl++;
}
}
m_ge_par.clear();
m_ge_par.resize(cl+1);
for (index_type k = 0; k < size(); k++)
m_ge_par[levGE[k]].push_back(k);
//for (std::size_t k = 0; k < m_ge_par.size(); k++)
// printf("%d %d\n", (int) k, (int) m_ge_par[k].size());
}
index_type m_size;
parray<index_type, base::Np1> ilu_rows;
std::size_t m_ILUp;
};
} // namespace plib

View File

@ -7,6 +7,7 @@
#include "pfmtlog.h"
#include "palloc.h"
#include "pstrutil.h"
#include <algorithm>
#include <array>

View File

@ -8,6 +8,7 @@
#include "pfunction.h"
#include "pexception.h"
#include "pfmtlog.h"
#include "pstrutil.h"
#include "putil.h"
#include <cmath>
@ -71,7 +72,7 @@ void pfunction::compile_postfix(const std::vector<pstring> &inputs,
if (rc.m_cmd != PUSH_INPUT)
{
rc.m_cmd = PUSH_CONST;
bool err;
bool err(false);
rc.m_param = plib::pstonum_ne<decltype(rc.m_param), true>(cmd, err);
if (err)
throw plib::pexception(plib::pfmt("pfunction: unknown/misformatted token <{1}> in <{2}>")(cmd)(expr));

View File

@ -7,6 +7,7 @@
#include "poptions.h"
#include "pexception.h"
#include "pstrutil.h"
#include "ptypes.h"
namespace plib {

View File

@ -184,7 +184,7 @@ public:
protected:
int parse(const pstring &argument) override
{
bool err;
bool err(false);
m_val = pstonum_ne<T, true>(argument, err);
return (err ? 1 : (m_val < m_min || m_val > m_max));
}

View File

@ -120,7 +120,7 @@ double ptokenizer::get_number_double()
{
error(pfmt("Expected a number, got <{1}>")(tok.str()) );
}
bool err;
bool err(false);
auto ret = plib::pstonum_ne<double, true>(tok.str(), err);
if (err)
error(pfmt("Expected a number, got <{1}>")(tok.str()) );
@ -134,7 +134,7 @@ long ptokenizer::get_number_long()
{
error(pfmt("Expected a long int, got <{1}>")(tok.str()) );
}
bool err;
bool err(false);
auto ret = plib::pstonum_ne<long, true>(tok.str(), err);
if (err)
error(pfmt("Expected a long int, got <{1}>")(tok.str()) );

View File

@ -13,6 +13,7 @@
#include "pexception.h"
#include "pfmtlog.h"
#include "pstring.h"
#include "pstrutil.h"
#include <array>
#include <fstream>

View File

@ -492,193 +492,6 @@ using putf8string = pstring_t<putf8_traits>;
using pu16string = pstring_t<putf16_traits>;
using pwstring = pstring_t<pwchar_traits>;
namespace plib
{
template<class T>
struct string_info
{
using mem_t = typename T::mem_t;
};
template<>
struct string_info<std::string>
{
using mem_t = char;
};
template<typename T>
pstring to_string(const T &v)
{
return pstring(std::to_string(v));
}
template<typename T>
pwstring to_wstring(const T &v)
{
return pwstring(std::to_wstring(v));
}
template<typename T>
typename T::size_type find_first_not_of(const T &str, const T &no)
{
typename T::size_type pos = 0;
for (auto it = str.begin(); it != str.end(); ++it, ++pos)
{
bool f = true;
for (typename T::value_type const jt : no)
{
if (*it == jt)
{
f = false;
break;
}
}
if (f)
return pos;
}
return T::npos;
}
template<typename T>
typename T::size_type find_last_not_of(const T &str, const T &no)
{
/* FIXME: reverse iterator */
typename T::size_type last_found = T::npos;
typename T::size_type pos = 0;
for (auto it = str.begin(); it != str.end(); ++it, ++pos)
{
bool f = true;
for (typename T::value_type const jt : no)
{
if (*it == jt)
{
f = false;
break;
}
}
if (f)
last_found = pos;
}
return last_found;
}
template<typename T>
T ltrim(const T &str, const T &ws = T(" \t\n\r"))
{
auto f = find_first_not_of(str, ws);
return (f == T::npos) ? T() : str.substr(f);
}
template<typename T>
T rtrim(const T &str, const T &ws = T(" \t\n\r"))
{
auto f = find_last_not_of(str, ws);
return (f == T::npos) ? T() : str.substr(0, f + 1);
}
template<typename T>
T trim(const T &str, const T &ws = T(" \t\n\r"))
{
return rtrim(ltrim(str, ws), ws);
}
template<typename T>
T left(const T &str, typename T::size_type len)
{
return str.substr(0, len);
}
template<typename T>
T right(const T &str, typename T::size_type nlen)
{
return nlen >= str.length() ? str : str.substr(str.length() - nlen, nlen);
}
template<typename T>
bool startsWith(const T &str, const T &arg)
{
return (arg == left(str, arg.length()));
}
template<typename T>
bool endsWith(const T &str, const T &arg)
{
return (right(str, arg.length()) == arg);
}
template<typename T, typename TA>
bool startsWith(const T &str, const TA &arg)
{
return startsWith(str, static_cast<pstring>(arg));
}
template<typename T, typename TA>
bool endsWith(const T &str, const TA &arg)
{
return endsWith(str, static_cast<pstring>(arg));
}
template<typename T>
std::size_t strlen(const T *str)
{
const T *p = str;
while (*p)
p++;
return static_cast<std::size_t>(p - str);
}
template<typename T>
T ucase(const T &str)
{
T ret;
for (const auto &c : str)
if (c >= 'a' && c <= 'z')
ret += (c - 'a' + 'A');
else
ret += c;
return ret;
}
template<typename T>
T rpad(const T &str, const T &ws, const typename T::size_type cnt)
{
// FIXME: pstringbuffer ret(*this);
T ret(str);
typename T::size_type wsl = ws.length();
for (auto i = ret.length(); i < cnt; i+=wsl)
ret += ws;
return ret;
}
template<typename T>
T replace_all(const T &str, const T &search, const T &replace)
{
T ret;
const typename T::size_type slen = search.length();
typename T::size_type last_s = 0;
typename T::size_type s = str.find(search, last_s);
while (s != T::npos)
{
ret += str.substr(last_s, s - last_s);
ret += replace;
last_s = s + slen;
s = str.find(search, last_s);
}
ret += str.substr(last_s);
return ret;
}
template<typename T, typename T1, typename T2>
T replace_all(const T &str, const T1 &search, const T2 &replace)
{
return replace_all(str, static_cast<T>(search), static_cast<T>(replace));
}
} // namespace plib
// custom specialization of std::hash can be injected in namespace std
namespace std
{

View File

@ -0,0 +1,206 @@
// license:GPL-2.0+
// copyright-holders:Couriersud
/*
* pstrutil.h
*/
#ifndef PSTRUTIL_H_
#define PSTRUTIL_H_
#include "pstring.h"
#include "ptypes.h"
#include <exception>
#include <iterator>
#include <limits>
#include <stdexcept>
#include <string>
#include <type_traits>
namespace plib
{
template<class T>
struct string_info
{
using mem_t = typename T::mem_t;
};
template<>
struct string_info<std::string>
{
using mem_t = char;
};
template<typename T>
pstring to_string(const T &v)
{
return pstring(std::to_string(v));
}
template<typename T>
pwstring to_wstring(const T &v)
{
return pwstring(std::to_wstring(v));
}
template<typename T>
typename T::size_type find_first_not_of(const T &str, const T &no)
{
typename T::size_type pos = 0;
for (auto it = str.begin(); it != str.end(); ++it, ++pos)
{
bool f = true;
for (typename T::value_type const jt : no)
{
if (*it == jt)
{
f = false;
break;
}
}
if (f)
return pos;
}
return T::npos;
}
template<typename T>
typename T::size_type find_last_not_of(const T &str, const T &no)
{
/* FIXME: reverse iterator */
typename T::size_type last_found = T::npos;
typename T::size_type pos = 0;
for (auto it = str.begin(); it != str.end(); ++it, ++pos)
{
bool f = true;
for (typename T::value_type const jt : no)
{
if (*it == jt)
{
f = false;
break;
}
}
if (f)
last_found = pos;
}
return last_found;
}
template<typename T>
T ltrim(const T &str, const T &ws = T(" \t\n\r"))
{
auto f = find_first_not_of(str, ws);
return (f == T::npos) ? T() : str.substr(f);
}
template<typename T>
T rtrim(const T &str, const T &ws = T(" \t\n\r"))
{
auto f = find_last_not_of(str, ws);
return (f == T::npos) ? T() : str.substr(0, f + 1);
}
template<typename T>
T trim(const T &str, const T &ws = T(" \t\n\r"))
{
return rtrim(ltrim(str, ws), ws);
}
template<typename T>
T left(const T &str, typename T::size_type len)
{
return str.substr(0, len);
}
template<typename T>
T right(const T &str, typename T::size_type nlen)
{
return nlen >= str.length() ? str : str.substr(str.length() - nlen, nlen);
}
template<typename T>
bool startsWith(const T &str, const T &arg)
{
return (arg == left(str, arg.length()));
}
template<typename T>
bool endsWith(const T &str, const T &arg)
{
return (right(str, arg.length()) == arg);
}
template<typename T, typename TA>
bool startsWith(const T &str, const TA &arg)
{
return startsWith(str, static_cast<pstring>(arg));
}
template<typename T, typename TA>
bool endsWith(const T &str, const TA &arg)
{
return endsWith(str, static_cast<pstring>(arg));
}
template<typename T>
std::size_t strlen(const T *str)
{
const T *p = str;
while (*p)
p++;
return static_cast<std::size_t>(p - str);
}
template<typename T>
T ucase(const T &str)
{
T ret;
for (const auto &c : str)
if (c >= 'a' && c <= 'z')
ret += (c - 'a' + 'A');
else
ret += c;
return ret;
}
template<typename T>
T rpad(const T &str, const T &ws, const typename T::size_type cnt)
{
// FIXME: pstringbuffer ret(*this);
T ret(str);
typename T::size_type wsl = ws.length();
for (auto i = ret.length(); i < cnt; i+=wsl)
ret += ws;
return ret;
}
template<typename T>
T replace_all(const T &str, const T &search, const T &replace)
{
T ret;
const typename T::size_type slen = search.length();
typename T::size_type last_s = 0;
typename T::size_type s = str.find(search, last_s);
while (s != T::npos)
{
ret += str.substr(last_s, s - last_s);
ret += replace;
last_s = s + slen;
s = str.find(search, last_s);
}
ret += str.substr(last_s);
return ret;
}
template<typename T, typename T1, typename T2>
T replace_all(const T &str, const T1 &search, const T2 &replace)
{
return replace_all(str, static_cast<T>(search), static_cast<T>(replace));
}
} // namespace plib
#endif /* PSTRUTIL_H_ */

View File

@ -3,11 +3,10 @@
#include "putil.h"
#include "plists.h"
#include "pstrutil.h"
#include "ptypes.h"
#include <algorithm>
//#include <cstdlib>
//#include <cstring>
#include <initializer_list>
namespace plib

View File

@ -291,7 +291,7 @@ struct input_t
: m_value(0.0)
{
std::array<char, 400> buf; // NOLINT(cppcoreguidelines-pro-type-member-init)
double t;
double t(0);
// NOLINTNEXTLINE(cppcoreguidelines-pro-type-vararg)
int e = std::sscanf(line.c_str(), "%lf,%[^,],%lf", &t, buf.data(), &m_value);
if (e != 3)

View File

@ -418,7 +418,7 @@ namespace devices
++m_stat_vsolver_calls;
if (has_dynamic_devices())
{
std::size_t this_resched;
std::size_t this_resched(0);
std::size_t newton_loops = 0;
do
{
@ -559,10 +559,14 @@ namespace devices
if (m_params.m_dynamic_ts)
{
#if 0
for (std::size_t k = 0, iN=m_terms.size(); k < iN; k++)
{
terms_for_net_t *t = m_terms[k].get();
#else
for (auto &t : m_terms)
{
#endif
//const nl_double DD_n = (n->Q_Analog() - t->m_last_V);
// avoid floating point exceptions
const nl_double DD_n = std::max(-1e100, std::min(1e100,(t->getV() - t->m_last_V)));
@ -570,7 +574,7 @@ namespace devices
//printf("%g %g %g %g\n", DD_n, hn, t->m_DD_n_m_1, t->m_h_n_m_1);
nl_double DD2 = (DD_n / hn - t->m_DD_n_m_1 / t->m_h_n_m_1) / (hn + t->m_h_n_m_1);
nl_double new_net_timestep;
nl_double new_net_timestep(0);
t->m_h_n_m_1 = hn;
t->m_DD_n_m_1 = DD_n;

View File

@ -10,6 +10,7 @@
#include "netlist/nl_base.h"
#include "netlist/nl_errstr.h"
#include "plib/mat_cr.h"
#include "plib/palloc.h"
#include "plib/pmatrix2d.h"
#include "plib/putil.h"
@ -307,12 +308,14 @@ namespace devices
template <typename T, typename M>
void log_fill(const T &fill, M &mat)
{
/* FIXME: move this to the cr matrix class and use computed
* parallel ordering once it makes sense.
*/
const std::size_t iN = fill.size();
// FIXME: Not yet working, mat_cr.h needs some more work
#if 0
auto mat_GE = dynamic_cast<plib::pGEmatrix_cr_t<typename M::base> *>(&mat);
#else
plib::unused_var(mat);
#endif
std::vector<unsigned> levL(iN, 0);
std::vector<unsigned> levU(iN, 0);
@ -346,7 +349,13 @@ namespace devices
if (fill[k][j] > fm)
fm = fill[k][j];
}
this->log().verbose("{1:4} {2} {3:4} {4:4} {5:4} {6:4}", k, ml, levL[k], levU[k], mat.get_parallel_level(k), fm);
#if 0
this->log().verbose("{1:4} {2} {3:4} {4:4} {5:4} {6:4}", k, ml,
levL[k], levU[k], mat_GE ? mat_GE->get_parallel_level(k) : 0, fm);
#else
this->log().verbose("{1:4} {2} {3:4} {4:4} {5:4} {6:4}", k, ml,
levL[k], levU[k], 0, fm);
#endif
}
}

View File

@ -10,8 +10,8 @@
#include "nld_matrix_solver.h"
#include "nld_solver.h"
#include "plib/vector_ops.h"
#include "plib/parray.h"
#include "plib/vector_ops.h"
#include <algorithm>
#include <cmath>
@ -95,10 +95,10 @@ namespace devices
const auto &nzrd = m_terms[i]->m_nzrd;
const auto &nzbd = m_terms[i]->m_nzbd;
for (std::size_t j : nzbd)
for (const std::size_t j : nzbd)
{
const FT f1 = -f * A(j, i);
for (std::size_t k : nzrd)
for (const std::size_t k : nzrd)
A(j, k) += A(i, k) * f1;
//RHS(j) += RHS(i) * f1;
}

View File

@ -30,7 +30,7 @@ namespace devices
{
public:
using mat_type = plib::matrix_compressed_rows_t<FT, SIZE>;
using mat_type = plib::pGEmatrix_cr_t<plib::pmatrix_cr_t<FT, SIZE>>;
// FIXME: dirty hack to make this compile
static constexpr const std::size_t storage_N = 100;
@ -54,7 +54,7 @@ namespace devices
private:
using mat_index_type = typename plib::matrix_compressed_rows_t<FT, SIZE>::index_type;
using mat_index_type = typename plib::pmatrix_cr_t<FT, SIZE>::index_type;
void generate_code(plib::putf8_fmt_writer &strm);

View File

@ -48,7 +48,7 @@ namespace devices
private:
using mattype = typename plib::matrix_compressed_rows_t<FT, SIZE>::index_type;
using mattype = typename plib::pmatrix_cr_t<FT, SIZE>::index_type;
//plib::mat_precondition_none<FT, SIZE> m_ops;
plib::mat_precondition_ILU<FT, SIZE> m_ops;

View File

@ -227,7 +227,7 @@ namespace devices
switch (net_count)
{
#if 1
#if 0
case 1:
ms = plib::make_unique<matrix_solver_direct1_t<double>>(state(), sname, &m_params);
break;

View File

@ -337,7 +337,7 @@ void nl_convert_spice_t::process_line(const pstring &line)
*/
pstring model;
pstring pins ="CBE";
bool err;
bool err(false);
auto nval = plib::pstonum_ne<long, true>(tt[4], err);
plib::unused_var(nval);

View File

@ -7,7 +7,7 @@
#endif
#ifndef NLTOOL_VERSION
#define USE_FRONTIERS 1
#define USE_FRONTIERS 0
#else
#define USE_FRONTIERS 0
#endif
@ -361,7 +361,7 @@ NETLIST_END()
NETLIST_START(kidniki)
#if (0 || USE_FRONTIERS)
#if (1 || USE_FRONTIERS)
SOLVER(Solver, 48000)
PARAM(Solver.ACCURACY, 1e-7)
PARAM(Solver.NR_LOOPS, 300)