Commit 9b4cba65 authored by webmanue's avatar webmanue
Browse files

Merge remote-tracking branch 'origin/master' into 89-add-eigenpair-extraction

parents 86ef6203 10d3ad77
Pipeline #138646 passed with stages
in 8 minutes and 21 seconds
......@@ -13,5 +13,10 @@
],
"dockerComposeFile": "../docker-compose.yml",
"service": "dev-real",
"workspaceFolder": "/mnt/io"
}
"workspaceFolder": "/mnt/io",
"remoteEnv": {
"OMPI_ALLOW_RUN_AS_ROOT": "1",
"OMPI_ALLOW_RUN_AS_ROOT_CONFIRM": "1",
"OMPI_MCA_btl_vader_single_copy_mechanism": "none"
}
}
\ No newline at end of file
......@@ -65,6 +65,7 @@ build-image-docu:
check-code-format:
stage: analyze
needs: ["build-image-dev-real"]
image: $CI_REGISTRY_IMAGE/dev-real:$CI_COMMIT_REF_SLUG
script:
- >
......@@ -77,10 +78,11 @@ check-code-format:
check-license-header:
stage: analyze
needs: ["build-image-dev-real"]
image: $CI_REGISTRY_IMAGE/dev-real:$CI_COMMIT_REF_SLUG
script:
- >
! grep \
grep \
--recursive \
--exclude-dir='.git' \
--exclude='*.json' \
......@@ -91,6 +93,7 @@ check-license-header:
check-script:
stage: analyze
needs: ["build-image-dev-real"]
image: $CI_REGISTRY_IMAGE/dev-real:$CI_COMMIT_REF_SLUG
script:
- find . -name "*.py" -print0 | xargs -0 python3 -m black --check
......@@ -100,6 +103,7 @@ check-script:
.build-library: &build-library
stage: build
coverage: '/^lines:\s*\d+\.\d+\%/'
script:
- useradd developer
- mkdir build
......@@ -115,16 +119,20 @@ check-script:
'
artifacts:
reports:
cobertura: build/coverage.xml
coverage_report:
coverage_format: cobertura
path: build/coverage.xml
junit: build/*/test/gtest-results.xml
when: always
build-library-real:
image: $CI_REGISTRY_IMAGE/dev-real:$CI_COMMIT_REF_SLUG
needs: ["build-image-dev-real"]
<<: *build-library
build-library-complex:
image: $CI_REGISTRY_IMAGE/dev-complex:$CI_COMMIT_REF_SLUG
needs: ["build-image-dev-complex"]
<<: *build-library
.install-library: &install-library
......@@ -147,14 +155,17 @@ build-library-complex:
install-library-real:
image: $CI_REGISTRY_IMAGE/dev-real:$CI_COMMIT_REF_SLUG
needs: ["build-image-dev-real"]
<<: *install-library
install-library-complex:
image: $CI_REGISTRY_IMAGE/dev-complex:$CI_COMMIT_REF_SLUG
needs: ["build-image-dev-complex"]
<<: *install-library
build-documentation:
stage: build
needs: ["build-image-docu"]
image: $CI_REGISTRY_IMAGE/docu:$CI_COMMIT_REF_SLUG
script:
- mkdir build
......
......@@ -56,9 +56,9 @@ The project uses [CMake](https://cmake.org) as its build system generator. The f
- [Google Test](https://github.com/google/googletest): version 1.8.1
- [HDF5](https://www.hdfgroup.org/solutions/hdf5/): version 1.10
- [MPI](https://cmake.org/cmake/help/latest/module/FindMPI.html): version 3.1
- [PETSc](https://www.mcs.anl.gov/petsc/): version 3.12
- [Range-v3](https://github.com/ericniebler/range-v3): version 0.10
- [SLEPc](https://slepc.upv.es/): version 3.12
- [PETSc](https://www.mcs.anl.gov/petsc/): version 3.15
- [Range-v3](https://github.com/ericniebler/range-v3): version 0.11
- [SLEPc](https://slepc.upv.es/): version 3.15
Of course, these libraries are covered by their own license terms. Since PETSc and SLEPc do not provide a CMake configuration file, these libraries are found using the provided find modules in ```cmake/modules/```, which in turn are based on ```pkg-config```.
......
......@@ -66,7 +66,7 @@ public:
/**
* @param args Are used to construct the plugins.
*/
template <class... Args> explicit Assembler(Args &&... args);
template <class... Args> explicit Assembler(Args &&...args);
using ElementView = typename cpppetsc::LocalElementView<mesh_type>;
......@@ -82,7 +82,7 @@ public:
* constructor.
*/
template <class... Args>
void emplaceElement(ElementView view, Args &&... constructorArguments);
void emplaceElement(ElementView view, Args &&...constructorArguments);
/**
* @brief Returns a range of iterators pointing to a struct with two methods:
......@@ -191,13 +191,13 @@ namespace assembly {
template <class Element, class Plugins, class Policy>
template <class... Args>
Assembler<Element, Plugins, Policy>::Assembler(Args &&... args)
Assembler<Element, Plugins, Policy>::Assembler(Args &&...args)
: PluginBase(std::forward<Args>(args)...) {}
template <class Element, class Plugins, class Policy>
template <class... Args>
void Assembler<Element, Plugins, Policy>::emplaceElement(
ElementView view, Args &&... constructorArguments) {
ElementView view, Args &&...constructorArguments) {
_elements.emplace_back(std::move(view),
std::forward<Args>(constructorArguments)...);
}
......@@ -220,7 +220,7 @@ public:
* @param args Are used to construct in-place the element instance.
*/
template <class... Args>
explicit AnnotatedElement(ElementView view, Args &&... args)
explicit AnnotatedElement(ElementView view, Args &&...args)
: _meshView(std::move(view)), _instance(std::forward<Args>(args)...) {}
const ElementView &meshView() const { return _meshView; }
......
......@@ -38,7 +38,7 @@ class AssemblerGroup : public DerivePluginsUniquely<
typename ConcatenatePlugins<typename PluginTypeTrait<
SingleElementAssemblers>::type...>::type> {
public:
template <class... Args> explicit AssemblerGroup(Args &&... args);
template <class... Args> explicit AssemblerGroup(Args &&...args);
/**
* @brief Get the Nth member assembler.
......@@ -158,7 +158,7 @@ namespace assembly {
template <class... SingleElementAssemblers>
template <class... Args>
AssemblerGroup<SingleElementAssemblers...>::AssemblerGroup(Args &&... args)
AssemblerGroup<SingleElementAssemblers...>::AssemblerGroup(Args &&...args)
: _assemblers(std::forward<Args>(args)...) {}
} // namespace assembly
......
......@@ -21,7 +21,7 @@
#include <utility>
#define DEFINE_ASSEMBLER_METHOD_BASE(PluginName, methodName, cvQualifiers) \
template <class... Args> void methodName(Args &&... args) cvQualifiers { \
template <class... Args> void methodName(Args &&...args) cvQualifiers { \
if constexpr (!::ae108::assembly::IsGroupTypeTrait< \
typename PluginName::assembler_type>::value) { \
execute(std::forward<Args>(args)...); \
......@@ -33,7 +33,7 @@
} \
}; \
std::apply( \
[&](auto &... as) { (..., std::invoke(conditionalCall, as)); }, \
[&](auto &...as) { (..., std::invoke(conditionalCall, as)); }, \
this->assembler().assemblers()); \
} \
}
......
......@@ -29,7 +29,7 @@ DEFINE_CONST_ASSEMBLER_PLUGIN(AssembleConsistentMassMatrixPlugin,
std::vector<value_type> outputBuffer;
for (const auto &meshElement : this->assembler().meshElements()) {
const auto matrix = meshElement.instance().computeConsistentMassMatrix();
const auto matrix = compute_consistent_mass_matrix(meshElement.instance());
outputBuffer.resize(matrix.rows() * matrix.cols());
utilities::serialize(matrix, outputBuffer.begin());
......
......@@ -29,7 +29,7 @@ DEFINE_CONST_ASSEMBLER_PLUGIN(AssembleLumpedMassMatrixPlugin,
std::vector<value_type> outputBuffer;
for (const auto &meshElement : this->assembler().meshElements()) {
const auto matrix = meshElement.instance().computeLumpedMassMatrix();
const auto matrix = compute_lumped_mass_matrix(meshElement.instance());
outputBuffer.resize(matrix.rows() * matrix.cols());
utilities::serialize(matrix, outputBuffer.begin());
......
......@@ -50,6 +50,17 @@ template <class ValueType_> struct Element {
MOCK_METHOD2_T(updateInternalVariables,
void(const NodalDisplacements &, const double));
};
template <class ValueType_>
auto compute_consistent_mass_matrix(const Element<ValueType_> &element) {
return element.computeConsistentMassMatrix();
}
template <class ValueType_>
auto compute_lumped_mass_matrix(const Element<ValueType_> &element) {
return element.computeLumpedMassMatrix();
}
} // namespace test
} // namespace assembly
} // namespace ae108
......@@ -20,9 +20,9 @@ find_dependency(Boost 1.67 COMPONENTS program_options REQUIRED)
find_dependency(Eigen3 3.3 CONFIG REQUIRED)
list(APPEND CMAKE_MODULE_PATH "${CMAKE_CURRENT_LIST_DIR}/modules")
find_dependency(AE108_PETSc MODULE 3.12)
find_dependency(AE108_PETSc MODULE 3.15)
find_dependency(range-v3 0.10.0 CONFIG REQUIRED)
find_dependency(range-v3 0.11.0 CONFIG REQUIRED)
foreach(AE108_LIBRARY elements cpppetsc cppslepc assembly solve cmdline)
include("${CMAKE_CURRENT_LIST_DIR}/ae108-${AE108_LIBRARY}-export.cmake")
......
# © 2021 ETH Zurich, Mechanics and Materials Lab
#
# Licensed under the Apache License, Version 2.0 (the "License");
# you may not use this file except in compliance with the License.
# You may obtain a copy of the License at
#
# http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
find_package(GTest 1.8.1 CONFIG)
if (NOT GTest_FOUND)
set(FETCHCONTENT_SOURCE_DIR_AE108_GOOGLETEST "/usr/src/googletest")
if (IS_DIRECTORY "${FETCHCONTENT_SOURCE_DIR_AE108_GOOGLETEST}")
message(STATUS
"GTest was not found, using source at ${FETCHCONTENT_SOURCE_DIR_AE108_GOOGLETEST} instead."
)
include(FetchContent)
FetchContent_Declare("ae108_googletest")
FetchContent_GetProperties(
ae108_googletest
SOURCE_DIR ae108_googletest_SOURCE_DIR
BINARY_DIR ae108_googletest_BINARY_DIR
POPULATED ae108_googletest_POPULATED
)
if (NOT ae108_googletest_POPULATED)
FetchContent_Populate("ae108_googletest")
add_subdirectory(
"${ae108_googletest_SOURCE_DIR}"
"${ae108_googletest_BINARY_DIR}"
EXCLUDE_FROM_ALL
)
endif()
foreach(LIBRARY_NAME gtest gmock gtest_main gmock_main)
add_library("GTest::${LIBRARY_NAME}" ALIAS "${LIBRARY_NAME}")
endforeach()
set(GTest_FOUND TRUE)
endif()
endif()
find_package_handle_standard_args(AE108_GTest
REQUIRED_VARS GTest_FOUND
)
\ No newline at end of file
......@@ -51,6 +51,7 @@ if(AE108_PETSc_FOUND AND NOT TARGET ae108::external::petsc)
endif()
try_compile_with_petsc("-DAE108_PETSC_COMPLEX")
set(AE108_PETSC_COMPLEX ${AE108_PETSc_COMPILE_RESULT} CACHE BOOL "ae108: PETSc with complex scalar type")
if(AE108_PETSc_COMPILE_RESULT)
target_compile_definitions(ae108::external::petsc
INTERFACE AE108_PETSC_COMPLEX=1
......
......@@ -14,8 +14,8 @@
# limitations under the License.
include(FindPkgConfig)
pkg_check_modules(SLEPc SLEPc IMPORTED_TARGET)
find_package(AE108_PETSc 1.12 MODULE)
pkg_check_modules(SLEPc slepc IMPORTED_TARGET)
find_package(AE108_PETSc 3.15 MODULE)
find_package_handle_standard_args(AE108_SLEPc
REQUIRED_VARS SLEPc_FOUND AE108_PETSc_FOUND
......@@ -28,4 +28,4 @@ if(AE108_SLEPc_FOUND AND NOT TARGET ae108::external::slepc)
INTERFACE PkgConfig::SLEPc
INTERFACE ae108::external::petsc
)
endif()
\ No newline at end of file
endif()
......@@ -19,7 +19,7 @@ add_executable(${PROJECT_NAME}-CmdLineTests
./CommandLineOptionParser_Test.cc
)
find_package(AE108_GTest MODULE REQUIRED)
find_package(GTest 1.8.1 CONFIG REQUIRED)
target_link_libraries(${PROJECT_NAME}-CmdLineTests
PRIVATE ${PROJECT_NAME}::cmdline
PRIVATE GTest::gmock_main
......
......@@ -76,8 +76,8 @@ target_include_directories(${PROJECT_NAME}-cpppetsc PUBLIC
$<BUILD_INTERFACE:${CMAKE_CURRENT_SOURCE_DIR}/include>
)
find_package(AE108_PETSc 3.12 MODULE REQUIRED)
find_package(range-v3 0.10.0 CONFIG REQUIRED)
find_package(AE108_PETSc 3.15 MODULE REQUIRED)
find_package(range-v3 0.11.0 CONFIG REQUIRED)
target_link_libraries(${PROJECT_NAME}-cpppetsc
PUBLIC ${PROJECT_NAME}::external::petsc
PUBLIC range-v3::range-v3
......
......@@ -37,10 +37,6 @@ template <> PetscErrorCode callDestructor(Vec *const ptr) noexcept {
return VecDestroy(ptr);
}
template <> PetscErrorCode callDestructor(VecScatter *const ptr) noexcept {
return VecScatterDestroy(ptr);
}
template <> PetscErrorCode callDestructor(PetscSF *const ptr) noexcept {
return PetscSFDestroy(ptr);
}
......
......@@ -189,9 +189,15 @@ private:
static UniqueEntity<DM> createDM();
/**
* @brief Adds the default uniform section to the mesh.
* @brief Creates an uniform section with the provided degrees of freedom.
*/
void addSection(const size_type dofPerVertex, const size_type dofPerElement);
UniqueEntity<PetscSection> createSection(const size_type dofPerVertex,
const size_type dofPerElement) const;
/**
* @brief Sets the default section of the mesh.
*/
void setSection(UniqueEntity<PetscSection> section);
/**
* @brief Adds the chart as described by the connectivity to the mesh.
......@@ -237,6 +243,9 @@ private:
void assertCorrectBaseMesh(const vector_type &vector) const;
UniqueEntity<DM> _mesh;
UniqueEntity<PetscSF> _migration;
UniqueEntity<PetscSection> _section;
UniqueEntity<PetscSF> _globalToNatural;
size_type _totalNumberOfElements = 0;
size_type _totalNumberOfVertices = 0;
};
......@@ -304,9 +313,12 @@ Mesh<Policy> Mesh<Policy>::fromConnectivity(const size_type dimension,
// calculate the strata
Policy::handleError(DMPlexStratify(mesh._mesh.get()));
// set unknown cell type
Policy::handleError(DMCreateLabel(mesh._mesh.get(), "celltype"));
distributeMesh(&mesh);
mesh.addSection(dofPerVertex, dofPerElement);
mesh.setSection(mesh.createSection(dofPerVertex, dofPerElement));
mesh.setGlobalToNaturalSF(mesh.createReorderingSF(
mesh.canonicalRowIndices(dofPerVertex, dofPerElement)));
......@@ -383,8 +395,9 @@ UniqueEntity<PetscSF> Mesh<Policy>::createReorderingSF(
template <class Policy>
void Mesh<Policy>::setGlobalToNaturalSF(UniqueEntity<PetscSF> globalToNatural) {
_globalToNatural = std::move(globalToNatural);
Policy::handleError(
DMPlexSetGlobalToNaturalSF(_mesh.get(), globalToNatural.release()));
DMPlexSetGlobalToNaturalSF(_mesh.get(), _globalToNatural.get()));
}
template <class Policy>
......@@ -402,13 +415,17 @@ Mesh<Policy> Mesh<Policy>::cloneWithDofs(const size_type dofPerVertex,
auto migration = PetscSF();
Policy::handleError(DMPlexGetMigrationSF(_mesh.get(), &migration));
if (migration) {
auto cloned = PetscSF();
mesh._migration = [&]() {
auto sf = PetscSF();
Policy::handleError(
PetscSFDuplicate(migration, PETSCSF_DUPLICATE_GRAPH, &sf));
return makeUniqueEntity<Policy>(sf);
}();
Policy::handleError(
PetscSFDuplicate(migration, PETSCSF_DUPLICATE_GRAPH, &cloned));
Policy::handleError(DMPlexSetMigrationSF(mesh._mesh.get(), cloned));
DMPlexSetMigrationSF(mesh._mesh.get(), mesh._migration.get()));
}
mesh.addSection(dofPerVertex, dofPerElement);
mesh.setSection(mesh.createSection(dofPerVertex, dofPerElement));
mesh.setGlobalToNaturalSF(mesh.createReorderingSF(
mesh.canonicalRowIndices(dofPerVertex, dofPerElement)));
......@@ -476,15 +493,17 @@ void Mesh<Policy>::addChart(const Container &elementVertexIDs,
}
template <class Policy>
void Mesh<Policy>::addSection(const size_type dofPerVertex,
const size_type dofPerElement) {
auto section = PetscSection();
UniqueEntity<PetscSection>
Mesh<Policy>::createSection(const size_type dofPerVertex,
const size_type dofPerElement) const {
const std::array<size_type, 1> numberOfComponents = {{1}};
std::vector<size_type> numberOfDofsPerDim(2);
numberOfDofsPerDim.front() = dofPerVertex;
numberOfDofsPerDim.back() = dofPerElement;
Policy::handleError(DMSetNumFields(_mesh.get(), numberOfComponents.size()));
auto section = PetscSection();
Policy::handleError(DMPlexCreateSection(
_mesh.get(), nullptr /* label */, numberOfComponents.data(),
numberOfDofsPerDim.data(), 0 /* number of boundary conditions */,
......@@ -492,11 +511,13 @@ void Mesh<Policy>::addSection(const size_type dofPerVertex,
nullptr /* boundary condition components */,
nullptr /* boundary condition points */, nullptr /* permutation */,
&section));
return makeUniqueEntity<Policy>(section);
}
Policy::handleError(DMSetSection(_mesh.get(), section));
// create default global section
Policy::handleError(DMGetGlobalSection(_mesh.get(), &section));
template <class Policy>
void Mesh<Policy>::setSection(UniqueEntity<PetscSection> section) {
_section = std::move(section);
Policy::handleError(DMSetSection(_mesh.get(), _section.get()));
}
template <class Policy> void Mesh<Policy>::distributeMesh(Mesh *const mesh) {
......
// © 2020 ETH Zurich, Mechanics and Materials Lab
// © 2020, 2022 ETH Zurich, Mechanics and Materials Lab
// © 2020 California Institute of Technology
//
// Licensed under the Apache License, Version 2.0 (the "License");
......@@ -25,6 +25,7 @@
#include "ae108/cpppetsc/UniqueEntity.h"
#include "ae108/cpppetsc/Vector.h"
#include <functional>
#include <optional>
#include <petscmat.h>
#include <petscmath.h>
#include <petscsys.h>
......@@ -71,12 +72,17 @@ public:
gpcg,
nm,
pounders,
brgn,
lcl,
ssils,
ssfls,
asils,
asfls,
ipm
ipm,
pdipm,
shell,
admm,
almm,
};
/**
......@@ -132,6 +138,29 @@ public:
void setBounds(const distributed<vector_type> &lowerBounds,
const distributed<vector_type> &upperBounds);
struct EqualityConstraints {
/**
* @brief The solver attempts to reach a residual of zero.
*/
GradientFunctor residual;
/**
* @brief The jacobian of the residual.
*/
HessianFunctor jacobian;
/**
* @brief Used to store residuals and jacobians internally.
*/
std::pair<distributed<vector_type>, matrix_type> buffer;
};
/**
* @brief Set the equality constraints of the problem.
* @note Not all solvers consider these constraints.
*/
void setConstraints(EqualityConstraints constraints);
/**
* @brief Minimizes the objective function using optimization.
*
......@@ -164,6 +193,8 @@ private:
UniqueEntity<Tao> _tao;
matrix_type _buffer_matrix;
std::optional<EqualityConstraints> _constraints;
};
template <class Policy>
......@@ -288,6 +319,10 @@ template <class Policy> void TAOSolver<Policy>::setType(const Type type) {
taoType = TAOPOUNDERS;
break;
}
case Type::brgn: {
taoType = TAOBRGN;
break;
}
case Type::lcl: {
taoType = TAOLCL;
break;
......@@ -312,6 +347,22 @@ template <class Policy> void TAOSolver<Policy>::setType(const Type type) {
taoType = TAOIPM;
break;
}
case Type::pdipm: {
taoType = TAOPDIPM;
break;
}
case Type::shell: {
taoType = TAOSHELL;
break;
}
case Type::admm: {
taoType = TAOADMM;
break;
}
case Type::almm: {
taoType = TAOALMM;
break;
}
}
Policy::handleError(TaoSetType(_tao.get(), taoType));
}
......@@ -323,6 +374,19 @@ void TAOSolver<Policy>::setBounds(const distributed<vector_type> &lowerBounds,
_tao.get(), lowerBounds.unwrap().data(), upperBounds.unwrap().data()));
}
template <class Policy>
void TAOSolver<Policy>::setConstraints(EqualityConstraints constraints) {
_constraints.emplace(std::move(constraints));
Policy::handleError(TaoSetEqualityConstraintsRoutine(
_tao.get(), _constraints->buffer.first.unwrap().data(),
&TAOSolver::gradientAdapter, &_constraints->residual));
Policy::handleError(TaoSetJacobianEqualityRoutine(
_tao.get(), _constraints->buffer.second.data(),
_constraints->buffer.second.data(), &TAOSolver::hessianAdapter,
&_constraints->jacobian));
}
template <class Policy> void TAOSolver<Policy>::checkConvergence() const {
const auto errorCode = [this]() {
auto reason = TaoConvergedReason{};
......
......@@ -33,21 +33,22 @@ public:
/**
* @brief Constructs the contained value with the given arguments.
*/
template <class... Args> explicit constexpr TaggedEntity(Args &&... args);
template <class... Args> explicit constexpr TaggedEntity(Args &&...args);
/**
* @brief Calls operator() of the contained value.
*/
template <class... Args>
auto operator()(Args &&... args)
auto operator()(Args &&...args)
-> decltype(std::declval<value_type>()(std::forward<Args>(args)...));
/**
* @brief Calls operator() of the contained value.
*/
template <class... Args>
constexpr auto operator()(Args &&... args) const -> decltype(
std::declval<const value_type>()(std::forward<Args>(args)...));
constexpr auto operator()(Args &&...args) const
-> decltype(std::declval<const value_type>()(