diff --git a/python/py_akantu.cc b/python/py_akantu.cc index 071814a72..245d02bd6 100644 --- a/python/py_akantu.cc +++ b/python/py_akantu.cc @@ -1,169 +1,169 @@ /** * Copyright (©) 2018-2023 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * This file is part of Akantu * * Akantu 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. * * Akantu 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 Akantu. If not, see . */ /* -------------------------------------------------------------------------- */ #include "aka_config.hh" // for NLSNotConvergedException #include "non_linear_solver.hh" /* -------------------------------------------------------------------------- */ #include "py_aka_common.hh" #include "py_aka_error.hh" #include "py_boundary_conditions.hh" #include "py_constitutive_law.hh" #include "py_constitutive_law_selector.hh" #include "py_dof_manager.hh" #include "py_dumpable.hh" #include "py_fe_engine.hh" #include "py_group_manager.hh" #include "py_integration_scheme.hh" #include "py_mesh.hh" #include "py_model.hh" #include "py_parser.hh" #include "py_solver.hh" #if defined(AKANTU_SOLID_MECHANICS) #include "py_material.hh" #include "py_solid_mechanics_model.hh" #endif #if defined(AKANTU_DIFFUSION) #include "py_heat_transfer_model.hh" #endif #if defined(AKANTU_COHESIVE_ELEMENT) #include "py_fragment_manager.hh" #include "py_solid_mechanics_model_cohesive.hh" #endif #if defined(AKANTU_CONTACT_MECHANICS) #include "py_contact_mechanics_model.hh" #include "py_model_couplers.hh" #endif #if defined(AKANTU_PHASE_FIELD) #include "py_phase_field_model.hh" #endif #if defined(AKANTU_STRUCTURAL_MECHANICS) #include "py_structural_mechanics_model.hh" #endif /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ namespace py = pybind11; namespace akantu { void register_all(pybind11::module & mod) { register_initialize(mod); register_enums(mod); register_error(mod); register_functions(mod); register_parser(mod); register_solvers(mod); - register_group_manager(mod); register_dumpable(mod); + register_group_manager(mod); register_mesh(mod); register_fe_engine(mod); register_integration_schemes(mod); register_dof_manager(mod); register_boundary_conditions(mod); register_model(mod); register_constitutive_law_selector(mod); register_constitutive_law_internal_handler(mod); #if defined(AKANTU_DIFFUSION) register_heat_transfer_model(mod); #endif #if defined(AKANTU_SOLID_MECHANICS) register_solid_mechanics_model(mod); register_material(mod); #endif #if defined(AKANTU_COHESIVE_ELEMENT) register_solid_mechanics_model_cohesive(mod); register_fragment_manager(mod); #endif #if defined(AKANTU_STRUCTURAL_MECHANICS) register_structural_mechanics_model(mod); #endif #if defined(AKANTU_CONTACT_MECHANICS) register_contact_mechanics_model(mod); register_model_couplers(mod); #endif #if defined(AKANTU_PHASE_FIELD) register_phase_field_model(mod); register_phase_field_coupler(mod); #endif } } // namespace akantu /* -------------------------------------------------------------------------- */ /* -------------------------------------------------------------------------- */ PYBIND11_MODULE(py11_akantu, mod) { mod.doc() = "Akantu python interface"; static py::exception akantu_exception(mod, "Exception"); static py::exception akantu_exception_nls_not_converged(mod, "NLSNotConvergedException"); py::register_exception_translator([](std::exception_ptr ptr) { try { if (ptr) { std::rethrow_exception(ptr); } } catch (akantu::debug::NLSNotConvergedException & e) { akantu_exception_nls_not_converged(e.info().c_str()); } catch (akantu::debug::Exception & e) { if (akantu::debug::debugger.printBacktrace()) { akantu::debug::printBacktrace(); } akantu_exception(e.info().c_str()); } }); akantu::register_all(mod); mod.def("has_mpi", []() { #if defined(AKANTU_USE_MPI) return true; #else return false; #endif }) .def("getVersion", &akantu::getVersion); } // Module akantu diff --git a/python/py_fe_engine.cc b/python/py_fe_engine.cc index 4de1f9f0c..71cbcaeab 100644 --- a/python/py_fe_engine.cc +++ b/python/py_fe_engine.cc @@ -1,139 +1,211 @@ /** * Copyright (©) 2019-2023 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * This file is part of Akantu * * Akantu 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. * * Akantu 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 Akantu. If not, see . */ /* -------------------------------------------------------------------------- */ #include "py_aka_array.hh" #include "py_aka_common.hh" /* -------------------------------------------------------------------------- */ #include #include #include /* -------------------------------------------------------------------------- */ #include #include #include /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ namespace py = pybind11; /* -------------------------------------------------------------------------- */ namespace akantu { void register_fe_engine(py::module & mod) { py::class_(mod, "Element") .def(py::init([](ElementType type, Int id, GhostType ghost_type) { return std::make_unique(Element{type, id, ghost_type}); }), py::arg("type"), py::arg("ghost_type"), py::arg("ghost_type") = _not_ghost) .def("__lt__", [](Element & self, const Element & other) { return (self < other); }) .def("__repr__", [](Element & self) { return std::to_string(self); }); mod.attr("ElementNull") = ElementNull; py::class_(mod, "FEEngine") .def( "getNbIntegrationPoints", [](FEEngine & fem, ElementType type, GhostType ghost_type) { return fem.getNbIntegrationPoints(type, ghost_type); }, py::arg("type"), py::arg("ghost_type") = _not_ghost) .def("initShapeFunctions", &FEEngine::initShapeFunctions, py::arg("ghost_type") = _not_ghost) + + .def( + "integrate", + [](FEEngine & self, const Array & f, Array & intf, + Int nb_degree_of_freedom, ElementType type, + GhostType ghost_type = _not_ghost, + const Array & filter_elements = empty_filter) { + self.integrate(f, intf, nb_degree_of_freedom, type, ghost_type, + filter_elements); + }, + py::arg("f"), py::arg("intf"), py::arg("nb_degree_of_freedom"), + py::arg("type"), py::arg("ghost_type") = _not_ghost, + py::arg("filter_elements") = empty_filter) + .def( + "integrate", + [](FEEngine & self, const Array & f, ElementType type, + GhostType ghost_type = _not_ghost, + const Array & filter_elements = empty_filter) -> Real { + return self.integrate(f, type, ghost_type, filter_elements); + }, + py::arg("f"), py::arg("type"), py::arg("ghost_type") = _not_ghost, + py::arg("filter_elements") = empty_filter) .def( "gradientOnIntegrationPoints", [](FEEngine & fem, const Array & u, Array & nablauq, const Int nb_degree_of_freedom, ElementType type, GhostType ghost_type, const Array & filter_elements) { fem.gradientOnIntegrationPoints(u, nablauq, nb_degree_of_freedom, type, ghost_type, filter_elements); }, py::arg("u"), py::arg("nablauq"), py::arg("nb_degree_of_freedom"), py::arg("type"), py::arg("ghost_type") = _not_ghost, py::arg("filter_elements") = empty_filter) .def( "interpolateOnIntegrationPoints", [](FEEngine & self, const Array & u, Array & uq, Int nb_degree_of_freedom, ElementType type, GhostType ghost_type, const Array & filter_elements) { self.interpolateOnIntegrationPoints( u, uq, nb_degree_of_freedom, type, ghost_type, filter_elements); }, py::arg("u"), py::arg("uq"), py::arg("nb_degree_of_freedom"), py::arg("type"), py::arg("ghost_type") = _not_ghost, py::arg("filter_elements") = empty_filter) .def( "interpolateOnIntegrationPoints", [](FEEngine & self, const Array & u, std::shared_ptr> uq, std::shared_ptr> filter_elements) { self.interpolateOnIntegrationPoints(u, *uq, filter_elements.get()); }, py::arg("u"), py::arg("uq"), py::arg("filter_elements") = nullptr) + + .def( + "computeBtD", + [](FEEngine & self, const Array & Ds, Array & BtDs, + ElementType type, GhostType ghost_type = _not_ghost, + const Array & filter_elements = empty_filter) { + self.computeBtD(Ds, BtDs, type, ghost_type, filter_elements); + }, + py::arg("Ds"), py::arg("BtDs"), py::arg("type"), + py::arg("ghost_type") = _not_ghost, + py::arg("filter_elements") = empty_filter) + + .def( + "computeBtDB", + [](FEEngine & self, const Array & Ds, Array & BtDBs, + Int order_d, ElementType type, GhostType ghost_type = _not_ghost, + const Array & filter_elements = empty_filter) { + self.computeBtDB(Ds, BtDBs, order_d, type, ghost_type, + filter_elements); + }, + py::arg("Ds"), py::arg("BtDBs"), py::arg("order_d"), py::arg("type"), + py::arg("ghost_type") = _not_ghost, + py::arg("filter_elements") = empty_filter) + + .def( + "computeNtb", + [](FEEngine & self, const Array & bs, Array & Ntbs, + ElementType type, GhostType ghost_type = _not_ghost, + const Array & filter_elements = empty_filter) { + self.computeNtb(bs, Ntbs, type, ghost_type, filter_elements); + }, + py::arg("bs"), py::arg("Ntbs"), py::arg("type"), + py::arg("ghost_type") = _not_ghost, + py::arg("filter_elements") = empty_filter) + .def( + "computeNtbN", + [](FEEngine & self, const Array & bs, Array & NtbNs, + ElementType type, GhostType ghost_type = _not_ghost, + const Array & filter_elements = empty_filter) { + self.computeNtbN(bs, NtbNs, type, ghost_type, filter_elements); + }, + py::arg("bs"), py::arg("NtbNs"), py::arg("type"), + py::arg("ghost_type") = _not_ghost, + py::arg("filter_elements") = empty_filter) .def( "computeIntegrationPointsCoordinates", [](FEEngine & self, std::shared_ptr> coordinates, std::shared_ptr> filter_elements) -> decltype(auto) { return self.computeIntegrationPointsCoordinates( *coordinates, filter_elements.get()); }, py::arg("coordinates"), py::arg("filter_elements") = nullptr) .def( "assembleFieldLumped", [](FEEngine & fem, const std::function &, const Element &)> & field_funct, const ID & matrix_id, const ID & dof_id, DOFManager & dof_manager, ElementType type, GhostType ghost_type) { fem.assembleFieldLumped(field_funct, matrix_id, dof_id, dof_manager, type, ghost_type); }, py::arg("field_funct"), py::arg("matrix_id"), py::arg("dof_id"), py::arg("dof_manager"), py::arg("type"), py::arg("ghost_type") = _not_ghost) .def( "assembleFieldMatrix", [](FEEngine & fem, const std::function &, const Element &)> & field_funct, const ID & matrix_id, const ID & dof_id, DOFManager & dof_manager, ElementType type, GhostType ghost_type = _not_ghost) { fem.assembleFieldMatrix(field_funct, matrix_id, dof_id, dof_manager, type, ghost_type); }, py::arg("field_funct"), py::arg("matrix_id"), py::arg("dof_id"), py::arg("dof_manager"), py::arg("type"), py::arg("ghost_type") = _not_ghost) .def("getElementInradius", [](FEEngine & self, const Element & element) { return self.getElementInradius(element); }) .def("getNormalsOnIntegrationPoints", &FEEngine::getNormalsOnIntegrationPoints, py::arg("type"), py::arg("ghost_type") = _not_ghost, - py::return_value_policy::reference); + py::return_value_policy::reference) + .def("getShapes", &FEEngine::getShapes, py::arg("type"), + py::arg("ghost_type") = _not_ghost, py::arg("id") = 0, + py::return_value_policy::reference) + .def("getShapesDerivatives", &FEEngine::getShapesDerivatives, + py::arg("type"), py::arg("ghost_type") = _not_ghost, + py::arg("id") = 0, py::return_value_policy::reference); py::class_(mod, "IntegrationPoint"); } } // namespace akantu diff --git a/python/py_group_manager.cc b/python/py_group_manager.cc index a4d05bb73..5e43fef7a 100644 --- a/python/py_group_manager.cc +++ b/python/py_group_manager.cc @@ -1,172 +1,176 @@ /** * Copyright (©) 2019-2023 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * This file is part of Akantu * * Akantu 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. * * Akantu 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 Akantu. If not, see . */ /* -------------------------------------------------------------------------- */ #include "py_aka_array.hh" /* -------------------------------------------------------------------------- */ +#include #include +#include #include /* -------------------------------------------------------------------------- */ +#include +/* -------------------------------------------------------------------------- */ #include #include /* -------------------------------------------------------------------------- */ namespace py = pybind11; /* -------------------------------------------------------------------------- */ namespace akantu { /* -------------------------------------------------------------------------- */ void register_group_manager(py::module & mod) { /* ------------------------------------------------------------------------ */ py::class_(mod, "NodeGroup") .def( "getNodes", [](NodeGroup & self) -> decltype(auto) { return self.getNodes(); }, py::return_value_policy::reference) .def("__len__", &NodeGroup::size) .def( "__iter__", [](const NodeGroup & self) { return py::make_iterator(self.begin(), self.end()); }, py::keep_alive<0, 1>()) .def("__contains__", [](const NodeGroup & self, UInt node) { return self.find(node) != -1; }) .def("getName", &NodeGroup::getName) .def("clear", &NodeGroup::clear) .def("empty", &NodeGroup::empty) .def("append", &NodeGroup::append) .def( "add", [](NodeGroup & self, UInt node, bool check_for_duplicate) { auto && it = self.add(node, check_for_duplicate); return *it; }, py::arg("node"), py::arg("check_for_duplicate") = true) .def("remove", &NodeGroup::add); /* ------------------------------------------------------------------------ */ - py::class_(mod, "ElementGroup") + py::class_(mod, "ElementGroup") .def( "getNodeGroup", [](ElementGroup & self) -> decltype(auto) { return self.getNodeGroup(); }, py::return_value_policy::reference) .def("getName", &ElementGroup::getName) .def( "getElements", [](ElementGroup & self) -> decltype(auto) { return self.getElements(); }, py::return_value_policy::reference) .def( "getNodeGroup", [](ElementGroup & self) -> decltype(auto) { return self.getNodeGroup(); }, py::return_value_policy::reference) .def("__len__", [](const ElementGroup & self) { return self.size(); }) .def("clear", [](ElementGroup & self) { self.clear(); }) .def("empty", &ElementGroup::empty) .def("append", &ElementGroup::append) .def("optimize", &ElementGroup::optimize) .def( "add", [](ElementGroup & self, const Element & element, bool add_nodes, bool check_for_duplicate) { self.add(element, add_nodes, check_for_duplicate); }, py::arg("element"), py::arg("add_nodes") = false, py::arg("check_for_duplicate") = true) .def("fillFromNodeGroup", &ElementGroup::fillFromNodeGroup) .def("addDimension", &ElementGroup::addDimension); /* ------------------------------------------------------------------------ */ py::class_(mod, "GroupManager") .def( "getElementGroup", [](GroupManager & self, const std::string & name) -> decltype(auto) { return self.getElementGroup(name); }, py::return_value_policy::reference) .def("iterateElementGroups", [](GroupManager & self) -> decltype(auto) { std::vector> groups; for (auto & group : self.iterateElementGroups()) { groups.emplace_back(group); } return groups; }) .def("iterateNodeGroups", [](GroupManager & self) -> decltype(auto) { std::vector> groups; for (auto & group : self.iterateNodeGroups()) { groups.emplace_back(group); } return groups; }) .def("createNodeGroup", &GroupManager::createNodeGroup, py::return_value_policy::reference) .def( "createElementGroup", [](GroupManager & self, const std::string & id, Int spatial_dimension, bool replace_group) -> decltype(auto) { return self.createElementGroup(id, spatial_dimension, replace_group); }, py::arg("id"), py::arg("spatial_dimension"), py::arg("replace_group") = false, py::return_value_policy::reference) .def("createGroupsFromMeshDataUInt", &GroupManager::createGroupsFromMeshData) .def("createElementGroupFromNodeGroup", &GroupManager::createElementGroupFromNodeGroup, py::arg("name"), py::arg("node_group"), py::arg("dimension") = _all_dimensions) .def( "getNodeGroup", [](GroupManager & self, const std::string & name) -> decltype(auto) { return self.getNodeGroup(name); }, py::return_value_policy::reference) .def( "nodeGroups", [](GroupManager & self) { std::vector groups; for (auto & g : self.iterateNodeGroups()) { groups.push_back(&g); } return groups; }, py::return_value_policy::reference) .def( "elementGroups", [](GroupManager & self) { std::vector groups; for (auto & g : self.iterateElementGroups()) { groups.push_back(&g); } return groups; }, py::return_value_policy::reference) .def("createBoundaryGroupFromGeometry", &GroupManager::createBoundaryGroupFromGeometry); } } // namespace akantu diff --git a/python/py_model.cc b/python/py_model.cc index 21c9a7544..8ea74cd06 100644 --- a/python/py_model.cc +++ b/python/py_model.cc @@ -1,202 +1,238 @@ /** * Copyright (©) 2019-2023 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * This file is part of Akantu * * Akantu 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. * * Akantu 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 Akantu. If not, see . */ /* -------------------------------------------------------------------------- */ #include "py_aka_array.hh" /* -------------------------------------------------------------------------- */ #include #include #include #include #include #include /* -------------------------------------------------------------------------- */ #include #include #include /* -------------------------------------------------------------------------- */ namespace py = pybind11; /* -------------------------------------------------------------------------- */ namespace akantu { template void register_constitutive_law_handler(py::module & mod, const std::string & name) { py::class_, ModelType>(mod, name.c_str(), py::multiple_inheritance()) .def( "getConstitutiveLaw", [](ConstitutiveLawsHandler & self, UInt constitutive_law_id) -> decltype(auto) { return self.getConstitutiveLaw(constitutive_law_id); }, py::arg("constitutive_law_id"), py::return_value_policy::reference) .def( "getConstitutiveLaw", [](ConstitutiveLawsHandler & self, const ID & constitutive_law_name) -> decltype(auto) { return self.getConstitutiveLaw(constitutive_law_name); }, py::arg("constitutive_law_name"), py::return_value_policy::reference) .def( "getConstitutiveLaw", [](ConstitutiveLawsHandler & self, const Element & element) -> decltype(auto) { return self.getConstitutiveLaw(element); }, py::arg("element"), py::return_value_policy::reference) .def("getNbConstitutiveLaws", &ConstitutiveLawsHandler::getNbConstitutiveLaws) .def("getConstitutiveLawIndex", &ConstitutiveLawsHandler::getConstitutiveLawIndex) .def("setConstitutiveLawSelector", [](ConstitutiveLawsHandler & self, std::shared_ptr constitutive_law_selector) { self.setConstitutiveLawSelector(constitutive_law_selector); }) .def("getConstitutiveLawSelector", &ConstitutiveLawsHandler::getConstitutiveLawSelector) .def( "getConstitutiveLawByElement", [](const ConstitutiveLawsHandler & self) -> decltype(auto) { return self.getConstitutiveLawByElement(); }, py::return_value_policy::reference, py::keep_alive<0, 1>()) .def("reassignConstitutiveLaw", &ConstitutiveLawsHandler::reassignConstitutiveLaw) .def( "registerNewConstitutiveLaw", [](ConstitutiveLawsHandler & self, const ID & cl_name, const ID & cl_type, const ID & opt_param) -> decltype(auto) { return self.registerNewConstitutiveLaw(cl_name, cl_type, opt_param); }, py::arg("constitutive_law_name"), py::arg("constitutive_law_type"), py::arg("option") = "", py::return_value_policy::reference) .def("initConstitutiveLaws", &ConstitutiveLawsHandler::initConstitutiveLaws) .def("flattenInternal", &ConstitutiveLawsHandler::flattenInternal, py::return_value_policy::reference) .def("inflateInternal", &ConstitutiveLawsHandler::inflateInternal, py::return_value_policy::reference); } /* -------------------------------------------------------------------------- */ void register_model(py::module & mod) { py::class_(mod, "ModelSolver", py::multiple_inheritance()) .def( "getNonLinearSolver", [](ModelSolver & self, const ID & solver_id) -> NonLinearSolver & { return self.getNonLinearSolver(solver_id); }, py::arg("solver_id") = "", py::return_value_policy::reference) .def( "getTimeStepSolver", [](ModelSolver & self, const ID & solver_id) -> TimeStepSolver & { return self.getTimeStepSolver(solver_id); }, py::arg("solver_id") = "", py::return_value_policy::reference) .def( "solveStep", [](ModelSolver & self, const ID & solver_id) { self.solveStep(solver_id); }, py::arg("solver_id") = "") .def( "solveStep", [](ModelSolver & self, SolverCallback & callback, const ID & solver_id) { self.solveStep(callback, solver_id); }, py::arg("callback"), py::arg("solver_id") = ""); py::class_(mod, "Model", py::multiple_inheritance()) .def("setBaseName", &Model::setBaseName) .def("setDirectory", &Model::setDirectory) .def("getFEEngine", &Model::getFEEngine, py::arg("name") = "", py::return_value_policy::reference) .def("getFEEngineBoundary", &Model::getFEEngine, py::arg("name") = "", py::return_value_policy::reference) .def("addDumpFieldVector", &Model::addDumpFieldVector) .def("addDumpField", &Model::addDumpField) .def("setBaseNameToDumper", &Model::setBaseNameToDumper) .def("addDumpFieldVectorToDumper", &Model::addDumpFieldVectorToDumper) .def("addDumpFieldToDumper", &Model::addDumpFieldToDumper) .def("dump", [](Model & self) { self.dump(); }) .def( "dump", [](Model & self, UInt step) { self.dump(step); }, py::arg("step")) .def( "dump", [](Model & self, Real time, UInt step) { self.dump(time, step); }, py::arg("time"), py::arg("step")) .def( "dump", [](Model & self, const std::string & dumper) { self.dump(dumper); }, py::arg("dumper_name")) .def( "dump", [](Model & self, const std::string & dumper, UInt step) { self.dump(dumper, step); }, py::arg("dumper_name"), py::arg("step")) .def( "dump", [](Model & self, const std::string & dumper, Real time, UInt step) { self.dump(dumper, time, step); }, py::arg("dumper_name"), py::arg("time"), py::arg("step")) + + .def( + "dumpGroup", + [](Model & self, const std::string & group_name) { + self.dumpGroup(group_name); + }, + py::arg("group_name")) + .def( + "setGroupDirectory", + [](Model & self, const std::string & directory, + const std::string & group_name) { + self.setGroupDirectory(directory, group_name); + }, + py::arg("directory"), py::arg("group_name")) + .def( + "setGroupBaseName", + [](Model & self, const std::string & basename, + const std::string & group_name) { + self.setGroupBaseName(basename, group_name); + }, + py::arg("basename"), py::arg("group_name")) + .def( + "addDumpGroupField", + [](Model & self, const std::string & field_id, + const std::string & group_name) { + self.addDumpGroupField(field_id, group_name); + }, + py::arg("field_id"), py::arg("group_name")) + + .def( + "addDumpGroupFieldVector", + [](Model & self, const std::string & field_id, + const std::string & group_name) { + self.addDumpGroupFieldVector(field_id, group_name); + }, + py::arg("field_id"), py::arg("group_name")) .def("initNewSolver", &Model::initNewSolver) .def( "getNewSolver", [](Model & self, const std::string id, const TimeStepSolverType & time, const NonLinearSolverType & type) { self.getNewSolver(id, time, type); }, py::return_value_policy::reference) .def( "setIntegrationScheme", [](Model & self, const std::string id, const std::string primal, const IntegrationSchemeType & scheme_type, IntegrationScheme::SolutionType solution_type) { self.setIntegrationScheme(id, primal, scheme_type, solution_type); }, py::arg("id"), py::arg("primal"), py::arg("scheme_type"), py::arg("solution_type") = IntegrationScheme::SolutionType::_not_defined) .def("getDOFManager", &Model::getDOFManager, py::return_value_policy::reference) .def("assembleMatrix", &Model::assembleMatrix); } } // namespace akantu diff --git a/src/fe_engine/fe_engine.hh b/src/fe_engine/fe_engine.hh index 24bc43b8c..ee6fa7f90 100644 --- a/src/fe_engine/fe_engine.hh +++ b/src/fe_engine/fe_engine.hh @@ -1,370 +1,363 @@ /** * Copyright (©) 2010-2023 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * This file is part of Akantu * * Akantu 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. * * Akantu 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 Akantu. If not, see . */ /* -------------------------------------------------------------------------- */ #include "element_type_map.hh" #include "mesh_events.hh" /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ #ifndef AKANTU_FE_ENGINE_HH_ #define AKANTU_FE_ENGINE_HH_ namespace akantu { class Mesh; class Integrator; class ShapeFunctions; class DOFManager; class Element; } // namespace akantu /* -------------------------------------------------------------------------- */ namespace akantu { /* -------------------------------------------------------------------------- */ /** * The generic FEEngine class derived in a FEEngineTemplate class * containing the * shape functions and the integration method */ class FEEngine : public MeshEventHandler { /* ------------------------------------------------------------------------ */ /* Constructors/Destructors */ /* ------------------------------------------------------------------------ */ public: FEEngine(Mesh & mesh, Int element_dimension = _all_dimensions, const ID & id = "fem"); ~FEEngine() override; /* ------------------------------------------------------------------------ */ /* Methods */ /* ------------------------------------------------------------------------ */ public: /// pre-compute all the shape functions, their derivatives and the jacobians virtual void initShapeFunctions(GhostType ghost_type = _not_ghost) = 0; /// extract the nodal values and store them per element template static void extractNodalToElementField(const Mesh & mesh, const Array & nodal_f, Array & elemental_f, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter); /// filter a field template static void filterElementalData(const Mesh & mesh, const Array & quad_f, Array & filtered_f, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter); /* ------------------------------------------------------------------------ */ /* Integration method bridges */ /* ------------------------------------------------------------------------ */ /// integrate f for all elements of type "type" virtual void integrate(const Array & f, Array & intf, Int nb_degree_of_freedom, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter) const = 0; /// integrate a scalar value f on all elements of type "type" [[nodiscard]] virtual Real integrate(const Array & f, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter) const = 0; - /// integrate f for all integration points of type "type" but don't sum over - /// all integration points - virtual void integrateOnIntegrationPoints( - const Array & f, Array & intf, Int nb_degree_of_freedom, - ElementType type, GhostType ghost_type = _not_ghost, - const Array & filter_elements = empty_filter) const = 0; - /// integrate one element scalar value on all elements of type "type" [[nodiscard]] Real integrate(const Ref f, const Element & element) const { return integrate(f, element.type, element.element, element.ghost_type); } private: [[nodiscard]] virtual Real integrate(const Ref f, ElementType type, Idx index, GhostType ghost_type = _not_ghost) const = 0; /* ------------------------------------------------------------------------ */ /* compatibility with old FEEngine fashion */ /* ------------------------------------------------------------------------ */ public: /// get the number of integration points [[nodiscard]] virtual Int getNbIntegrationPoints(ElementType type, GhostType ghost_type = _not_ghost) const = 0; /// get the precomputed shapes [[nodiscard]] const virtual Array & getShapes(ElementType type, GhostType ghost_type = _not_ghost, Idx id = 0) const = 0; /// get the derivatives of shapes [[nodiscard]] virtual const Array & getShapesDerivatives(ElementType type, GhostType ghost_type = _not_ghost, Idx id = 0) const = 0; /// get integration points [[nodiscard]] virtual const MatrixXr & getIntegrationPoints(ElementType type, GhostType ghost_type = _not_ghost) const = 0; /* ------------------------------------------------------------------------ */ /* Shape method bridges */ /* ------------------------------------------------------------------------ */ /// Compute the gradient nablauq on the integration points of an element type /// from nodal values u virtual void gradientOnIntegrationPoints( const Array & u, Array & nablauq, Int nb_degree_of_freedom, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter) const = 0; /// Interpolate a nodal field u at the integration points of an element type /// -> uq virtual void interpolateOnIntegrationPoints( const Array & u, Array & uq, Int nb_degree_of_freedom, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter) const = 0; /// Interpolate a nodal field u at the integration points of many element /// types -> uq virtual void interpolateOnIntegrationPoints( const Array & u, ElementTypeMapArray & uq, const ElementTypeMapArray * filter_elements = nullptr) const = 0; /// pre multiplies a tensor by the shapes derivaties virtual void computeBtD(const Array & Ds, Array & BtDs, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter) const = 0; /// left and right multiplies a tensor by the shapes derivaties virtual void computeBtDB(const Array & Ds, Array & BtDBs, Int order_d, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter) const = 0; /// left multiples a vector by the shape functions virtual void computeNtb(const Array & bs, Array & Ntbs, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter) const = 0; /// left and right multiplies a tensor by the shapes virtual void computeNtbN(const Array & bs, Array & NtbNs, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter) const = 0; /// Compute the interpolation point position in the global coordinates for /// many element types virtual void computeIntegrationPointsCoordinates( ElementTypeMapArray & integration_points_coordinates, const ElementTypeMapArray * filter_elements = nullptr) const = 0; /// Compute the interpolation point position in the global coordinates for an /// element type virtual void computeIntegrationPointsCoordinates( Array & integration_points_coordinates, ElementType type, GhostType ghost_type = _not_ghost, const Array & filter_elements = empty_filter) const = 0; /// Build pre-computed matrices for interpolation of field form integration /// points at other given positions (interpolation_points) virtual void initElementalFieldInterpolationFromIntegrationPoints( const ElementTypeMapArray & interpolation_points_coordinates, ElementTypeMapArray & interpolation_points_coordinates_matrices, ElementTypeMapArray & integration_points_coordinates_inv_matrices, const ElementTypeMapArray * element_filter) const = 0; /// interpolate field at given position (interpolation_points) from given /// values of this field at integration points (field) virtual void interpolateElementalFieldFromIntegrationPoints( const ElementTypeMapArray & field, const ElementTypeMapArray & interpolation_points_coordinates, ElementTypeMapArray & result, const GhostType ghost_type, const ElementTypeMapArray * element_filter) const = 0; /// Interpolate field at given position from given values of this field at /// integration points (field) /// using matrices precomputed with /// initElementalFieldInterplationFromIntegrationPoints virtual void interpolateElementalFieldFromIntegrationPoints( const ElementTypeMapArray & field, const ElementTypeMapArray & interpolation_points_coordinates_matrices, const ElementTypeMapArray & integration_points_coordinates_inv_matrices, ElementTypeMapArray & result, const GhostType ghost_type, const ElementTypeMapArray * element_filter) const = 0; /// interpolate on a phyiscal point inside an element virtual void interpolate(const Ref real_coords, const Ref nodal_values, Ref interpolated, const Element & element) const = 0; /// compute the shape on a provided point virtual void computeShapes(const Ref real_coords, Int elem, ElementType type, Ref shapes, GhostType ghost_type = _not_ghost) const = 0; /// compute the shape derivatives on a provided point virtual void computeShapeDerivatives(const Ref real_coords, Int element, ElementType type, Ref shape_derivatives, GhostType ghost_type = _not_ghost) const = 0; /// assembles the lumped version of @f[ \int N^t rho N @f] virtual void assembleFieldLumped( const std::function &, const Element &)> & field_funct, const ID & matrix_id, const ID & dof_id, DOFManager & dof_manager, ElementType type, GhostType ghost_type = _not_ghost) const = 0; /// assembles the matrix @f[ \int N^t rho N @f] virtual void assembleFieldMatrix( const std::function &, const Element &)> & field_funct, const ID & matrix_id, const ID & dof_id, DOFManager & dof_manager, ElementType type, GhostType ghost_type = _not_ghost) const = 0; /* ------------------------------------------------------------------------ */ /* Other methods */ /* ------------------------------------------------------------------------ */ /// pre-compute normals on integration points virtual void computeNormalsOnIntegrationPoints(GhostType ghost_type = _not_ghost) = 0; /// pre-compute normals on integration points virtual void computeNormalsOnIntegrationPoints(const Array & /*field*/, GhostType /*ghost_type*/ = _not_ghost) { AKANTU_TO_IMPLEMENT(); } /// pre-compute normals on integration points virtual void computeNormalsOnIntegrationPoints( const Array & /*field*/, Array & /*normal*/, ElementType /*type*/, GhostType /*ghost_type*/ = _not_ghost) const { AKANTU_TO_IMPLEMENT(); } /// function to print the containt of the class virtual void printself(std::ostream & stream, int indent = 0) const; private: /// initialise the class void init(); /* ------------------------------------------------------------------------ */ /* Accessors */ /* ------------------------------------------------------------------------ */ public: using ElementTypesIteratorHelper = ElementTypeMapArray::ElementTypesIteratorHelper; [[nodiscard]] ElementTypesIteratorHelper elementTypes(Int dim = _all_dimensions, GhostType ghost_type = _not_ghost, ElementKind kind = _ek_regular) const; /// get the dimension of the element handeled by this fe_engine object AKANTU_GET_MACRO_AUTO(ElementDimension, element_dimension); /// get the mesh contained in the fem object AKANTU_GET_MACRO_AUTO(Mesh, mesh); /// get the mesh contained in the fem object AKANTU_GET_MACRO_NOT_CONST(Mesh, mesh, Mesh &); /// get the in-radius of an element template [[nodiscard]] static inline Real getElementInradius(const Eigen::MatrixBase & coord, ElementType type); [[nodiscard]] inline Real getElementInradius(const Element & element) const; /// get the normals on integration points AKANTU_GET_MACRO_BY_ELEMENT_TYPE_CONST(NormalsOnIntegrationPoints, normals_on_integration_points, Real); /// get cohesive element type for a given facet type [[nodiscard]] static inline constexpr ElementType getCohesiveElementType(ElementType type_facet); /// get igfem element type for a given regular type [[nodiscard]] static inline Vector getIGFEMElementTypes(ElementType type); /// get the interpolation element associated to an element type [[nodiscard]] static inline constexpr auto getInterpolationType(ElementType el_type); /// get the shape function class (probably useless: see getShapeFunction in /// fe_engine_template.hh) [[nodiscard]] virtual const ShapeFunctions & getShapeFunctionsInterface() const = 0; /// get the integrator class (probably useless: see getIntegrator in /// fe_engine_template.hh) [[nodiscard]] virtual const Integrator & getIntegratorInterface() const = 0; AKANTU_GET_MACRO(ID, id, ID); /* ------------------------------------------------------------------------ */ /* Class Members */ /* ------------------------------------------------------------------------ */ protected: ID id; /// spatial dimension of the problem Int element_dimension; /// the mesh on which all computation are made Mesh & mesh; /// normals at integration points ElementTypeMapArray normals_on_integration_points; }; /* -------------------------------------------------------------------------- */ /* inline functions */ /* -------------------------------------------------------------------------- */ /// standard output stream operator inline std::ostream & operator<<(std::ostream & stream, const FEEngine & _this) { _this.printself(stream); return stream; } } // namespace akantu #include "fe_engine_inline_impl.hh" #include "fe_engine_template.hh" #endif /* AKANTU_FE_ENGINE_HH_ */ diff --git a/src/fe_engine/fe_engine_template.hh b/src/fe_engine/fe_engine_template.hh index b53502f77..a442e0ee7 100644 --- a/src/fe_engine/fe_engine_template.hh +++ b/src/fe_engine/fe_engine_template.hh @@ -1,481 +1,474 @@ /** * Copyright (©) 2010-2023 EPFL (Ecole Polytechnique Fédérale de Lausanne) * Laboratory (LSMS - Laboratoire de Simulation en Mécanique des Solides) * * This file is part of Akantu * * Akantu 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. * * Akantu 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 Akantu. If not, see . */ /* -------------------------------------------------------------------------- */ #include "fe_engine.hh" /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ /* -------------------------------------------------------------------------- */ #ifndef AKANTU_FE_ENGINE_TEMPLATE_HH_ #define AKANTU_FE_ENGINE_TEMPLATE_HH_ namespace akantu { class Integrator; class ShapeFunctions; } // namespace akantu namespace akantu { class DOFManager; namespace fe_engine { namespace details { template struct AssembleLumpedTemplateHelper; template struct AssembleFieldMatrixHelper; } // namespace details } // namespace fe_engine template struct AssembleFieldMatrixStructHelper; struct DefaultIntegrationOrderFunctor { template static inline constexpr int getOrder() { return ElementClassProperty::polynomial_degree; } }; /* -------------------------------------------------------------------------- */ template