diff --git a/src/gui/add_duct_dialog.cpp b/src/gui/add_duct_dialog.cpp index 2d4ce90..e481aee 100644 --- a/src/gui/add_duct_dialog.cpp +++ b/src/gui/add_duct_dialog.cpp @@ -141,8 +141,8 @@ void AddDuctDialog::accept(){ try { // If duct can be built without exceptions, everything is OK - TaSystem sys(TaSystem::testSystem()); - Duct duct (sys,0,_duct); + std::unique_ptr sys(TaSystem::testSystem()); + Duct duct (*sys,0,_duct); } catch(TaSMETError& e) { @@ -220,10 +220,10 @@ void AddDuctDialog::changed(){ PreviewShow pshow = (PreviewShow) _dialog->previewshow->currentIndex(); - TaSystem sys = TaSystem::testSystem(); + std::unique_ptr sys(TaSystem::testSystem()); std::unique_ptr duct; try { - duct = std::unique_ptr(new Duct(sys,0,_duct)); + duct = std::unique_ptr(new Duct(*sys,0,_duct)); } catch(TaSMETError& e) { return; diff --git a/src/gui/mainwindow.cpp b/src/gui/mainwindow.cpp index 86bfdb2..77d49d1 100644 --- a/src/gui/mainwindow.cpp +++ b/src/gui/mainwindow.cpp @@ -417,14 +417,43 @@ void TaSMETMainWindow::on_actionSolve_triggered() { TRACE(15,"actionSolve()"); SolverDialog *d; + + std::unique_ptr sys; try { - d = new SolverDialog(this,_system,*_model.mutable_sparams()); + sys = std::unique_ptr(new TaSystem(_model.system())); + + if(_model.solution_size()>0) { + vd solution(_model.solution_size()); + for(us i=0;i<_model.solution_size();i++) { + solution(i) = _model.solution(i); + } + + // Throws in case the solution size does not match. + try { + sys->updateSolution(solution); + } + catch(...) {} + } + d = new SolverDialog(this,*sys.get(),*_model.mutable_sparams()); } + catch(TaSMETError &e) { e.show_user("Solver failed to initialize"); return; } - d->exec(); + + if(d->exec()) { + + // On succes, we copy the solution to the model + const vd& sol = d->getSolution(); + + // Clear old solution + _model.clear_solution(); + for(us i=0;i - &System type: + S&ystem type: systemtype @@ -151,7 +151,7 @@ - &Fundamental frequency: + F&undamental frequency: freq @@ -348,7 +348,7 @@ 0 0 683 - 21 + 18 @@ -372,11 +372,11 @@ - Solver + So&lver - - + + @@ -446,17 +446,17 @@ - Solve... + &Solve... - + - Reinitialize solver + &Reinitialize solution - + - Posprocess model... + Create postprocessing file diff --git a/src/gui/solver_dialog.cpp b/src/gui/solver_dialog.cpp index 1e15de6..8f7849b 100644 --- a/src/gui/solver_dialog.cpp +++ b/src/gui/solver_dialog.cpp @@ -22,7 +22,7 @@ #include SolverDialog::SolverDialog(QWidget* parent, - pb::System& sys, + const GradientNonlinearSystem& sys, pb::SolverParams& sparams): QDialog(parent), _sys(sys), @@ -125,8 +125,13 @@ void SolverDialog::solver_progress(const SolverProgress& progress){ TRACE(15,"SolverDialog::solver_progress()"); - // VARTRACE(15,progress.fun_err); - // VARTRACE(15,progress.iteration); + if(progress.error) { + QMessageBox::warning(this, + "Solver error", + QString::fromStdString(progress.err_msg)); + + } + d funtol = _sparams.funtol(); d reltol = _sparams.reltol(); @@ -159,17 +164,18 @@ void SolverDialog::on_solve_clicked() { _funtol->setData(empty,empty); _reltol->setData(empty,empty); - assert(!_solver_worker); qRegisterMetaType(); _solver_worker = new SolverWorker(_sys,_sparams); + QThread* thread = new QThread; + if(!_solver_worker || !thread) throw TaSMETBadAlloc(); + // Move the _solver_worker to its own thread _solver_worker->moveToThread(thread); - connect(thread, &QThread::started, _solver_worker, &SolverWorker::solver_start); @@ -192,12 +198,9 @@ void SolverDialog::on_solve_clicked() { } -void SolverDialog::on_singleiteration_clicked() { - -} + void SolverDialog::setEnabled(bool enabled){ _dialog->solve->setEnabled(enabled); - _dialog->singleiteration->setEnabled(enabled); _dialog->stop->setEnabled(!enabled); _dialog->solvertype->setEnabled(enabled); @@ -210,6 +213,9 @@ void SolverDialog::solver_stopped(bool converged) { // stop the solver and delete it if(_solver_worker!=nullptr) { _solver_worker->solver_stop(); + if(converged) { + _solution = _solver_worker->getSolution(); + } _solver_worker = nullptr; } diff --git a/src/gui/solver_dialog.h b/src/gui/solver_dialog.h index 87b57f1..b8b075a 100644 --- a/src/gui/solver_dialog.h +++ b/src/gui/solver_dialog.h @@ -10,7 +10,7 @@ #define SOLVER_DIALOG_H #include #include - +#include "tasmet_types.h" #include "solver.pb.h" namespace pb { @@ -20,7 +20,7 @@ namespace pb { namespace Ui { class solver_dialog; } - +class GradientNonlinearSystem; class QCustomPlot; class QCPGraph; class SolverWorker; @@ -29,7 +29,6 @@ class SolverProgress; class SolverDialog: public QDialog { Q_OBJECT - pb::System& _sys; // Reference to system pb::SolverParams& _sparams; Ui::solver_dialog* _dialog; @@ -41,14 +40,21 @@ class SolverDialog: public QDialog { SolverWorker* _solver_worker = nullptr; + const GradientNonlinearSystem& _sys; + + /// Place where the final solution will be stored + vd _solution; + public: SolverDialog(QWidget* parent, - pb::System& sys, + const GradientNonlinearSystem& sys, pb::SolverParams& sparams); ~SolverDialog(); + const vd& getSolution() const { return _solution;}; + public slots: void solver_progress(const SolverProgress&); private slots: @@ -57,17 +63,23 @@ private slots: void solver_stopped(bool converged); void on_solve_clicked(); - void on_singleiteration_clicked(); void on_stop_clicked(); void on_funtol_textChanged() { changed();} void on_reltol_textChanged() { changed();} void on_solvertype_currentIndexChanged(int) { changed();} + + void on_buttons_accepted() {accept();} + void on_buttons_rejected() {reject();} + private: + // Called whenever the user changes input values void changed(); void setEnabled(bool); void set(const pb::SolverParams&); + + }; diff --git a/src/gui/solver_dialog.ui b/src/gui/solver_dialog.ui index a452900..81dd7f3 100644 --- a/src/gui/solver_dialog.ui +++ b/src/gui/solver_dialog.ui @@ -13,143 +13,149 @@ TaSMET Solver - + - - - - - - 0 - 0 - - - - Solver Settings - - - - - - Solver type - - - - - - - - - - - 0 - 0 - - - - Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter - - - - - - - Solution tolerance: - - - - - - - Residual tolerance: - - - - - - - - 0 - 0 - - - - Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter - - - - - - - - - - - 0 - 0 - - - - Solve! - - - - - - Solve - - - - - - - Stop - - - - - - - Single iteration - - - - - - - - - - Qt::Vertical - - - - 20 - 40 - - - - - - - - - - - 0 - 0 - - - - Progress - - + + - + + + + + + 0 + 0 + + + + Solver Settings + + + + + + Solver type + + + + + + + + + + + 0 + 0 + + + + Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter + + + + + + + Solution tolerance: + + + + + + + Residual tolerance: + + + + + + + + 0 + 0 + + + + Qt::AlignRight|Qt::AlignTrailing|Qt::AlignVCenter + + + + + + + + + + + 0 + 0 + + + + Solve! + + + + + + Solve + + + + + + + Stop + + + + + + + + + + Qt::Vertical + + + + 20 + 40 + + + + + + + + + + + 0 + 0 + + + + Progress + + + + + + + + + + + QDialogButtonBox::Cancel|QDialogButtonBox::Save + + + diff --git a/src/gui/solver_worker.cpp b/src/gui/solver_worker.cpp index ba5f3f7..3c960bd 100644 --- a/src/gui/solver_worker.cpp +++ b/src/gui/solver_worker.cpp @@ -9,22 +9,38 @@ #include "solver_worker.h" #include #include "tasmet_tracer.h" - +#include "solver.h" #include "system.pb.h" #include "solver.pb.h" #include +// Solvers available +#include "newton_raphson.h" +#include "tasystem.h" -SolverWorker::SolverWorker(const pb::System& sys,const pb::SolverParams& sparams): +SolverWorker::SolverWorker(const GradientNonlinearSystem& sys, + const pb::SolverParams& sparams): _run(false), - _reltol(sparams.reltol()), - _funtol(sparams.funtol()) + _funtol(sparams.funtol()), + _reltol(sparams.reltol()) { + TRACE(15,"SolverWorker"); + switch (sparams.solvertype()) { + case pb::NewtonRaphson: { + _solver = new NewtonRaphson(sys); + break; + } + default: + tasmet_assert(false,"Not implemented solver type"); + break; + } + } SolverWorker::~SolverWorker(){ TRACE(15,"~SolverWorker"); + if(_solver!=nullptr) delete _solver; } void SolverWorker::solver_stop() { _run = false; @@ -39,23 +55,10 @@ void SolverWorker::solver_start() { progress_callback callback = std::bind(&SolverWorker::pg_callback, this,_1); - SolverProgress p; - // For testing purposes + + tasmet_assert(_solver!=nullptr,"Solver not initialized"); - SolverAction action; - while(true) { - TRACE(15,"Solver start virtual iteration"); - - - SolverAction action = callback(p); - if(action != Continue) break; - sleep(1); - - p.fun_err/=10; - p.rel_err/=10; - p.iteration++; - - } + _solver->start(&callback); emit solver_stopped(_converged); } @@ -66,11 +69,15 @@ SolverAction SolverWorker::pg_callback(SolverProgress pg) { emit progress(pg); + if(pg.error) { + _converged = false; + return Stop; + } if(pg.fun_err <= _funtol && pg.rel_err <= _reltol) { _converged = true; return Stop; } - + return Continue; diff --git a/src/gui/solver_worker.h b/src/gui/solver_worker.h index 0a9a5fb..94ce5e5 100644 --- a/src/gui/solver_worker.h +++ b/src/gui/solver_worker.h @@ -19,18 +19,28 @@ namespace pb{ class SolverParams; } +class TaSystem; + Q_DECLARE_METATYPE(SolverProgress); class SolverWorker: public QObject { Q_OBJECT std::atomic _run; - Solver* _solver; + Solver* _solver = nullptr; bool _converged = false; d _funtol,_reltol; public: - SolverWorker(const pb::System& sys,const pb::SolverParams& sparams); + SolverWorker(const GradientNonlinearSystem& sys, + const pb::SolverParams& sparams); + ~SolverWorker(); + + vd getSolution() const { return _solver ? _solver->getSolution() : zeros(0);} + + /** + * Stop the solver. Called when the user presses the `stop' button + */ void solver_stop(); public slots: diff --git a/src/main.cpp b/src/main.cpp index 9ab1d95..080f1dc 100644 --- a/src/main.cpp +++ b/src/main.cpp @@ -40,7 +40,7 @@ int main(int argc, char *argv[]) { // This can be used to store files in the application. Very useful // when our application grows bigger // Q_INIT_RESOURCE(application); - INITTRACE(15); + INITTRACE(16); QApplication app(argc, argv); diff --git a/src/solver/newton_raphson.cpp b/src/solver/newton_raphson.cpp index b9d5863..82aba6d 100644 --- a/src/solver/newton_raphson.cpp +++ b/src/solver/newton_raphson.cpp @@ -9,7 +9,7 @@ #include "newton_raphson.h" #include "tasmet_tracer.h" -#define DEBUG_TASMET_SYSTEM +// #define DEBUG_TASMET_SYSTEM void NewtonRaphson::start_implementation(GradientNonlinearSystem& system, progress_callback* callback) { @@ -48,7 +48,15 @@ void NewtonRaphson::start_implementation(GradientNonlinearSystem& system, #endif // DEBUG_TASMET_SYSTEM TRACE(15,"Solving system of eqs"); - dx = -1*_dampfac*arma::spsolve(resjac.jacobian,resjac.residual,"superlu"); + try { + dx = -1*_dampfac*arma::spsolve(resjac.jacobian,resjac.residual,"superlu"); + } + catch(...) { + progress.error = true; + progress.err_msg = "Failed to solve linear system of equations"; + (*callback)(progress); + return; + } TRACE(15,"Solving system of eqs done"); progress.rel_err = norm(dx); diff --git a/src/solver/solver.h b/src/solver/solver.h index f459ce9..6114629 100644 --- a/src/solver/solver.h +++ b/src/solver/solver.h @@ -32,6 +32,8 @@ struct SolverProgress Solver stops as an action of itself. Probably due to an internal error. */ + string err_msg; + }; diff --git a/src/sys/tasystem.cpp b/src/sys/tasystem.cpp index 0a5cc10..e7b8183 100644 --- a/src/sys/tasystem.cpp +++ b/src/sys/tasystem.cpp @@ -71,8 +71,7 @@ TaSystem::TaSystem(const pb::System& sys): } TaSystem::TaSystem(const TaSystem& o): GlobalConf(o), // Share a ptr to the Global conf - _gas(o._gas->copy()), - _solution(o._solution) + _gas(o._gas->copy()) { TRACE(25,"TaSystem::TaSystem(TaSystem&) copy"); @@ -84,7 +83,11 @@ TaSystem::TaSystem(const TaSystem& o): if(!seg.second) throw TaSMETBadAlloc(); } - initSolRes(); + initSolRes(); // This empties the + // solution. Therefore we copy it + // after this phase. + + _solution = o._solution; } void TaSystem::initSolRes() { @@ -272,7 +275,6 @@ void TaSystem::residualJac(ResidualJac& resjac) const { } - } @@ -295,6 +297,7 @@ void TaSystem::residualJac(ResidualJac& resjac) const { _residual(arbitrateMassEq)=mass - _mass; } + // Store the Jacobian and the residual resjac.jacobian = sdmat(triplets); resjac.residual = _residual; @@ -302,9 +305,20 @@ void TaSystem::residualJac(ResidualJac& resjac) const { vd TaSystem::getSolution() const { return _solution; } +void TaSystem::updateSolution(const vd& solution) { + + if(_solution.size() == solution.size()) { + _solution = solution; + } + else { + throw TaSMETError("Solution vector size does not match."); + } + +} const SegPositionMapper& TaSystem::getSolution(const us seg_id) const { return _solution_dofs.at(seg_id); } + vus TaSystem::getNDofs() const { TRACE(0,"TaSystem::getNDofs()"); vus Ndofs(_segs.size()); @@ -481,7 +495,7 @@ void TaSystem::exportHDF5(const string& filename) const { seg_.second->exportHDF5(grp_id); - + // Close the group H5Gclose(grp_id); } diff --git a/src/sys/tasystem.h b/src/sys/tasystem.h index 2ba4474..d95387f 100644 --- a/src/sys/tasystem.h +++ b/src/sys/tasystem.h @@ -64,14 +64,24 @@ private: * */ void initSolRes(); -public: + TaSystem(const TaSystem& o); +public: + TaSystem(const pb::System&); - static TaSystem testSystem() { - return TaSystem(pb::System::default_instance()); + /** + * Returns an empty TaSystem with default GlobalConf + * parameters. The caller of this function owns the object and is + * responsible for its destruction. + * + * + * @return the TaSystem instance pointer. + */ + static TaSystem* testSystem() { + return new TaSystem(pb::System::default_instance()); } const Gas& gas() const {return *_gas;} @@ -96,7 +106,7 @@ public: // Obtain the solution vector for the Segment with given id const SegPositionMapper& getSolution(const us seg_id) const; - virtual void updateSolution(const vd& sol) {_solution = sol; } // Update the solution + virtual void updateSolution(const vd& sol); // Change Nf in the system, while keeping the results. void updateNf(us);