diff --git a/python/wrap/solvers.cpp b/python/wrap/solvers.cpp index cf843fe..769718d 100644 --- a/python/wrap/solvers.cpp +++ b/python/wrap/solvers.cpp @@ -1,223 +1,238 @@ /* * SPDX-License-Indentifier: AGPL-3.0-or-later * * Copyright (©) 2016-2022 EPFL (École Polytechnique Fédérale de Lausanne), * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * This program is free software: you can redistribute it and/or modify * it under the terms of the GNU Affero General Public License as published * by the Free Software Foundation, either version 3 of the License, or * (at your option) any later version. * * This program is distributed in the hope that it will be useful, * but WITHOUT ANY WARRANTY; without even the implied warranty of * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the * GNU Affero General Public License for more details. * * You should have received a copy of the GNU Affero General Public License * along with this program. If not, see . * */ /* -------------------------------------------------------------------------- */ #include "beck_teboulle.hh" #include "condat.hh" #include "contact_solver.hh" #include "dfsane_solver.hh" #include "ep_solver.hh" #include "epic.hh" #include "kato.hh" #include "kato_saturated.hh" #include "polonsky_keer_rey.hh" #include "polonsky_keer_tan.hh" #include "wrap.hh" #include /* -------------------------------------------------------------------------- */ namespace tamaas { namespace wrap { /* -------------------------------------------------------------------------- */ using namespace py::literals; class PyEPSolver : public EPSolver { public: using EPSolver::EPSolver; // NOLINTNEXTLINE(readability-else-after-return) void solve() override { PYBIND11_OVERLOAD_PURE(void, EPSolver, solve); } void updateState() override { // NOLINTNEXTLINE(readability-else-after-return) PYBIND11_OVERLOAD(void, EPSolver, updateState); } }; +class PyContactSolver : public ContactSolver { +public: + using ContactSolver::ContactSolver; + + Real solve(std::vector load) override { + PYBIND11_OVERLOAD(Real, ContactSolver, solve, load); + } + + Real solve(Real load) override { + PYBIND11_OVERLOAD(Real, ContactSolver, solve, load); + } +}; + /* -------------------------------------------------------------------------- */ void wrapSolvers(py::module& mod) { - py::class_(mod, "ContactSolver") + py::class_(mod, "ContactSolver") + .def(py::init&, Real>(), + py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def( "setMaxIterations", [](ContactSolver& m, UInt iter) { TAMAAS_DEPRECATE("setMaxIterations()", "the max_iter property"); m.setMaxIterations(iter); }, "max_iter"_a) .def( "setDumpFrequency", [](ContactSolver& m, UInt freq) { TAMAAS_DEPRECATE("setDumpFrequency()", "the dump_freq property"); m.setDumpFrequency(freq); }, "dump_freq"_a) .def_property("tolerance", &ContactSolver::getTolerance, &ContactSolver::setTolerance, "Solver tolerance") .def_property("max_iter", &ContactSolver::getMaxIterations, &ContactSolver::setMaxIterations, "Maximum number of iterations") .def_property("dump_freq", &ContactSolver::getDumpFrequency, &ContactSolver::setDumpFrequency, "Frequency of displaying solver info") .def_property_readonly("model", &ContactSolver::getModel) .def_property_readonly("surface", &ContactSolver::getSurface) .def_property_readonly("functional", &ContactSolver::getFunctional) .def("addFunctionalTerm", &ContactSolver::addFunctionalTerm, "Add a term to the contact functional to minimize") .def("solve", py::overload_cast>(&ContactSolver::solve), "target_force"_a, py::call_guard(), "Solve the contact for a mean traction/gap vector") .def("solve", py::overload_cast(&ContactSolver::solve), "target_normal_pressure"_a, py::call_guard(), "Solve the contact for a mean normal pressure/gap"); py::class_ pkr( mod, "PolonskyKeerRey", "Main solver class for normal elastic contact problems. Its functional " "can be customized to add an adhesion term, and its primal variable can " "be set to either the gap or the pressure."); // Need to export enum values before defining PKR constructor py::enum_(pkr, "type") .value("gap", PolonskyKeerRey::gap) .value("pressure", PolonskyKeerRey::pressure) .export_values(); pkr.def(py::init&, Real, PolonskyKeerRey::type, PolonskyKeerRey::type>(), "model"_a, "surface"_a, "tolerance"_a, "primal_type"_a = PolonskyKeerRey::type::pressure, "constraint_type"_a = PolonskyKeerRey::type::pressure, py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def("computeError", &PolonskyKeerRey::computeError); py::class_( mod, "KatoSaturated", "Solver for pseudo-plasticity problems where the normal pressure is " "constrained above by a saturation pressure \"pmax\"") .def(py::init&, Real, Real>(), "model"_a, "surface"_a, "tolerance"_a, "pmax"_a, py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def_property("pmax", &KatoSaturated::getPMax, &KatoSaturated::setPMax, "Saturation normal pressure"); py::class_ kato(mod, "Kato"); kato.def(py::init&, Real, Real>(), "model"_a, "surface"_a, "tolerance"_a, "mu"_a, py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def("solve", &Kato::solve, "p0"_a, "proj_iter"_a = 50) .def("solveRelaxed", &Kato::solveRelaxed, "g0"_a) .def("solveRegularized", &Kato::solveRegularized, "p0"_a, "r"_a = 0.01) .def("computeCost", &Kato::computeCost, "use_tresca"_a = false); py::class_ bt(mod, "BeckTeboulle"); bt.def(py::init&, Real, Real>(), "model"_a, "surface"_a, "tolerance"_a, "mu"_a, py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def("solve", &BeckTeboulle::solve, "p0"_a) .def("computeCost", &BeckTeboulle::computeCost); py::class_ cd( mod, "Condat", "Main solver for frictional contact problems. It has no restraint on the " "material properties or friction coefficient values, but solves an " "associated version of the Coulomb friction law, which differs from the " "traditional Coulomb friction in that the normal and tangential slip " "components are coupled."); cd.def(py::init&, Real, Real>(), "model"_a, "surface"_a, "tolerance"_a, "mu"_a, py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def("solve", &Condat::solve, "p0"_a, "grad_step"_a = 0.9) .def("computeCost", &Condat::computeCost); py::class_ pkt(mod, "PolonskyKeerTan"); pkt.def(py::init&, Real, Real>(), "model"_a, "surface"_a, "tolerance"_a, "mu"_a, py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def("solve", &PolonskyKeerTan::solve, "p0"_a) .def("solveTresca", &PolonskyKeerTan::solveTresca, "p0"_a) .def("computeCost", &PolonskyKeerTan::computeCost, "use_tresca"_a = false); py::class_( mod, "_tolerance_manager", "Manager object for the tolereance of nonlinear plasticity solvers. " "Decreases the solver tolerance by geometric progression.") .def(py::init(), "start_tol"_a, "end_tol"_a, "rate"_a); py::class_( mod, "EPSolver", "Mother class for nonlinear plasticity solvers") .def(py::init(), "residual"_a, py::keep_alive<1, 2>()) .def("solve", &EPSolver::solve) .def("getStrainIncrement", &EPSolver::getStrainIncrement, py::return_value_policy::reference_internal) .def("getResidual", &EPSolver::getResidual, py::return_value_policy::reference_internal) .def("updateState", &EPSolver::updateState) .def_property("tolerance", &EPSolver::getTolerance, &EPSolver::setTolerance) .def("setToleranceManager", &EPSolver::setToleranceManager) .def("beforeSolve", &EPSolver::beforeSolve); py::class_(mod, "_DFSANESolver") .def(py::init(), "residual"_a, py::keep_alive<1, 2>()) .def(py::init([](Residual& res, Model&) { return std::make_unique(res); }), "residual"_a, "model"_a, py::keep_alive<1, 2>(), "Deprecated constructor: only here for compatibility reasons"); py::class_( mod, "EPICSolver", "Main solver class for elastic-plastic contact problems") .def(py::init(), "contact_solver"_a, "elasto_plastic_solver"_a, "tolerance"_a = 1e-10, "relaxation"_a = 0.3, py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def( "solve", [](EPICSolver& solver, Real pressure) { return solver.solve(std::vector{pressure}); }, "normal_pressure"_a, py::call_guard(), "Solves the EP contact with a relaxed fixed-point scheme. Adjust " "the relaxation parameter to help convergence.") .def( "acceleratedSolve", [](EPICSolver& solver, Real pressure) { return solver.acceleratedSolve(std::vector{pressure}); }, "normal_pressure"_a, py::call_guard(), "Solves the EP contact with an accelerated fixed-point scheme. May " "not converge!"); } } // namespace wrap } // namespace tamaas