From 8334ef7c61f700cc9c0469386f1ea69cbe7a5ad4 Mon Sep 17 00:00:00 2001 From: DONNOT Benjamin Date: Fri, 8 Dec 2023 16:54:17 +0100 Subject: [PATCH] python calls updated, now need to solve the segfault [skip ci] --- lightsim2grid/lightSimBackend.py | 11 +++-- src/BaseSolver.cpp | 1 + src/ChooseSolver.cpp | 77 ++++++++++++++++++++++++++++++++ src/ChooseSolver.h | 16 ++++++- src/DataDCLine.cpp | 2 - src/GridModel.cpp | 69 +++++++++++++++++++++++----- src/GridModel.h | 4 +- 7 files changed, 159 insertions(+), 21 deletions(-) diff --git a/lightsim2grid/lightSimBackend.py b/lightsim2grid/lightSimBackend.py index cf81459..a8c66c5 100644 --- a/lightsim2grid/lightSimBackend.py +++ b/lightsim2grid/lightSimBackend.py @@ -735,7 +735,7 @@ def _aux_finish_setup_after_reading(self): self.storage_theta = np.full(cls.n_storage, dtype=dt_float, fill_value=np.NaN).reshape(-1) self._count_object_per_bus() - self._grid.tell_topo_changed() + self._grid.tell_solver_need_reset() self.__me_at_init = self._grid.copy() self.__init_topo_vect = np.ones(cls.dim_topo, dtype=dt_int) self.__init_topo_vect[:] = self.topo_vect @@ -894,7 +894,9 @@ def runpf(self, is_dc=False): self._grid.deactivate_result_computation() # if I init with dc values, it should depends on previous state self.V[:] = self._grid.get_init_vm_pu() # see issue 30 + print("before dc pf") Vdc = self._grid.dc_pf(copy.deepcopy(self.V), self.max_it, self.tol) + print("after dc pf") self._grid.reactivate_result_computation() if Vdc.shape[0] == 0: raise BackendError(f"Divergence of DC powerflow (non connected grid) at the initialization of AC powerflow. Detailed error: {self._grid.get_dc_solver().get_error()}") @@ -903,6 +905,7 @@ def runpf(self, is_dc=False): V_init = copy.deepcopy(self.V) tick = time.perf_counter() self._timer_preproc += tick - beg_preproc + print("before ac pf") V = self._grid.ac_pf(V_init, self.max_it, self.tol) self._timer_solver += time.perf_counter() - tick if V.shape[0] == 0: @@ -968,11 +971,11 @@ def runpf(self, is_dc=False): if (self.line_or_theta >= 1e6).any() or (self.line_ex_theta >= 1e6).any(): raise BackendError(f"Some theta are above 1e6 which should not be happening !") res = True - self._grid.unset_topo_changed() + self._grid.unset_changes() self._timer_postproc += time.perf_counter() - beg_postroc except Exception as exc_: # of the powerflow has not converged, results are Nan - self._grid.tell_topo_changed() + self._grid.unset_changes() self._fill_nans() res = False my_exc_ = exc_ @@ -1190,7 +1193,7 @@ def get_current_solver_type(self): def reset(self, grid_path, grid_filename=None): self._fill_nans() self._grid = self.__me_at_init.copy() - self._grid.tell_topo_changed() + self._grid.unset_changes() self._grid.change_solver(self.__current_solver_type) self._handle_turnedoff_pv() self.topo_vect[:] = self.__init_topo_vect diff --git a/src/BaseSolver.cpp b/src/BaseSolver.cpp index c2df65e..61b9d49 100644 --- a/src/BaseSolver.cpp +++ b/src/BaseSolver.cpp @@ -23,6 +23,7 @@ void BaseSolver::reset(){ err_ = ErrorType::NotInitError; //error message: _solver_control = SolverControl(); + _solver_control.tell_all_changed(); } diff --git a/src/ChooseSolver.cpp b/src/ChooseSolver.cpp index 165f787..bfd87dd 100644 --- a/src/ChooseSolver.cpp +++ b/src/ChooseSolver.cpp @@ -7,3 +7,80 @@ // This file is part of LightSim2grid, LightSim2grid implements a c++ backend targeting the Grid2Op platform. #include "ChooseSolver.h" + +std::ostream& operator<<(std::ostream& out, const SolverType& solver_type) +{ + switch (solver_type) + { + case SolverType::SparseLU: + out << "SparseLU"; + break; + case SolverType::KLU: + out << "KLU"; + break; + case SolverType::GaussSeidel: + out << "GaussSeidel"; + break; + case SolverType::DC: + out << "DC"; + break; + case SolverType::GaussSeidelSynch: + out << "GaussSeidelSynch"; + break; + case SolverType::NICSLU: + out << "NICSLU"; + break; + case SolverType::SparseLUSingleSlack: + out << "SparseLUSingleSlack"; + break; + case SolverType::KLUSingleSlack: + out << "KLUSingleSlack"; + break; + case SolverType::NICSLUSingleSlack: + out << "NICSLUSingleSlack"; + break; + case SolverType::KLUDC: + out << "KLUDC"; + break; + case SolverType::NICSLUDC: + out << "NICSLUDC"; + break; + case SolverType::CKTSO: + out << "CKTSO"; + break; + case SolverType::CKTSOSingleSlack: + out << "CKTSOSingleSlack"; + break; + case SolverType::CKTSODC: + out << "CKTSODC"; + break; + case SolverType::FDPF_XB_SparseLU: + out << "FDPF_XB_SparseLU"; + break; + case SolverType::FDPF_BX_SparseLU: + out << "FDPF_BX_SparseLU"; + break; + case SolverType::FDPF_XB_KLU: + out << "FDPF_XB_KLU"; + break; + case SolverType::FDPF_BX_KLU: + out << "FDPF_BX_KLU"; + break; + case SolverType::FDPF_XB_NICSLU: + out << "FDPF_XB_NICSLU"; + break; + case SolverType::FDPF_BX_NICSLU: + out << "FDPF_BX_NICSLU"; + break; + case SolverType::FDPF_XB_CKTSO: + out << "FDPF_XB_CKTSO"; + break; + case SolverType::FDPF_BX_CKTSO: + out << "FDPF_BX_CKTSO"; + break; + default: + out << "(unknown)"; + break; + } + return out; +} diff --git a/src/ChooseSolver.h b/src/ChooseSolver.h index 08fb0cd..2961d8a 100644 --- a/src/ChooseSolver.h +++ b/src/ChooseSolver.h @@ -27,6 +27,8 @@ enum class SolverType {SparseLU, KLU, GaussSeidel, DC, GaussSeidelSynch, NICSLU, FDPF_XB_CKTSO, FDPF_BX_CKTSO // from 0.7.5 }; + +std::ostream& operator<<(std::ostream& out, const SolverType& solver_type); // TODO define a template class instead of these weird stuff !!! // TODO export all methods from base class ! @@ -278,7 +280,7 @@ class ChooseSolver } void tell_solver_control(const SolverControl & solver_control){ - auto p_solver = get_prt_solver("get_ptdf", true); + auto p_solver = get_prt_solver("tell_solver_control", false); p_solver -> tell_solver_control(solver_control); } @@ -312,7 +314,17 @@ class ChooseSolver private: void check_right_solver(const std::string & error_msg) const { - if(_solver_type != _type_used_for_nr) throw std::runtime_error("ChooseSolver: Solver mismatch when calling '"+error_msg+"': current solver is not the last solver used to perform a powerflow"); + if(_solver_type != _type_used_for_nr){ + std::ostringstream exc_; + exc_ << "ChooseSolver: Solver mismatch when calling '"; + exc_ << error_msg; + exc_ << ": current solver ("; + exc_ << _solver_type; + exc_ << ") is not the one used to perform a powerflow ("; + exc_ << _type_used_for_nr; + exc_ << ")."; + throw std::runtime_error(exc_.str()); + } #ifndef KLU_SOLVER_AVAILABLE if(_solver_type == SolverType::KLU){ diff --git a/src/DataDCLine.cpp b/src/DataDCLine.cpp index e7c51c4..c74301b 100644 --- a/src/DataDCLine.cpp +++ b/src/DataDCLine.cpp @@ -87,11 +87,9 @@ void DataDCLine::disconnect_if_not_in_main_component(std::vector & busbar_ auto bus_or = bus_or_id(i); auto bus_ex = bus_ex_id(i); if(!busbar_in_main_component[bus_or]) { - bool tmp = false; from_gen_.deactivate(i, unused_solver_control); } if(!busbar_in_main_component[bus_ex]) { - bool tmp = false; to_gen_.deactivate(i, unused_solver_control); } // if(!busbar_in_main_component[bus_or] || !busbar_in_main_component[bus_ex]){ diff --git a/src/GridModel.cpp b/src/GridModel.cpp index 2aee652..f4768f2 100644 --- a/src/GridModel.cpp +++ b/src/GridModel.cpp @@ -414,25 +414,67 @@ CplxVect GridModel::pre_process_solver(const CplxVect & Vinit, const SolverControl & solver_control) { // TODO get rid of the "is_ac" argument: this info is available in the _solver already - if(is_ac) _solver.tell_solver_control(solver_control); - else _dc_solver.tell_solver_control(solver_control); - - if (solver_control.has_slack_participate_changed()) slack_bus_id_ = generators_.get_slack_bus_id(); - if (solver_control.ybus_change_sparsity_pattern() || solver_control.has_dimension_changed()) init_Ybus(Ybus, id_me_to_solver, id_solver_to_me); - if (solver_control.ybus_change_sparsity_pattern() || solver_control.has_dimension_changed() || solver_control.need_recompute_ybus()) fillYbus(Ybus, is_ac, id_me_to_solver); - if (solver_control.has_dimension_changed()) init_Sbus(Sbus_, id_me_to_solver, id_solver_to_me, slack_bus_id_solver); - if (solver_control.has_slack_participate_changed() || solver_control.has_pv_changed() || solver_control.has_pq_changed()) fillpv_pq(id_me_to_solver, id_solver_to_me, slack_bus_id_solver, solver_control); + if(is_ac){ + _solver.tell_solver_control(solver_control); + if(solver_control.need_reset_solver()) _solver.reset(); + } else { + _dc_solver.tell_solver_control(solver_control); + if(solver_control.need_reset_solver()){ + std::cout << "_dc_solver.reset();" << std::endl; + _dc_solver.reset(); + } + } + + if (solver_control.need_reset_solver() || + solver_control.has_slack_participate_changed()){ + std::cout << "get_slack_bus_id;" << std::endl; + slack_bus_id_ = generators_.get_slack_bus_id(); + } + if (solver_control.need_reset_solver() || + solver_control.ybus_change_sparsity_pattern() || + solver_control.has_dimension_changed()){ + init_Ybus(Ybus, id_me_to_solver, id_solver_to_me); + std::cout << "init_Ybus;" << std::endl; + } + if (solver_control.need_reset_solver() || + solver_control.ybus_change_sparsity_pattern() || + solver_control.has_dimension_changed() || + solver_control.need_recompute_ybus()){ + fillYbus(Ybus, is_ac, id_me_to_solver); + std::cout << "fillYbus;" << std::endl; + } + if (solver_control.need_reset_solver() || + solver_control.has_dimension_changed()) { + init_Sbus(Sbus_, id_me_to_solver, id_solver_to_me, slack_bus_id_solver); + std::cout << "init_Sbus;" << std::endl; + } + if (solver_control.need_reset_solver() || + solver_control.has_slack_participate_changed() || + solver_control.has_pv_changed() || + solver_control.has_pq_changed()) { + fillpv_pq(id_me_to_solver, id_solver_to_me, slack_bus_id_solver, solver_control); + std::cout << "fillpv_pq;" << std::endl; + } - if (solver_control.has_dimension_changed() || solver_control.need_recompute_sbus() && is_ac){ + if (is_ac && (solver_control.need_reset_solver() || + solver_control.has_dimension_changed() || + solver_control.need_recompute_sbus())){ int nb_bus_total = static_cast(bus_vn_kv_.size()); total_q_min_per_bus_ = RealVect::Constant(nb_bus_total, 0.); total_q_max_per_bus_ = RealVect::Constant(nb_bus_total, 0.); total_gen_per_bus_ = Eigen::VectorXi::Constant(nb_bus_total, 0); generators_.init_q_vector(nb_bus_total, total_gen_per_bus_, total_q_min_per_bus_, total_q_max_per_bus_); dc_lines_.init_q_vector(nb_bus_total, total_gen_per_bus_, total_q_min_per_bus_, total_q_max_per_bus_); - fillSbus_me(Sbus_, is_ac, id_me_to_solver); + std::cout << "total_gen_per_bus_;" << std::endl; } + if (solver_control.need_reset_solver() || + solver_control.has_slack_participate_changed() || + solver_control.has_pq_changed()) { + fillSbus_me(Sbus_, is_ac, id_me_to_solver); + std::cout << "fillSbus_me;" << std::endl; + } + const int nb_bus_solver = static_cast(id_solver_to_me.size()); CplxVect V = CplxVect::Constant(nb_bus_solver, init_vm_pu_); for(int bus_solver_id = 0; bus_solver_id < nb_bus_solver; ++bus_solver_id){ @@ -442,6 +484,7 @@ CplxVect GridModel::pre_process_solver(const CplxVect & Vinit, } generators_.set_vm(V, id_me_to_solver); dc_lines_.set_vm(V, id_me_to_solver); + std::cout << nb_bus_solver << std::endl; return V; } @@ -701,15 +744,17 @@ CplxVect GridModel::dc_pf(const CplxVect & Vinit, CplxVect V = pre_process_solver(Vinit, Ybus_dc_, id_me_to_dc_solver_, id_dc_solver_to_me_, slack_bus_id_dc_solver_, is_ac, solver_control_); - + std::cout << "after pre proces\n"; // start the solver if(solver_control_.has_slack_participate_changed() || solver_control_.has_pv_changed() || solver_control_.has_slack_weight_changed()) slack_weights_ = generators_.get_slack_weights(Ybus_dc_.rows(), id_me_to_dc_solver_); + std::cout << "slack_weights\n"; conv = _dc_solver.compute_pf(Ybus_dc_, V, dcSbus_, slack_bus_id_dc_solver_, slack_weights_, bus_pv_, bus_pq_, max_iter, tol); - + std::cout << "after compute_pf\n"; // store results (fase -> because I am in dc mode) process_results(conv, res, Vinit, false, id_me_to_dc_solver_); + std::cout << "after compute_pf\n"; return res; } diff --git a/src/GridModel.h b/src/GridModel.h index 8295534..0030c4d 100644 --- a/src/GridModel.h +++ b/src/GridModel.h @@ -76,8 +76,11 @@ class GridModel : public DataGeneric solver_control_(), init_vm_pu_(1.04), sn_mva_(1.0){ + _solver.change_solver(SolverType::SparseLU); _dc_solver.change_solver(SolverType::DC); _solver.set_gridmodel(this); + _dc_solver.set_gridmodel(this); + solver_control_.tell_all_changed(); } GridModel(const GridModel & other); GridModel copy() const{ @@ -166,7 +169,6 @@ class GridModel : public DataGeneric const RealVect & trafo_x, const CplxVect & trafo_b, const RealVect & trafo_tap_step_pct, - // const RealVect & trafo_tap_step_degree, const RealVect & trafo_tap_pos, const RealVect & trafo_shift_degree, const std::vector & trafo_tap_hv, // is tap on high voltage (true) or low voltate