diff options
Diffstat (limited to 'src/lib/netlist/solver/nld_ms_sor.h')
-rw-r--r-- | src/lib/netlist/solver/nld_ms_sor.h | 112 |
1 files changed, 55 insertions, 57 deletions
diff --git a/src/lib/netlist/solver/nld_ms_sor.h b/src/lib/netlist/solver/nld_ms_sor.h index c31aaa6d46a..eea692d6c47 100644 --- a/src/lib/netlist/solver/nld_ms_sor.h +++ b/src/lib/netlist/solver/nld_ms_sor.h @@ -12,42 +12,33 @@ #ifndef NLD_MS_SOR_H_ #define NLD_MS_SOR_H_ +#include <algorithm> + #include "nld_ms_direct.h" #include "nld_solver.h" -#include <algorithm> - namespace netlist { namespace devices -{ - -template <typename FT, int SIZE> -class matrix_solver_SOR_t: public matrix_solver_direct_t<FT, SIZE> + { +template <std::size_t m_N, std::size_t storage_N> +class matrix_solver_SOR_t: public matrix_solver_direct_t<m_N, storage_N> { public: - using float_type = FT; - - matrix_solver_SOR_t(netlist_state_t &anetlist, const pstring &name, const solver_parameters_t *params, const std::size_t size) - : matrix_solver_direct_t<FT, SIZE>(anetlist, name, matrix_solver_t::ASCENDING, params, size) + matrix_solver_SOR_t(netlist_t &anetlist, const pstring &name, const solver_parameters_t *params, const std::size_t size) + : matrix_solver_direct_t<m_N, storage_N>(anetlist, name, matrix_solver_t::ASCENDING, params, size) , m_lp_fact(*this, "m_lp_fact", 0) - , w(size, 0.0) - , one_m_w(size, 0.0) - , RHS(size, 0.0) - //, new_V(size, 0.0) { } - void vsetup(analog_net_t::list_t &nets) override; - unsigned vsolve_non_dynamic(const bool newton_raphson) override; + virtual ~matrix_solver_SOR_t() override {} + + virtual void vsetup(analog_net_t::list_t &nets) override; + virtual unsigned vsolve_non_dynamic(const bool newton_raphson) override; private: - state_var<float_type> m_lp_fact; - std::vector<float_type> w; - std::vector<float_type> one_m_w; - std::vector<float_type> RHS; - //std::vector<float_type> new_V; + state_var<nl_double> m_lp_fact; }; // ---------------------------------------------------------------------------------------- @@ -55,16 +46,16 @@ private: // ---------------------------------------------------------------------------------------- -template <typename FT, int SIZE> -void matrix_solver_SOR_t<FT, SIZE>::vsetup(analog_net_t::list_t &nets) +template <std::size_t m_N, std::size_t storage_N> +void matrix_solver_SOR_t<m_N, storage_N>::vsetup(analog_net_t::list_t &nets) { - matrix_solver_direct_t<FT, SIZE>::vsetup(nets); + matrix_solver_direct_t<m_N, storage_N>::vsetup(nets); } -template <typename FT, int SIZE> -unsigned matrix_solver_SOR_t<FT, SIZE>::vsolve_non_dynamic(const bool newton_raphson) +template <std::size_t m_N, std::size_t storage_N> +unsigned matrix_solver_SOR_t<m_N, storage_N>::vsolve_non_dynamic(const bool newton_raphson) { - const std::size_t iN = this->size(); + const std::size_t iN = this->N(); bool resched = false; unsigned resched_cnt = 0; @@ -76,21 +67,26 @@ unsigned matrix_solver_SOR_t<FT, SIZE>::vsolve_non_dynamic(const bool newton_rap * omega = 2.0 / (1.0 + std::sqrt(1-rho)) */ - const float_type ws = this->m_params.m_gs_sor; + const nl_double ws = this->m_params.m_gs_sor; + + nl_double w[storage_N]; + nl_double one_m_w[storage_N]; + nl_double RHS[storage_N]; + nl_double new_V[storage_N]; 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_double gtot_t = 0.0; + nl_double gabs_t = 0.0; + nl_double RHS_t = 0.0; const std::size_t term_count = this->m_terms[k]->count(); - const float_type * const gt = this->m_gtn[k]; - const float_type * const go = this->m_gonn[k]; - const float_type * const Idr = this->m_Idrn[k]; - auto other_cur_analog = this->m_connected_net_Vn[k]; + const nl_double * const RESTRICT gt = this->m_terms[k]->gt(); + const nl_double * const RESTRICT go = this->m_terms[k]->go(); + const nl_double * const RESTRICT Idr = this->m_terms[k]->Idr(); + const nl_double * const *other_cur_analog = this->m_terms[k]->connected_net_V(); - this->m_new_V[k] = this->m_nets[k]->Q_Analog(); + new_V[k] = this->m_nets[k]->Q_Analog(); for (std::size_t i = 0; i < term_count; i++) { @@ -99,60 +95,61 @@ unsigned matrix_solver_SOR_t<FT, SIZE>::vsolve_non_dynamic(const bool newton_rap } for (std::size_t i = this->m_terms[k]->m_railstart; i < term_count; i++) - RHS_t = RHS_t - go[i] * *other_cur_analog[i]; + RHS_t = RHS_t + go[i] * *other_cur_analog[i]; RHS[k] = RHS_t; - if (this->m_params.m_use_gabs) + if (USE_GABS) { for (std::size_t i = 0; i < term_count; i++) gabs_t = gabs_t + std::abs(go[i]); - gabs_t *= plib::constants<nl_double>::cast(0.5); // derived by try and error + gabs_t *= NL_FCONST(0.5); // derived by try and error if (gabs_t <= gtot_t) { w[k] = ws / gtot_t; - one_m_w[k] = plib::constants<FT>::one() - ws; + one_m_w[k] = NL_FCONST(1.0) - 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] = NL_FCONST(1.0) / (gtot_t + gabs_t); + one_m_w[k] = NL_FCONST(1.0) - NL_FCONST(1.0) * gtot_t / (gtot_t + gabs_t); } } else { w[k] = ws / gtot_t; - one_m_w[k] = plib::constants<FT>::one() - ws; + one_m_w[k] = NL_FCONST(1.0) - ws; } } - const float_type accuracy = this->m_params.m_accuracy; + const nl_double accuracy = this->m_params.m_accuracy; do { resched = false; - float_type err = 0; + nl_double err = 0; for (std::size_t k = 0; k < iN; k++) { - const int * net_other = this->m_terms[k]->m_connected_net_idx.data(); + const int * RESTRICT net_other = this->m_terms[k]->connected_net_idx(); const std::size_t railstart = this->m_terms[k]->m_railstart; - const float_type * go = this->m_gonn[k]; + const nl_double * RESTRICT go = this->m_terms[k]->go(); - float_type Idrive = 0.0; + nl_double Idrive = 0.0; 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 + go[i] * new_V[net_other[i]]; - const float_type new_val = this->m_new_V[k] * one_m_w[k] + (Idrive + RHS[k]) * w[k]; + const nl_double new_val = new_V[k] * one_m_w[k] + (Idrive + RHS[k]) * w[k]; - err = std::max(std::abs(new_val - this->m_new_V[k]), err); - this->m_new_V[k] = new_val; + err = std::max(std::abs(new_val - new_V[k]), err); + new_V[k] = new_val; } if (err > accuracy) resched = true; resched_cnt++; - } while (resched && (resched_cnt < this->m_params.m_gs_loops)); + //} while (resched && (resched_cnt < this->m_params.m_gs_loops)); + } while (resched && ((resched_cnt < this->m_params.m_gs_loops))); this->m_iterative_total += resched_cnt; this->m_stat_calculations++; @@ -161,12 +158,13 @@ unsigned matrix_solver_SOR_t<FT, SIZE>::vsolve_non_dynamic(const bool newton_rap { // Fallback to direct solver ... this->m_iterative_fail++; - return matrix_solver_direct_t<FT, SIZE>::vsolve_non_dynamic(newton_raphson); + return matrix_solver_direct_t<m_N, storage_N>::vsolve_non_dynamic(newton_raphson); } - const float_type err = (newton_raphson ? this->delta(this->m_new_V) : 0.0); - this->store(this->m_new_V); - return (err > this->m_params.m_accuracy) ? 2 : 1; + for (std::size_t k = 0; k < iN; k++) + this->m_nets[k]->set_Q_Analog(new_V[k]); + + return resched_cnt; } } //namespace devices |