diff --git a/python/wrap/solvers.cpp b/python/wrap/solvers.cpp index a5dbb2e..6ef48a7 100644 --- a/python/wrap/solvers.cpp +++ b/python/wrap/solvers.cpp @@ -1,73 +1,82 @@ /** * @file * * @author Lucas Frérot * * @section LICENSE * * Copyright (©) 2017 EPFL (Ecole Polytechnique Fédérale de * Lausanne) Laboratory (LSMS - Laboratoire de Simulation en Mécanique des * Solides) * * Tamaas is free software: you can redistribute it and/or modify it under the * terms of the GNU Lesser General Public License as published by the Free * Software Foundation, either version 3 of the License, or (at your option) any * later version. * * Tamaas 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 Lesser General Public License for more * details. * * You should have received a copy of the GNU Lesser General Public License * along with Tamaas. If not, see . * */ /* -------------------------------------------------------------------------- */ #include "wrap.hh" #include "contact_solver.hh" #include "polonsky_keer_rey.hh" #include "kato.hh" +#include "beck_teboulle.hh" /* -------------------------------------------------------------------------- */ __BEGIN_TAMAAS__ namespace wrap { /* -------------------------------------------------------------------------- */ using namespace py::literals; /* -------------------------------------------------------------------------- */ void wrapSolvers(py::module& mod) { py::class_(mod, "ContactSolver") .def(py::init&, Real>(), "model"_a, "surface"_a, "tolerance"_a) .def("setMaxIterations", &ContactSolver::setMaxIterations, "max_iter"_a) .def("setDumpFrequency", &ContactSolver::setDumpFrequency, "dump_freq"_a) .def("addFunctionalTerm", &ContactSolver::addFunctionalTerm); py::class_ pkr(mod, "PolonskyKeerRey"); pkr.def(py::init&, Real, PolonskyKeerRey::type, PolonskyKeerRey::type>(), "model"_a, "surface"_a, "tolerance"_a, "primal_type"_a, "constraint_type"_a, py::keep_alive<1, 2>(), py::keep_alive<1, 3>()) .def("solve", &PolonskyKeerRey::solve, "target"_a) .def("computeError", &PolonskyKeerRey::computeError); py::enum_(pkr, "type") .value("gap", PolonskyKeerRey::gap) .value("pressure", PolonskyKeerRey::pressure) .export_values(); 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) .def("solveRelaxed", &Kato::solveRelaxed, "g0"_a) .def("computeCost", &Kato::computeCost); + + 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("solveRelaxed", &BeckTeboulle::solveRelaxed, "g0"_a) + .def("computeCost", &BeckTeboulle::computeCost); } } // namespace wrap __END_TAMAAS__ diff --git a/src/SConscript b/src/SConscript index b4e6fec..eef18ae 100644 --- a/src/SConscript +++ b/src/SConscript @@ -1,115 +1,116 @@ import os Import('main_env') def prepend(path, list): return [os.path.join(path, x) for x in list] env = main_env # Core core_list = """ fft_plan_manager.cpp fftransform.cpp fftransform_fftw.cpp grid.cpp grid_hermitian.cpp statistics.cpp surface.cpp tamaas.cpp types.cpp loop.cpp """.split() core_list = prepend('core', core_list) core_list.append('tamaas_info.cpp') # Lib roughcontact generator_list = """ surface_generator.cpp surface_generator_filter.cpp surface_generator_filter_fft.cpp isopowerlaw.cpp """.split() generator_list = prepend('surface', generator_list) # Lib PERCOLATION percolation_list = """ flood_fill.cpp """.split() percolation_list = prepend('percolation', percolation_list) # BEM PERCOLATION bem_list = """ bem_kato.cpp bem_polonski.cpp bem_gigi.cpp bem_gigipol.cpp bem_penalty.cpp bem_uzawa.cpp bem_fft_base.cpp bem_functional.cpp bem_meta_functional.cpp elastic_energy_functional.cpp exponential_adhesion_functional.cpp squared_exponential_adhesion_functional.cpp maugis_adhesion_functional.cpp complimentary_term_functional.cpp bem_grid.cpp bem_grid_polonski.cpp bem_grid_kato.cpp bem_grid_teboulle.cpp bem_grid_condat.cpp """.split() bem_list = prepend('bem', bem_list) # Model model_list = """ model.cpp model_factory.cpp model_type.cpp model_template.cpp westergaard.cpp elastic_functional.cpp meta_functional.cpp adhesion_functional.cpp elasto_plastic_model.cpp ve_engine.cpp kelvin.cpp mindlin.cpp """.split() model_list = prepend('model', model_list) # Solvers solvers_list = """ contact_solver.cpp polonsky_keer_rey.cpp kato.cpp +beck_teboulle.cpp """.split() solvers_list = prepend('solvers', solvers_list) # GPU API gpu_list = """ fftransform_cufft.cpp """.split() gpu_list = prepend('gpu', gpu_list) # Assembling total list rough_contact_list = \ core_list + generator_list + percolation_list + model_list + solvers_list + bem_list # Adding GPU if needed if env['backend'] == 'cuda': rough_contact_list += gpu_list # Generating dependency files # env.AppendUnique(CXXFLAGS=['-MMD']) # for file_name in rough_contact_list: # obj = file_name.replace('.cpp', '.os') # dep_file_name = file_name.replace('.cpp', '.d') # env.SideEffect(dep_file_name, obj) # env.ParseDepends(dep_file_name) libTamaas = env.SharedLibrary('Tamaas', rough_contact_list) Export('libTamaas') diff --git a/src/solvers/beck_teboulle.cpp b/src/solvers/beck_teboulle.cpp new file mode 100644 index 0000000..4881770 --- /dev/null +++ b/src/solvers/beck_teboulle.cpp @@ -0,0 +1,86 @@ +/** + * @file + * + * @author Son Pham-Ba + * + * @section LICENSE + * + * Copyright (©) 2016-2018 EPFL (Ecole Polytechnique Fédérale de + * Lausanne) Laboratory (LSMS - Laboratoire de Simulation en Mécanique des + * Solides) + * + * Tamaas is free software: you can redistribute it and/or modify it under the + * terms of the GNU Lesser General Public License as published by the Free + * Software Foundation, either version 3 of the License, or (at your option) any + * later version. + * + * Tamaas 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 Lesser General Public License for more + * details. + * + * You should have received a copy of the GNU Lesser General Public License + * along with Tamaas. If not, see . + * + */ + +/* -------------------------------------------------------------------------- */ +#include "beck_teboulle.hh" +#include +/* -------------------------------------------------------------------------- */ + +__BEGIN_TAMAAS__ + +BeckTeboulle::BeckTeboulle(Model& model, const GridBase& surface, + Real tolerance, Real mu) + : Kato(model, surface, tolerance, mu) {} + +/* -------------------------------------------------------------------------- */ + +Real BeckTeboulle::solve(GridBase& g0) { + Real cost = 0; + UInt n = 0; + + // Printing column headers + std::cout << std::setw(5) << "Iter" + << " " << std::setw(15) << "Cost_f" + << " " << std::setw(15) << "Error" << '\n' + << std::fixed; + + *pressure = 0; + + do { + engine.solveNeumann(*pressure, *gap); + (model.getType() == model_type::surface_1d) ? + addUniform<2>(*gap, g0) : + addUniform<3>(*gap, g0); + *gap -= surfaceComp; + *pressure -= *gap; + + (model.getType() == model_type::surface_1d) ? + enforcePressureCoulomb() : + enforcePressureCoulomb(); + cost = computeCost(); + printState(n, cost, cost); + } while (cost > this->tolerance && n++ < this->max_iterations); + + Real beta = 0; + switch (model.getType()) { + case model_type::surface_1d: + beta = computeBeta(); + computeFinalGap(beta); + break; + case model_type::surface_2d: + beta = computeBeta(); + computeFinalGap(beta); + break; + default: + break; + } + model.getDisplacement() = *gap; + + return cost; +} + +__END_TAMAAS__ +/* -------------------------------------------------------------------------- */ diff --git a/src/solvers/beck_teboulle.hh b/src/solvers/beck_teboulle.hh new file mode 100644 index 0000000..61ebdd2 --- /dev/null +++ b/src/solvers/beck_teboulle.hh @@ -0,0 +1,50 @@ +/** + * @file + * + * @author Son Pham-Ba + * + * @section LICENSE + * + * Copyright (©) 2016-2018 EPFL (Ecole Polytechnique Fédérale de + * Lausanne) Laboratory (LSMS - Laboratoire de Simulation en Mécanique des + * Solides) + * + * Tamaas is free software: you can redistribute it and/or modify it under the + * terms of the GNU Lesser General Public License as published by the Free + * Software Foundation, either version 3 of the License, or (at your option) any + * later version. + * + * Tamaas 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 Lesser General Public License for more + * details. + * + * You should have received a copy of the GNU Lesser General Public License + * along with Tamaas. If not, see . + * + */ +/* -------------------------------------------------------------------------- */ +#ifndef __BECK_TEBOULLE_HH__ +#define __BECK_TEBOULLE_HH__ +/* -------------------------------------------------------------------------- */ +#include "kato.hh" +/* -------------------------------------------------------------------------- */ +__BEGIN_TAMAAS__ + +class BeckTeboulle : public Kato { +public: + /// Constructor + BeckTeboulle(Model& model, const GridBase& surface, Real tolerance, + Real mu); + +public: + /// Solve + Real solve(GridBase& p0); + +private: + +}; + +__END_TAMAAS__ + +#endif // __BECK_TEBOULLE_HH__ diff --git a/src/solvers/kato.hh b/src/solvers/kato.hh index 78fce6a..b1571bb 100644 --- a/src/solvers/kato.hh +++ b/src/solvers/kato.hh @@ -1,95 +1,95 @@ /** * @file * * @author Son Pham-Ba * * @section LICENSE * * Copyright (©) 2016-2018 EPFL (Ecole Polytechnique Fédérale de * Lausanne) Laboratory (LSMS - Laboratoire de Simulation en Mécanique des * Solides) * * Tamaas is free software: you can redistribute it and/or modify it under the * terms of the GNU Lesser General Public License as published by the Free * Software Foundation, either version 3 of the License, or (at your option) any * later version. * * Tamaas 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 Lesser General Public License for more * details. * * You should have received a copy of the GNU Lesser General Public License * along with Tamaas. If not, see . * */ /* -------------------------------------------------------------------------- */ -#ifndef __KATO__ -#define __KATO__ +#ifndef __KATO_HH__ +#define __KATO_HH__ /* -------------------------------------------------------------------------- */ #include "contact_solver.hh" #include "meta_functional.hh" #include "model_type.hh" #include "static_types.hh" #include "tamaas.hh" /* -------------------------------------------------------------------------- */ __BEGIN_TAMAAS__ class Kato : public ContactSolver { public: /// Constructor Kato(Model& model, const GridBase& surface, Real tolerance, Real mu); public: /// Solve Real solve(GridBase& p0); /// Solve relaxed problem Real solveRelaxed(GridBase& g0); /// Compute cost function Real computeCost(); -private: +protected: /// Creates surfaceComp form surface template void initSurfaceWithComponents(); /// Project pressure on friction cone template void enforcePressureConstraints(GridBase& p0); /// Project on C template void enforcePressureMean(GridBase& p0); /// Project on D template void enforcePressureCoulomb(); /// Comupte mean value of field template void computeMean(GridBase& field, GridBase& mean); /// Add vector to each point of field template void addUniform(GridBase& field, GridBase& vec); /// Compute mean gap template Real computeBeta(); /// Compute grids of dual and primal variables template void computeValuesForCost(Real beta, GridBase& lambda, GridBase& eta, GridBase& p_N, GridBase& p_C); /// Compute total displacement template void computeFinalGap(Real beta); protected: BEEngine& engine; GridBase* gap = nullptr; GridBase* pressure = nullptr; Grid surfaceComp; Real mu = 0; UInt comp = 0; // number of components UInt N = 0; // number of points }; __END_TAMAAS__ -#endif // __KATO__ +#endif // __KATO_HH__