Code maintenance. In addition, added a parameter to select iterative

solver. (nw)
This commit is contained in:
couriersud 2015-06-24 19:13:20 +02:00
parent e5cd96af89
commit 5d49ee6fbd
8 changed files with 97 additions and 188 deletions

View File

@ -1021,7 +1021,7 @@ ATTR_COLD analog_output_t::analog_output_t()
this->set_net(m_my_net);
set_state(STATE_OUT);
net().as_analog().m_cur_Analog = 0.98;
net().as_analog().m_cur_Analog = NL_FCONST(0.99);
}
ATTR_COLD void analog_output_t::init_object(core_device_t &dev, const pstring &aname)
@ -1047,35 +1047,6 @@ ATTR_COLD param_t::param_t(const param_type_t atype)
{
}
ATTR_COLD param_double_t::param_double_t()
: param_t(DOUBLE)
, m_param(0.0)
{
}
ATTR_COLD param_int_t::param_int_t()
: param_t(INTEGER)
, m_param(0)
{
}
ATTR_COLD param_logic_t::param_logic_t()
: param_int_t()
{
}
ATTR_COLD param_str_t::param_str_t()
: param_t(STRING)
, m_param("")
{
}
ATTR_COLD param_model_t::param_model_t()
: param_t(MODEL)
, m_param("")
{
}
ATTR_COLD const pstring param_model_t::model_type() const
{
pstring tmp = this->Value();

View File

@ -914,88 +914,57 @@ namespace netlist
const param_type_t m_param_type;
};
class param_double_t : public param_t
template <class C, param_t::param_type_t T>
class param_template_t : public param_t
{
NETLIST_PREVENT_COPYING(param_double_t)
NETLIST_PREVENT_COPYING(param_template_t)
public:
ATTR_COLD param_double_t();
ATTR_COLD param_template_t()
: param_t(T)
, m_param(C(0))
{
}
ATTR_HOT void setTo(const nl_double param);
ATTR_COLD void initial(const nl_double val) { m_param = val; }
ATTR_HOT nl_double Value() const { return m_param; }
operator const C() const { return Value(); }
ATTR_HOT void setTo(const C &param);
ATTR_COLD void initial(const C &val) { m_param = val; }
ATTR_HOT C Value() const { return m_param; }
protected:
virtual void save_register()
{
save(NLNAME(m_param));
/* pstrings not yet supported, these need special logic */
if (T != param_t::STRING && T != param_t::MODEL)
save(NLNAME(m_param));
param_t::save_register();
}
C m_param;
private:
nl_double m_param;
};
class param_int_t : public param_t
{
NETLIST_PREVENT_COPYING(param_int_t)
public:
ATTR_COLD param_int_t();
ATTR_HOT void setTo(const int param);
ATTR_COLD void initial(const int val) { m_param = val; }
ATTR_HOT int Value() const { return m_param; }
protected:
virtual void save_register()
{
save(NLNAME(m_param));
param_t::save_register();
}
private:
int m_param;
};
typedef param_template_t<nl_double, param_t::DOUBLE> param_double_t;
typedef param_template_t<int, param_t::INTEGER> param_int_t;
typedef param_template_t<pstring, param_t::STRING> param_str_t;
class param_logic_t : public param_int_t
{
NETLIST_PREVENT_COPYING(param_logic_t)
public:
ATTR_COLD param_logic_t();
ATTR_COLD param_logic_t() : param_int_t() { };
};
class param_str_t : public param_t
{
NETLIST_PREVENT_COPYING(param_str_t)
public:
ATTR_COLD param_str_t();
ATTR_HOT void setTo(const pstring &param);
ATTR_COLD void initial(const pstring &val) { m_param = val; }
ATTR_HOT const pstring &Value() const { return m_param; }
private:
pstring m_param;
};
class param_model_t : public param_t
class param_model_t : public param_template_t<pstring, param_t::MODEL>
{
NETLIST_PREVENT_COPYING(param_model_t)
public:
ATTR_COLD param_model_t();
ATTR_COLD void initial(const pstring &val) { m_param = val; }
ATTR_HOT const pstring &Value() const { return m_param; }
ATTR_COLD param_model_t() : param_template_t<pstring, param_t::MODEL>() { }
/* these should be cached! */
ATTR_COLD nl_double model_value(const pstring &entity, const nl_double defval = 0.0) const;
ATTR_COLD const pstring model_value_str(const pstring &entity, const pstring defval = "") const;
ATTR_COLD const pstring model_type() const;
private:
pstring m_param;
};
// -----------------------------------------------------------------------------
@ -1296,22 +1265,8 @@ namespace netlist
PSTATE_INTERFACE(object_t, m_netlist, name())
ATTR_HOT inline void param_str_t::setTo(const pstring &param)
{
m_param = param;
netdev().update_param();
}
ATTR_HOT inline void param_int_t::setTo(const int param)
{
if (m_param != param)
{
m_param = param;
netdev().update_param();
}
}
ATTR_HOT inline void param_double_t::setTo(const nl_double param)
template <class C, param_t::param_type_t T>
ATTR_HOT inline void param_template_t<C, T>::setTo(const C &param)
{
if (m_param != param)
{

View File

@ -37,7 +37,7 @@ enum pstate_data_type_e {
DT_INT8,
DT_INT,
DT_BOOLEAN,
DT_FLOAT
DT_FLOAT,
};
template<typename _ItemType> struct nl_datatype

View File

@ -31,8 +31,6 @@ public:
, m_use_iLU_preconditioning(true)
, m_use_more_precise_stop_condition(false)
, m_accuracy_mult(1.0)
, m_gs_fail(0)
, m_gs_total(0)
{
int mr=this->N(); /* FIXME: maximum iterations locked in here */
@ -55,8 +53,6 @@ public:
delete[] m_v[i];
}
virtual void log_stats();
virtual void vsetup(analog_net_t::list_t &nets);
ATTR_HOT virtual int vsolve_non_dynamic(const bool newton_raphson);
protected:
@ -83,42 +79,17 @@ private:
double * RESTRICT m_v[_storage_N + 1]; /*(mr + 1), n */
//double m_y[_storage_N]; /* mr + 1 */
double m_accuracy_mult;
int m_gs_fail;
int m_gs_total;
double m_accuracy_mult; // FXIME: Save state
};
// ----------------------------------------------------------------------------------------
// matrix_solver - Gauss - Seidel
// matrix_solver - GMRES
// ----------------------------------------------------------------------------------------
template <unsigned m_N, unsigned _storage_N>
void matrix_solver_GMRES_t<m_N, _storage_N>::log_stats()
{
if (this->m_stat_calculations != 0 && this->m_params.m_log_stats)
{
this->netlist().log("==============================================");
this->netlist().log("Solver %s", this->name().cstr());
this->netlist().log(" ==> %d nets", this->N()); //, (*(*groups[i].first())->m_core_terms.first())->name().cstr());
this->netlist().log(" has %s elements", this->is_dynamic() ? "dynamic" : "no dynamic");
this->netlist().log(" has %s elements", this->is_timestep() ? "timestep" : "no timestep");
this->netlist().log(" %6.3f average newton raphson loops", (double) this->m_stat_newton_raphson / (double) this->m_stat_vsolver_calls);
this->netlist().log(" %10d invocations (%6d Hz) %10d gs fails (%6.2f%%) %6.3f average",
this->m_stat_calculations,
this->m_stat_calculations * 10 / (int) (this->netlist().time().as_double() * 10.0),
this->m_gs_fail,
100.0 * (double) this->m_gs_fail / (double) this->m_stat_calculations,
(double) this->m_gs_total / (double) this->m_stat_calculations);
}
}
template <unsigned m_N, unsigned _storage_N>
void matrix_solver_GMRES_t<m_N, _storage_N>::vsetup(analog_net_t::list_t &nets)
{
matrix_solver_direct_t<m_N, _storage_N>::vsetup(nets);
this->save(NLNAME(m_gs_fail));
this->save(NLNAME(m_gs_total));
int nz = 0;
const int iN = this->N();
@ -220,6 +191,7 @@ ATTR_HOT inline int matrix_solver_GMRES_t<m_N, _storage_N>::vsolve_non_dynamic(c
const nl_double accuracy = this->m_params.m_accuracy;
#if 1
int mr = std::min(iN-1,(int) sqrt(iN));
mr = std::min(mr, this->m_params.m_gs_loops);
int iter = 4;
int gsl = solve_ilu_gmres(new_V, RHS, iter, mr, accuracy);
int failed = mr * iter;
@ -228,7 +200,7 @@ ATTR_HOT inline int matrix_solver_GMRES_t<m_N, _storage_N>::vsolve_non_dynamic(c
//int gsl = tt_ilu_cr(new_V, RHS, failed, accuracy);
int gsl = tt_gs_cr(new_V, RHS, failed, accuracy);
#endif
m_gs_total += gsl;
this->m_iterative_total += gsl;
this->m_stat_calculations++;
if (gsl>=failed)
@ -236,7 +208,7 @@ ATTR_HOT inline int matrix_solver_GMRES_t<m_N, _storage_N>::vsolve_non_dynamic(c
//for (int k = 0; k < iN; k++)
// this->m_nets[k]->m_cur_Analog = new_V[k];
// Fallback to direct solver ...
this->m_gs_fail++;
this->m_iterative_fail++;
return matrix_solver_direct_t<m_N, _storage_N>::vsolve_non_dynamic(newton_raphson);
}

View File

@ -27,15 +27,11 @@ public:
matrix_solver_SOR_t(const solver_parameters_t *params, int size)
: matrix_solver_direct_t<m_N, _storage_N>(matrix_solver_t::GAUSS_SEIDEL, params, size)
, m_lp_fact(0)
, m_gs_fail(0)
, m_gs_total(0)
{
}
virtual ~matrix_solver_SOR_t() {}
virtual void log_stats();
virtual void vsetup(analog_net_t::list_t &nets);
ATTR_HOT virtual int vsolve_non_dynamic(const bool newton_raphson);
protected:
@ -43,41 +39,18 @@ protected:
private:
nl_double m_lp_fact;
int m_gs_fail;
int m_gs_total;
};
// ----------------------------------------------------------------------------------------
// matrix_solver - Gauss - Seidel
// ----------------------------------------------------------------------------------------
template <unsigned m_N, unsigned _storage_N>
void matrix_solver_SOR_t<m_N, _storage_N>::log_stats()
{
if (this->m_stat_calculations != 0 && this->m_params.m_log_stats)
{
this->netlist().log("==============================================");
this->netlist().log("Solver %s", this->name().cstr());
this->netlist().log(" ==> %d nets", this->N()); //, (*(*groups[i].first())->m_core_terms.first())->name().cstr());
this->netlist().log(" has %s elements", this->is_dynamic() ? "dynamic" : "no dynamic");
this->netlist().log(" has %s elements", this->is_timestep() ? "timestep" : "no timestep");
this->netlist().log(" %6.3f average newton raphson loops", (double) this->m_stat_newton_raphson / (double) this->m_stat_vsolver_calls);
this->netlist().log(" %10d invocations (%6d Hz) %10d gs fails (%6.2f%%) %6.3f average",
this->m_stat_calculations,
this->m_stat_calculations * 10 / (int) (this->netlist().time().as_double() * 10.0),
this->m_gs_fail,
100.0 * (double) this->m_gs_fail / (double) this->m_stat_calculations,
(double) this->m_gs_total / (double) this->m_stat_calculations);
}
}
template <unsigned m_N, unsigned _storage_N>
void matrix_solver_SOR_t<m_N, _storage_N>::vsetup(analog_net_t::list_t &nets)
{
matrix_solver_direct_t<m_N, _storage_N>::vsetup(nets);
this->save(NLNAME(m_lp_fact));
this->save(NLNAME(m_gs_fail));
this->save(NLNAME(m_gs_total));
}
template <unsigned m_N, unsigned _storage_N>
@ -194,13 +167,13 @@ ATTR_HOT inline int matrix_solver_SOR_t<m_N, _storage_N>::vsolve_non_dynamic(con
//} while (resched && (resched_cnt < this->m_params.m_gs_loops));
} while (resched && ((!interleaved_dynamic_updates && resched_cnt < this->m_params.m_gs_loops) || (interleaved_dynamic_updates && resched_cnt < 5 )));
this->m_gs_total += resched_cnt;
this->m_iterative_total += resched_cnt;
this->m_stat_calculations++;
if (resched && !interleaved_dynamic_updates)
{
// Fallback to direct solver ...
this->m_gs_fail++;
this->m_iterative_fail++;
return matrix_solver_direct_t<m_N, _storage_N>::vsolve_non_dynamic(newton_raphson);
}

View File

@ -87,6 +87,8 @@ ATTR_COLD matrix_solver_t::matrix_solver_t(const eSolverType type, const solver_
: m_stat_calculations(0),
m_stat_newton_raphson(0),
m_stat_vsolver_calls(0),
m_iterative_fail(0),
m_iterative_total(0),
m_params(*params),
m_cur_ts(0),
m_type(type)
@ -177,6 +179,7 @@ ATTR_COLD void matrix_solver_t::setup(analog_net_t::list_t &nets)
}
NL_VERBOSE_OUT(("added net with %" SIZETFMT " populated connections\n", net->m_core_terms.size()));
}
}
@ -203,6 +206,11 @@ ATTR_COLD void matrix_solver_t::start()
save(NLNAME(m_last_step));
save(NLNAME(m_cur_ts));
save(NLNAME(m_stat_calculations));
save(NLNAME(m_stat_newton_raphson));
save(NLNAME(m_stat_vsolver_calls));
save(NLNAME(m_iterative_fail));
save(NLNAME(m_iterative_total));
}
@ -319,6 +327,8 @@ NETLIB_START(solver)
register_param("FREQ", m_freq, 48000.0);
register_param("ITERATIVE", m_iterative_solver, "SOR");
register_param("ACCURACY", m_accuracy, 1e-7);
register_param("GS_LOOPS", m_gs_loops, 9); // Gauss-Seidel loops
register_param("GS_THRESHOLD", m_gs_threshold, 6); // below this value, gaussian elimination is used
@ -411,7 +421,7 @@ NETLIB_UPDATE(solver)
}
template <int m_N, int _storage_N>
matrix_solver_t * NETLIB_NAME(solver)::create_solver(int size, const int gs_threshold, const bool use_specific)
matrix_solver_t * NETLIB_NAME(solver)::create_solver(int size, const bool use_specific)
{
if (use_specific && m_N == 1)
return palloc(matrix_solver_direct1_t, &m_params);
@ -419,19 +429,25 @@ matrix_solver_t * NETLIB_NAME(solver)::create_solver(int size, const int gs_thre
return palloc(matrix_solver_direct2_t, &m_params);
else
{
if (size >= gs_threshold)
if (size >= m_gs_threshold)
{
if (USE_MATRIX_GS)
if (pstring("SOR_MAT").equals(m_iterative_solver))
{
typedef matrix_solver_SOR_mat_t<m_N,_storage_N> solver_mat;
return palloc(solver_mat, &m_params, size);
}
else
else if (pstring("SOR").equals(m_iterative_solver))
{
typedef matrix_solver_SOR_t<m_N,_storage_N> solver_GS;
//typedef matrix_solver_GMRES_t<m_N,_storage_N> solver_GS;
return palloc(solver_GS, &m_params, size);
}
else if (pstring("GMRES").equals(m_iterative_solver))
{
typedef matrix_solver_GMRES_t<m_N,_storage_N> solver_GMRES;
return palloc(solver_GMRES, &m_params, size);
}
else
netlist().error("Unknown solver type: %s\n", m_iterative_solver.Value().cstr());
}
else
{
@ -445,7 +461,6 @@ ATTR_COLD void NETLIB_NAME(solver)::post_start()
{
analog_net_t::list_t groups[256];
int cur_group = -1;
const int gs_threshold = m_gs_threshold.Value();
const bool use_specific = true;
m_params.m_accuracy = m_accuracy.Value();
@ -502,52 +517,52 @@ ATTR_COLD void NETLIB_NAME(solver)::post_start()
switch (net_count)
{
case 1:
ms = create_solver<1,1>(1, gs_threshold, use_specific);
ms = create_solver<1,1>(1, use_specific);
break;
case 2:
ms = create_solver<2,2>(2, gs_threshold, use_specific);
ms = create_solver<2,2>(2, use_specific);
break;
case 3:
ms = create_solver<3,3>(3, gs_threshold, use_specific);
ms = create_solver<3,3>(3, use_specific);
break;
case 4:
ms = create_solver<4,4>(4, gs_threshold, use_specific);
ms = create_solver<4,4>(4, use_specific);
break;
case 5:
ms = create_solver<5,5>(5, gs_threshold, use_specific);
ms = create_solver<5,5>(5, use_specific);
break;
case 6:
ms = create_solver<6,6>(6, gs_threshold, use_specific);
ms = create_solver<6,6>(6, use_specific);
break;
case 7:
ms = create_solver<7,7>(7, gs_threshold, use_specific);
ms = create_solver<7,7>(7, use_specific);
break;
case 8:
ms = create_solver<8,8>(8, gs_threshold, use_specific);
ms = create_solver<8,8>(8, use_specific);
break;
case 12:
ms = create_solver<12,12>(12, gs_threshold, use_specific);
ms = create_solver<12,12>(12, use_specific);
break;
case 87:
ms = create_solver<87,87>(87, gs_threshold, use_specific);
ms = create_solver<87,87>(87, use_specific);
break;
default:
if (net_count <= 16)
{
ms = create_solver<0,16>(net_count, gs_threshold, use_specific);
ms = create_solver<0,16>(net_count, use_specific);
}
else if (net_count <= 32)
{
ms = create_solver<0,32>(net_count, gs_threshold, use_specific);
ms = create_solver<0,32>(net_count, use_specific);
}
else if (net_count <= 64)
{
ms = create_solver<0,64>(net_count, gs_threshold, use_specific);
ms = create_solver<0,64>(net_count, use_specific);
}
else
if (net_count <= 128)
{
ms = create_solver<0,128>(net_count, gs_threshold, use_specific);
ms = create_solver<0,128>(net_count, use_specific);
}
else
{

View File

@ -131,10 +131,29 @@ public:
virtual void reset();
ATTR_COLD int get_net_idx(net_t *net);
virtual void log_stats() {};
inline eSolverType type() const { return m_type; }
virtual void log_stats()
{
if (this->m_stat_calculations != 0 && this->m_params.m_log_stats)
{
this->netlist().log("==============================================");
this->netlist().log("Solver %s", this->name().cstr());
this->netlist().log(" ==> %d nets", (unsigned) this->m_nets.size()); //, (*(*groups[i].first())->m_core_terms.first())->name().cstr());
this->netlist().log(" has %s elements", this->is_dynamic() ? "dynamic" : "no dynamic");
this->netlist().log(" has %s elements", this->is_timestep() ? "timestep" : "no timestep");
this->netlist().log(" %6.3f average newton raphson loops", (double) this->m_stat_newton_raphson / (double) this->m_stat_vsolver_calls);
this->netlist().log(" %10d invocations (%6d Hz) %10d gs fails (%6.2f%%) %6.3f average",
this->m_stat_calculations,
this->m_stat_calculations * 10 / (int) (this->netlist().time().as_double() * 10.0),
this->m_iterative_fail,
100.0 * (double) this->m_iterative_fail / (double) this->m_stat_calculations,
(double) this->m_iterative_total / (double) this->m_stat_calculations);
}
}
protected:
ATTR_COLD void setup(analog_net_t::list_t &nets);
@ -151,6 +170,8 @@ protected:
int m_stat_calculations;
int m_stat_newton_raphson;
int m_stat_vsolver_calls;
int m_iterative_fail;
int m_iterative_total;
const solver_parameters_t &m_params;
@ -205,6 +226,7 @@ protected:
param_logic_t m_dynamic;
param_double_t m_min_timestep;
param_str_t m_iterative_solver;
param_int_t m_nr_loops;
param_int_t m_gs_loops;
param_int_t m_gs_threshold;
@ -218,7 +240,7 @@ private:
solver_parameters_t m_params;
template <int m_N, int _storage_N>
matrix_solver_t *create_solver(int size, int gs_threshold, bool use_specific);
matrix_solver_t *create_solver(int size, bool use_specific);
};
NETLIB_NAMESPACE_DEVICES_END()

View File

@ -402,13 +402,14 @@ NETLIST_START(kidniki_interface)
PARAM(Solver.GS_LOOPS, 2)
#else
SOLVER(Solver, 12000)
PARAM(Solver.ACCURACY, 1e-8)
PARAM(Solver.ACCURACY, 1e-9)
PARAM(Solver.NR_LOOPS, 300)
PARAM(Solver.GS_LOOPS, 2)
PARAM(Solver.GS_LOOPS, 1)
PARAM(Solver.GS_THRESHOLD, 99)
PARAM(Solver.ITERATIVE, "SOR")
#endif
//PARAM(Solver.GS_THRESHOLD, 99) // Force Gaussian elimination here
PARAM(Solver.SOR_FACTOR, 1.05)
PARAM(Solver.SOR_FACTOR, 1.00)
//FIXME proper models!
NET_MODEL(".model 2SC945 NPN(Is=2.04f Xti=3 Eg=1.11 Vaf=6 Bf=400 Ikf=20m Xtb=1.5 Br=3.377 Rc=1 Cjc=1p Mjc=.3333 Vjc=.75 Fc=.5 Cje=25p Mje=.3333 Vje=.75 Tr=450n Tf=20n Itf=0 Vtf=0 Xtf=0 VCEO=45V ICrating=150M MFG=Toshiba)")
NET_MODEL(".model 1S1588 D(Is=2.52n Rs=.568 N=1.752 Cjo=4p M=.4 tt=20n Iave=200m Vpk=75 mfg=OnSemi type=silicon)")