diff --git a/python/akantu/__init__.py b/python/akantu/__init__.py index 15a981286..23be37500 100644 --- a/python/akantu/__init__.py +++ b/python/akantu/__init__.py @@ -1,19 +1,49 @@ +import scipy.sparse as _aka_sparse +import numpy as _aka_np import py11_akantu private_keys = set(dir(py11_akantu)) - set(dir()) for k in private_keys: globals()[k] = getattr(py11_akantu, k) def initialize(*args, **kwargs): raise RuntimeError("No need to call initialize," " use parseInput to read an input file") def finalize(*args, **kwargs): raise RuntimeError("No need to call finalize") +class AkantuSparseMatrix (_aka_sparse.coo_matrix): + + def __init__(self, aka_sparse): + + self.aka_sparse = aka_sparse + matrix_type = self.aka_sparse.getMatrixType() + sz = self.aka_sparse.size() + row = self.aka_sparse.getIRN()[:, 0] - 1 + col = self.aka_sparse.getJCN()[:, 0] - 1 + data = self.aka_sparse.getA()[:, 0] + + row = row.copy() + col = col.copy() + data = data.copy() + + if matrix_type == py11_akantu._symmetric: + non_diags = (row != col) + row_sup = col[non_diags] + col_sup = row[non_diags] + data_sup = data[non_diags] + col = _aka_np.concatenate((col, col_sup)) + row = _aka_np.concatenate((row, row_sup)) + data = _aka_np.concatenate((data, data_sup)) + + _aka_sparse.coo_matrix.__init__( + self, (data, (row, col)), shape=(sz, sz)) + + FromStress = py11_akantu.FromHigherDim FromTraction = py11_akantu.FromSameDim py11_akantu.__initialize() diff --git a/python/py_aka_common.cc b/python/py_aka_common.cc index db6e49ff4..9f75a7951 100644 --- a/python/py_aka_common.cc +++ b/python/py_aka_common.cc @@ -1,170 +1,175 @@ /* -------------------------------------------------------------------------- */ #include "py_aka_common.hh" #include "py_aka_boundary_conditions.hh" #include "py_aka_error.hh" #include "py_aka_fe_engine.hh" #include "py_aka_material.hh" #include "py_aka_mesh.hh" #include "py_aka_model.hh" #include "py_aka_parser.hh" #include "py_aka_solid_mechanics_model.hh" #include "py_aka_solid_mechanics_model_cohesive.hh" /* -------------------------------------------------------------------------- */ #include /* -------------------------------------------------------------------------- */ #include #include #include /* -------------------------------------------------------------------------- */ namespace py = pybind11; namespace akantu { /* -------------------------------------------------------------------------- */ void register_initialize(py::module & mod) { mod.def("__initialize", []() { int nb_args = 0; char ** null = nullptr; initialize(nb_args, null); }); } /* -------------------------------------------------------------------------- */ void register_enums(py::module & mod) { py::enum_(mod, "SpatialDirection") .value("_x", _x) .value("_y", _y) .value("_z", _z) .export_values(); py::enum_(mod, "AnalysisMethod") .value("_static", _static) .value("_implicit_dynamic", _implicit_dynamic) .value("_explicit_lumped_mass", _explicit_lumped_mass) .value("_explicit_lumped_capacity", _explicit_lumped_capacity) .value("_explicit_consistent_mass", _explicit_consistent_mass) .export_values(); py::enum_(mod, "NonLinearSolverType") .value("_nls_linear", _nls_linear) .value("_nls_newton_raphson", _nls_newton_raphson) .value("_nls_newton_raphson_modified", _nls_newton_raphson_modified) .value("_nls_lumped", _nls_lumped) .value("_nls_auto", _nls_auto) .export_values(); py::enum_(mod, "TimeStepSolverType") .value("_tsst_static", _tsst_static) .value("_tsst_dynamic", _tsst_dynamic) .value("_tsst_dynamic_lumped", _tsst_dynamic_lumped) .value("_tsst_not_defined", _tsst_not_defined) .export_values(); py::enum_(mod, "IntegrationSchemeType") .value("_ist_pseudo_time", _ist_pseudo_time) .value("_ist_forward_euler", _ist_forward_euler) .value("_ist_trapezoidal_rule_1", _ist_trapezoidal_rule_1) .value("_ist_backward_euler", _ist_backward_euler) .value("_ist_central_difference", _ist_central_difference) .value("_ist_fox_goodwin", _ist_fox_goodwin) .value("_ist_trapezoidal_rule_2", _ist_trapezoidal_rule_2) .value("_ist_linear_acceleration", _ist_linear_acceleration) .value("_ist_newmark_beta", _ist_newmark_beta) .value("_ist_generalized_trapezoidal", _ist_generalized_trapezoidal) .export_values(); py::enum_(mod, "SolveConvergenceCriteria") .value("_scc_residual", _scc_residual) .value("_scc_solution", _scc_solution) .value("_scc_residual_mass_wgh", _scc_residual_mass_wgh) .export_values(); py::enum_(mod, "CohesiveMethod") .value("_intrinsic", _intrinsic) .value("_extrinsic", _extrinsic) .export_values(); py::enum_(mod, "GhostType") .value("_not_ghost", _not_ghost) .value("_ghost", _ghost) .value("_casper", _casper) .export_values(); py::enum_(mod, "MeshIOType") .value("_miot_auto", _miot_auto) .value("_miot_gmsh", _miot_gmsh) .value("_miot_gmsh_struct", _miot_gmsh_struct) .value("_miot_diana", _miot_diana) .value("_miot_abaqus", _miot_abaqus) .export_values(); py::enum_(mod, "ModelType") .value("_model", ModelType::_model) .value("_solid_mechanics_model", ModelType::_solid_mechanics_model) .value("_solid_mechanics_model_cohesive", ModelType::_solid_mechanics_model_cohesive) .value("_heat_transfer_model", ModelType::_heat_transfer_model) .value("_structural_mechanics_model", ModelType::_structural_mechanics_model) .value("_embedded_model", ModelType::_embedded_model) .export_values(); py::enum_(mod, "ElementType") .value("_point_1", _point_1) .value("_segment_2", _segment_2) .value("_segment_3", _segment_3) .value("_triangle_3", _triangle_3) .value("_triangle_6", _triangle_6) .value("_quadrangle_4", _quadrangle_4) .value("_quadrangle_8", _quadrangle_8) .value("_tetrahedron_4", _tetrahedron_4) .value("_tetrahedron_10", _tetrahedron_10) .value("_pentahedron_6", _pentahedron_6) .value("_pentahedron_15", _pentahedron_15) .value("_hexahedron_8", _hexahedron_8) .value("_hexahedron_20", _hexahedron_20) .export_values(); py::enum_(mod, "ElementKind") .value("_ek_regular", _ek_regular) .export_values(); + + py::enum_(mod, "MatrixType") + .value("_unsymmetric", _unsymmetric) + .value("_symmetric", _symmetric) + .export_values(); } /* -------------------------------------------------------------------------- */ #define AKANTU_PP_STR_TO_TYPE2(s, data, elem) ({BOOST_PP_STRINGIZE(elem), elem}) void register_functions(py::module & mod) { mod.def("getElementTypes", []() { std::map element_types{ BOOST_PP_SEQ_FOR_EACH_I( AKANTU_PP_ENUM, BOOST_PP_SEQ_SIZE(AKANTU_ek_regular_ELEMENT_TYPE), BOOST_PP_SEQ_TRANSFORM(AKANTU_PP_STR_TO_TYPE2, akantu, AKANTU_ek_regular_ELEMENT_TYPE))}; return element_types; }); } #undef AKANTU_PP_STR_TO_TYPE2 /* -------------------------------------------------------------------------- */ __attribute__((visibility("default"))) void register_all(py::module & mod) { register_initialize(mod); register_enums(mod); register_error(mod); register_functions(mod); register_parser(mod); register_boundary_conditions(mod); register_fe_engine(mod); register_model(mod); register_solid_mechanics_model(mod); register_solid_mechanics_model_cohesive(mod); register_material(mod); register_mesh(mod); } } // namespace akantu diff --git a/python/py_aka_model.cc b/python/py_aka_model.cc index c70c290d4..0379db5af 100644 --- a/python/py_aka_model.cc +++ b/python/py_aka_model.cc @@ -1,63 +1,73 @@ /* -------------------------------------------------------------------------- */ #include "py_aka_array.hh" /* -------------------------------------------------------------------------- */ #include #include -#include +#include /* -------------------------------------------------------------------------- */ #include #include #include /* -------------------------------------------------------------------------- */ namespace py = pybind11; /* -------------------------------------------------------------------------- */ namespace akantu { /* -------------------------------------------------------------------------- */ [[gnu::visibility("default")]] void register_model(py::module & mod) { + py::class_(mod, "SparseMatrix") + .def("getMatrixType", &SparseMatrix::getMatrixType) + .def("size", &SparseMatrix::size); + + py::class_(mod, "SparseMatrixAIJ") + .def("getIRN", &SparseMatrixAIJ::getIRN) + .def("getJCN", &SparseMatrixAIJ::getJCN) + .def("getA", &SparseMatrixAIJ::getA); + py::class_(mod, "DOFManager") .def("getMatrix", - [](DOFManager::getMatrix & self, const std::string & name) { - return + [](DOFManager & self, const std::string & name) { + return dynamic_cast( + self.getMatrix(name)); }, py::return_value_policy::reference); py::class_(mod, "Model") .def("setBaseName", &Model::setBaseName) .def("getFEEngine", &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::dump) .def("initNewSolver", &Model::initNewSolver) .def("getDOFManager", &Model::getDOFManager, py::return_value_policy::reference); py::class_(mod, "NonLinearSolver") .def( "set", [](NonLinearSolver & self, const std::string & id, const Real & val) { if (id == "max_iterations") self.set(id, int(val)); else if (id == "convergence_type") self.set(id, akantu::SolveConvergenceCriteria(UInt(val))); else self.set(id, val); }) .def("set", &NonLinearSolver::set); py::class_(mod, "ModelSolver") .def("getNonLinearSolver", (NonLinearSolver & (ModelSolver::*)(const ID &)) & ModelSolver::getNonLinearSolver, py::arg("solver_id") = "", py::return_value_policy::reference) .def("solveStep", &ModelSolver::solveStep, py::arg("solver_id") = ""); } } // namespace akantu