diff options
Diffstat (limited to 'src/lib/netlist/solver/nld_ms_sor_mat.h')
-rw-r--r-- | src/lib/netlist/solver/nld_ms_sor_mat.h | 52 |
1 files changed, 29 insertions, 23 deletions
diff --git a/src/lib/netlist/solver/nld_ms_sor_mat.h b/src/lib/netlist/solver/nld_ms_sor_mat.h index 5eb4405600e..1892ed472f5 100644 --- a/src/lib/netlist/solver/nld_ms_sor_mat.h +++ b/src/lib/netlist/solver/nld_ms_sor_mat.h @@ -33,13 +33,10 @@ public: typedef FT float_type; matrix_solver_SOR_mat_t(netlist_base_t &anetlist, const pstring &name, const solver_parameters_t *params, std::size_t size) - : matrix_solver_direct_t<FT, SIZE>(anetlist, name, matrix_solver_t::DESCENDING, params, size) + : matrix_solver_direct_t<FT, SIZE>(anetlist, name, matrix_solver_t::ASCENDING, params, size) , m_Vdelta(*this, "m_Vdelta", std::vector<float_type>(size)) , m_omega(*this, "m_omega", params->m_gs_sor) , m_lp_fact(*this, "m_lp_fact", 0) - , m_gs_fail(*this, "m_gs_fail", 0) - , m_gs_total(*this, "m_gs_total", 0) - , new_V(size, 0.0) { } @@ -55,10 +52,7 @@ private: state_var<float_type> m_omega; state_var<float_type> m_lp_fact; - state_var<int> m_gs_fail; - state_var<int> m_gs_total; - std::vector<float_type> new_V; }; // ---------------------------------------------------------------------------------------- @@ -81,13 +75,13 @@ float_type matrix_solver_SOR_mat_t<m_N, storage_N>::vsolve() */ if (USE_LINEAR_PREDICTION) - for (unsigned k = 0; k < this->N(); k++) + for (unsigned k = 0; k < this->size(); k++) { this->m_last_V[k] = this->m_nets[k]->m_cur_Analog; this->m_nets[k]->m_cur_Analog = this->m_nets[k]->m_cur_Analog + this->m_Vdelta[k] * this->current_timestep() * m_lp_fact; } else - for (unsigned k = 0; k < this->N(); k++) + for (unsigned k = 0; k < this->size(); k++) { this->m_last_V[k] = this->m_nets[k]->m_cur_Analog; } @@ -99,7 +93,7 @@ float_type matrix_solver_SOR_mat_t<m_N, storage_N>::vsolve() float_type sq = 0; float_type sqo = 0; const float_type rez_cts = 1.0 / this->current_timestep(); - for (unsigned k = 0; k < this->N(); k++) + for (unsigned k = 0; k < this->size(); k++) { const analog_net_t *n = this->m_nets[k]; const float_type nv = (n->Q_Analog() - this->m_last_V[k]) * rez_cts ; @@ -129,7 +123,7 @@ unsigned matrix_solver_SOR_mat_t<FT, SIZE>::vsolve_non_dynamic(const bool newton */ - const std::size_t iN = this->N(); + const std::size_t iN = this->size(); this->build_LE_A(*this); this->build_LE_RHS(*this); @@ -176,7 +170,7 @@ unsigned matrix_solver_SOR_mat_t<FT, SIZE>::vsolve_non_dynamic(const bool newton #endif for (std::size_t k = 0; k < iN; k++) - new_V[k] = this->m_nets[k]->Q_Analog(); + this->m_new_V[k] = this->m_nets[k]->Q_Analog(); do { resched = false; @@ -190,11 +184,26 @@ unsigned matrix_solver_SOR_mat_t<FT, SIZE>::vsolve_non_dynamic(const bool newton const std::size_t e = this->m_terms[k]->m_nz.size(); for (std::size_t i = 0; i < e; i++) - Idrive = Idrive + this->A(k,p[i]) * new_V[p[i]]; - - const float_type delta = m_omega * (this->RHS(k) - Idrive) / this->A(k,k); + Idrive = Idrive + this->A(k,p[i]) * this->m_new_V[p[i]]; + + FT w = m_omega / this->A(k,k); + if (USE_GABS) + { + FT gabs_t = 0.0; + for (std::size_t i = 0; i < e; i++) + if (p[i] != k) + gabs_t = gabs_t + std::abs(this->A(k,p[i])); + + gabs_t *= NL_FCONST(1.0); // derived by try and error + if (gabs_t > this->A(k,k)) + { + w = NL_FCONST(1.0) / (this->A(k,k) + gabs_t); + } + } + + const float_type delta = w * (this->RHS(k) - Idrive) ; cerr = std::max(cerr, std::abs(delta)); - new_V[k] += delta; + this->m_new_V[k] += delta; } if (cerr > this->m_params.m_accuracy) @@ -206,20 +215,17 @@ unsigned matrix_solver_SOR_mat_t<FT, SIZE>::vsolve_non_dynamic(const bool newton this->m_stat_calculations++; this->m_iterative_total += resched_cnt; - this->m_gs_total += resched_cnt; if (resched) { this->m_iterative_fail++; //this->netlist().warning("Falling back to direct solver .. Consider increasing RESCHED_LOOPS"); - this->m_gs_fail++; - return matrix_solver_direct_t<FT, SIZE>::solve_non_dynamic(newton_raphson); } - else { - this->store(new_V.data()); - return resched_cnt; - } + + 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; } |