diff --git a/src/materials/CMakeLists.txt b/src/materials/CMakeLists.txt index 3cf4bca..b5a089e 100644 --- a/src/materials/CMakeLists.txt +++ b/src/materials/CMakeLists.txt @@ -1,56 +1,56 @@ # ============================================================================= # file CMakeLists.txt # # @author Till Junge # # @date 08 Jan 2018 # # @brief configuration for material laws # # @section LICENSE # # Copyright © 2018 Till Junge # # µSpectre 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, or (at # your option) any later version. # # µSpectre 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 # General Public License for more details. # # You should have received a copy of the GNU Lesser General Public License # along with µSpectre; see the file COPYING. If not, write to the # Free Software Foundation, Inc., 59 Temple Place - Suite 330, # Boston, MA 02111-1307, USA. # # Additional permission under GNU GPL version 3 section 7 # # If you modify this Program, or any covered work, by linking or combining it # with proprietary FFT implementations or numerical libraries, containing parts # covered by the terms of those libraries' licenses, the licensors of this # Program grant you additional permission to convey the resulting work. # ============================================================================= set (materials_SRC ${CMAKE_CURRENT_SOURCE_DIR}/material_base.cc ${CMAKE_CURRENT_SOURCE_DIR}/material_anisotropic.cc ${CMAKE_CURRENT_SOURCE_DIR}/material_orthotropic.cc ${CMAKE_CURRENT_SOURCE_DIR}/material_linear_elastic1.cc ${CMAKE_CURRENT_SOURCE_DIR}/material_linear_elastic2.cc ${CMAKE_CURRENT_SOURCE_DIR}/material_linear_elastic3.cc ${CMAKE_CURRENT_SOURCE_DIR}/material_linear_elastic4.cc ${CMAKE_CURRENT_SOURCE_DIR}/material_hyper_elasto_plastic1.cc - + ${CMAKE_CURRENT_SOURCE_DIR}/material_linear_elastic_generic.cc ) if(${SPLIT_CELL}) set (materials_SRC ${materials_SRC} ${CMAKE_CURRENT_SOURCE_DIR}/laminate_homogenisation.cc ) endif() target_sources(muSpectre PRIVATE ${materials_SRC}) diff --git a/src/materials/material_linear_elastic_generic.cc b/src/materials/material_linear_elastic_generic.cc new file mode 100644 index 0000000..54390f4 --- /dev/null +++ b/src/materials/material_linear_elastic_generic.cc @@ -0,0 +1,71 @@ +/** + * @file material_linear_elastic_generic.cc + * + * @author Till Junge + * + * @date 21 Sep 2018 + * + * @brief implementation for MaterialLinearElasticGeneric + * + * Copyright © 2018 Till Junge + * + * µSpectre 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, or (at + * your option) any later version. + * + * µSpectre 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 + * General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public License + * along with µSpectre; see the file COPYING. If not, write to the + * Free Software Foundation, Inc., 59 Temple Place - Suite 330, + * * Boston, MA 02111-1307, USA. + * + * Additional permission under GNU GPL version 3 section 7 + * + * If you modify this Program, or any covered work, by linking or combining it + * with proprietary FFT implementations or numerical libraries, containing parts + * covered by the terms of those libraries' licenses, the licensors of this + * Program grant you additional permission to convey the resulting work. + */ + +#include "materials/material_linear_elastic_generic.hh" +#include "common/voigt_conversion.hh" + +namespace muSpectre { + + /* ---------------------------------------------------------------------- */ + template + MaterialLinearElasticGeneric::MaterialLinearElasticGeneric( + const std::string & name, const CInput_t & C_voigt) + : Parent{name} { + using VC_t = VoigtConversion; + constexpr Dim_t VSize{vsize(DimM)}; + if (not(C_voigt.rows() == VSize) or not(C_voigt.cols() == VSize)) { + std::stringstream err_str{}; + err_str << "The stiffness tensor should be input as a " << VSize << " × " + << VSize << " Matrix in Voigt notation. You supplied" + << " a " << C_voigt.rows() << " × " << C_voigt.cols() + << " matrix"; + } + + for (int i{0}; i < DimM; ++i) { + for (int j{0}; j < DimM; ++j) { + for (int k{0}; k < DimM; ++k) { + for (int l{0}; l < DimM; ++l) { + get(this->C, i, j, k, l) = + C_voigt(VC_t::sym_mat(i, j), VC_t::sym_mat(k, l)); + } + } + } + } + } + + template class MaterialLinearElasticGeneric; + template class MaterialLinearElasticGeneric; + template class MaterialLinearElasticGeneric; + +} // namespace muSpectre diff --git a/src/materials/material_linear_elastic_generic.hh b/src/materials/material_linear_elastic_generic.hh new file mode 100644 index 0000000..31b078a --- /dev/null +++ b/src/materials/material_linear_elastic_generic.hh @@ -0,0 +1,187 @@ +/** + * @file material_linear_elastic_generic.hh + * + * @author Till Junge + * + * @date 21 Sep 2018 + * + * @brief Implementation fo a generic linear elastic material that + * stores the full elastic stiffness tensor. Convenient but not the + * most efficient + * + * Copyright © 2018 Till Junge + * + * µSpectre 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, or (at + * your option) any later version. + * + * µSpectre 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 + * General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public License + * along with µSpectre; see the file COPYING. If not, write to the + * Free Software Foundation, Inc., 59 Temple Place - Suite 330, + * * Boston, MA 02111-1307, USA. + * + * Additional permission under GNU GPL version 3 section 7 + * + * If you modify this Program, or any covered work, by linking or combining it + * with proprietary FFT implementations or numerical libraries, containing parts + * covered by the terms of those libraries' licenses, the licensors of this + * Program grant you additional permission to convey the resulting work. + */ + +#ifndef SRC_MATERIALS_MATERIAL_LINEAR_ELASTIC_GENERIC_HH_ +#define SRC_MATERIALS_MATERIAL_LINEAR_ELASTIC_GENERIC_HH_ + +#include "common/common.hh" +#include "common/T4_map_proxy.hh" +#include "materials/material_muSpectre_base.hh" +#include "common/tensor_algebra.hh" + +namespace muSpectre { + + /** + * forward declaration + */ + template + class MaterialLinearElasticGeneric; + + /** + * traits for use by MaterialMuSpectre for crtp + */ + template + struct MaterialMuSpectre_traits> { + //! global field collection + using GFieldCollection_t = + typename MaterialBase::GFieldCollection_t; + + //! expected map type for strain fields + using StrainMap_t = + MatrixFieldMap; + //! expected map type for stress fields + using StressMap_t = MatrixFieldMap; + //! expected map type for tangent stiffness fields + using TangentMap_t = T4MatrixFieldMap; + + //! declare what type of strain measure your law takes as input + constexpr static auto strain_measure{StrainMeasure::GreenLagrange}; + //! declare what type of stress measure your law yields as output + constexpr static auto stress_measure{StressMeasure::PK2}; + + //! elasticity without internal variables + using InternalVariables = std::tuple<>; + }; + + /** + * Linear elastic law defined by a full stiffness tensor. Very + * generic, but not most efficient + */ + template + class MaterialLinearElasticGeneric + : public MaterialMuSpectre, DimS, + DimM> { + public: + //! parent type + using Parent = + MaterialMuSpectre, DimS, DimM>; + //! generic input tolerant to python input + using CInput_t = + Eigen::Ref, 0, + Eigen::Stride>; + //! Default constructor + MaterialLinearElasticGeneric() = delete; + + /** + * Constructor by name and stiffness tensor. + * + * @param name unique material name + * @param C_voigt elastic tensor in Voigt notation + */ + MaterialLinearElasticGeneric(const std::string & name, + const CInput_t & C_voigt); + + //! Copy constructor + MaterialLinearElasticGeneric(const MaterialLinearElasticGeneric & other) = + delete; + + //! Move constructor + MaterialLinearElasticGeneric(MaterialLinearElasticGeneric && other) = + delete; + + //! Destructor + virtual ~MaterialLinearElasticGeneric() = default; + + //! Copy assignment operator + MaterialLinearElasticGeneric & + operator=(const MaterialLinearElasticGeneric & other) = delete; + + //! Move assignment operator + MaterialLinearElasticGeneric & + operator=(MaterialLinearElasticGeneric && other) = delete; + + //! see + //! http://eigen.tuxfamily.org/dox/group__TopicStructHavingEigenMembers.html + EIGEN_MAKE_ALIGNED_OPERATOR_NEW; + + /** + * evaluates second Piola-Kirchhoff stress given the Green-Lagrange + * strain (or Cauchy stress if called with a small strain tensor) + */ + template + inline decltype(auto) evaluate_stress(const Eigen::MatrixBase & E); + + /** + * evaluates both second Piola-Kirchhoff stress and stiffness given + * the Green-Lagrange strain (or Cauchy stress and stiffness if + * called with a small strain tensor) + */ + template + inline decltype(auto) evaluate_stress_tangent(s_t && E); + + /** + * return the empty internals tuple + */ + std::tuple<> & get_internals() { return this->internal_variables; } + + /** + * return a reference to teh stiffness tensor + */ + const T4Mat & get_C() const { return this->C; } + + protected: + T4Mat C{}; + //! empty tuple + std::tuple<> internal_variables{}; + + private: + }; + + /* ---------------------------------------------------------------------- */ + template + template + auto MaterialLinearElasticGeneric::evaluate_stress( + const Eigen::MatrixBase & E) -> decltype(auto) { + static_assert(Derived::ColsAtCompileTime == DimM, "wrong input size"); + static_assert(Derived::RowsAtCompileTime == DimM, "wrong input size"); + return Matrices::tensmult(this->C, E); + } + + /* ---------------------------------------------------------------------- */ + template + template + auto + MaterialLinearElasticGeneric::evaluate_stress_tangent(s_t && E) + -> decltype(auto) { + using Stress_t = decltype(this->evaluate_stress(std::forward(E))); + using Stiffness_t = Eigen::Map>; + using Ret_t = std::tuple; + return Ret_t{this->evaluate_stress(std::forward(E)), + Stiffness_t(this->C.data())}; + } +} // namespace muSpectre + +#endif // SRC_MATERIALS_MATERIAL_LINEAR_ELASTIC_GENERIC_HH_