From 704f77c036b719191f3a27be23ca07199f9692c1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?St=C3=A9phane=20Del=20Pino?= <stephane.delpino44@gmail.com> Date: Mon, 14 Feb 2022 16:32:23 +0100 Subject: [PATCH] Improve variable compatibility checking Improve the following - check that primal variables are defined on the same mesh - check that dual variables are defined on the same mesh - check that the dual mesh is the diamond dual mesh of the primal mesh --- src/scheme/ScalarDiamondScheme.cpp | 45 +++++++++++++++++------- src/scheme/VectorDiamondScheme.cpp | 56 +++++++++++++++++------------- 2 files changed, 63 insertions(+), 38 deletions(-) diff --git a/src/scheme/ScalarDiamondScheme.cpp b/src/scheme/ScalarDiamondScheme.cpp index 32eafcc75..07561e6f5 100644 --- a/src/scheme/ScalarDiamondScheme.cpp +++ b/src/scheme/ScalarDiamondScheme.cpp @@ -281,13 +281,10 @@ class ScalarDiamondSchemeHandler::ScalarDiamondScheme : public ScalarDiamondSche const std::shared_ptr<const DiscreteFunctionP0<Dimension, double>>& f, const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list) { - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mu->mesh()) { - throw NormalError("diffusion coefficient is not defined on the dual mesh!"); - } - - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mub->mesh()) { - throw NormalError("boundary diffusion coefficient is not defined on the dual mesh!"); - } + Assert(DualMeshManager::instance().getDiamondDualMesh(*mesh) == dual_mu->mesh(), + "diffusion coefficient is not defined on the dual mesh!"); + Assert(DualMeshManager::instance().getDiamondDualMesh(*mesh) == dual_mub->mesh(), + "boundary diffusion coefficient is not defined on the dual mesh!"); using BoundaryCondition = std::variant<DirichletBoundaryCondition, FourierBoundaryCondition, NeumannBoundaryCondition, SymmetryBoundaryCondition>; @@ -781,6 +778,13 @@ ScalarDiamondSchemeHandler::ScalarDiamondSchemeHandler( const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list) { const std::shared_ptr i_mesh = getCommonMesh({alpha, f}); + if (not i_mesh) { + throw NormalError("primal discrete functions are not defined on the same mesh"); + } + const std::shared_ptr i_dual_mesh = getCommonMesh({dual_mub, dual_mu}); + if (not i_dual_mesh) { + throw NormalError("dual discrete functions are not defined on the same mesh"); + } checkDiscretizationType({alpha, dual_mub, dual_mu, f}, DiscreteFunctionType::P0); switch (i_mesh->dimension()) { @@ -788,9 +792,14 @@ ScalarDiamondSchemeHandler::ScalarDiamondSchemeHandler( using MeshType = Mesh<Connectivity<1>>; using DiscreteScalarFunctionType = DiscreteFunctionP0<1, double>; + std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh); + + if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != i_dual_mesh) { + throw NormalError("dual variables are is not defined on the diamond dual of the primal mesh"); + } + m_scheme = - std::make_unique<ScalarDiamondScheme<1>>(std::dynamic_pointer_cast<const MeshType>(i_mesh), - std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(alpha), + std::make_unique<ScalarDiamondScheme<1>>(mesh, std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(alpha), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(dual_mub), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(dual_mu), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(f), @@ -801,9 +810,14 @@ ScalarDiamondSchemeHandler::ScalarDiamondSchemeHandler( using MeshType = Mesh<Connectivity<2>>; using DiscreteScalarFunctionType = DiscreteFunctionP0<2, double>; + std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh); + + if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != i_dual_mesh) { + throw NormalError("dual variables are is not defined on the diamond dual of the primal mesh"); + } + m_scheme = - std::make_unique<ScalarDiamondScheme<2>>(std::dynamic_pointer_cast<const MeshType>(i_mesh), - std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(alpha), + std::make_unique<ScalarDiamondScheme<2>>(mesh, std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(alpha), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(dual_mub), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(dual_mu), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(f), @@ -814,9 +828,14 @@ ScalarDiamondSchemeHandler::ScalarDiamondSchemeHandler( using MeshType = Mesh<Connectivity<3>>; using DiscreteScalarFunctionType = DiscreteFunctionP0<3, double>; + std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh); + + if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != i_dual_mesh) { + throw NormalError("dual variables are is not defined on the diamond dual of the primal mesh"); + } + m_scheme = - std::make_unique<ScalarDiamondScheme<3>>(std::dynamic_pointer_cast<const MeshType>(i_mesh), - std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(alpha), + std::make_unique<ScalarDiamondScheme<3>>(mesh, std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(alpha), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(dual_mub), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(dual_mu), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(f), diff --git a/src/scheme/VectorDiamondScheme.cpp b/src/scheme/VectorDiamondScheme.cpp index 32cfbdab5..368713589 100644 --- a/src/scheme/VectorDiamondScheme.cpp +++ b/src/scheme/VectorDiamondScheme.cpp @@ -289,9 +289,8 @@ class VectorDiamondSchemeHandler::VectorDiamondScheme : public VectorDiamondSche Assert(dual_lambda->mesh() == dual_mu->mesh()); Assert(dual_lambdab->mesh() == dual_mu->mesh()); Assert(dual_mub->mesh() == dual_mu->mesh()); - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mu->mesh()) { - throw NormalError("diffusion coefficient is not defined on the dual mesh!"); - } + Assert(DualMeshManager::instance().getDiamondDualMesh(*mesh) == dual_mu->mesh(), + "diffusion coefficient is not defined on the dual mesh!"); using MeshDataType = MeshData<Dimension>; @@ -1244,9 +1243,8 @@ class EnergyComputerHandler::EnergyComputer : public EnergyComputerHandler::IEne Assert(dual_lambdab->mesh() == dual_U->mesh()); Assert(U->mesh() == source->mesh()); Assert(dual_lambdab->mesh() == dual_mub->mesh()); - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mub->mesh()) { - throw NormalError("diffusion coefficient is not defined on the dual mesh!"); - } + Assert(DualMeshManager::instance().getDiamondDualMesh(*mesh) == dual_mub->mesh(), + "diffusion coefficient is not defined on the dual mesh!"); using MeshDataType = MeshData<Dimension>; @@ -1865,8 +1863,14 @@ VectorDiamondSchemeHandler::VectorDiamondSchemeHandler( const std::shared_ptr<const IDiscreteFunction>& f, const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list) { - const std::shared_ptr i_mesh = getCommonMesh({alpha, f}); - const std::shared_ptr i_dual_mesh = getCommonMesh({dual_lambda, dual_mu}); + const std::shared_ptr i_mesh = getCommonMesh({alpha, f}); + if (not i_mesh) { + throw NormalError("primal discrete functions are not defined on the same mesh"); + } + const std::shared_ptr i_dual_mesh = getCommonMesh({dual_lambda, dual_lambdab, dual_mu, dual_mub}); + if (not i_dual_mesh) { + throw NormalError("dual discrete functions are not defined on the same mesh"); + } checkDiscretizationType({alpha, dual_lambdab, dual_mub, dual_lambda, dual_mu, f}, DiscreteFunctionType::P0); switch (i_mesh->dimension()) { @@ -1877,10 +1881,6 @@ VectorDiamondSchemeHandler::VectorDiamondSchemeHandler( std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh); - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mu->mesh()) { - throw NormalError("mu_dual is not defined on the dual mesh"); - } - m_scheme = std::make_unique<VectorDiamondScheme<1>>(mesh, std::dynamic_pointer_cast<const DiscreteScalarFunctionType>(alpha), std::dynamic_pointer_cast<const DiscreteScalarFunctionType>( @@ -1899,8 +1899,8 @@ VectorDiamondSchemeHandler::VectorDiamondSchemeHandler( std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh); - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mu->mesh()) { - throw NormalError("mu_dual is not defined on the dual mesh"); + if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != i_dual_mesh) { + throw NormalError("dual variables are is not defined on the diamond dual of the primal mesh"); } m_scheme = @@ -1921,8 +1921,8 @@ VectorDiamondSchemeHandler::VectorDiamondSchemeHandler( std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh); - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mu->mesh()) { - throw NormalError("mu_dual is not defined on the dual mesh"); + if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != i_dual_mesh) { + throw NormalError("dual variables are is not defined on the diamond dual of the primal mesh"); } m_scheme = @@ -1952,9 +1952,15 @@ EnergyComputerHandler::EnergyComputerHandler( const std::shared_ptr<const IDiscreteFunction>& source, const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list) { - const std::shared_ptr i_mesh = getCommonMesh({U, source}); - const std::shared_ptr i_dual_mesh = getCommonMesh({dual_lambdab, dual_mub}); - checkDiscretizationType({dual_lambdab, dual_mub, U, source}, DiscreteFunctionType::P0); + const std::shared_ptr i_mesh = getCommonMesh({U, source}); + if (not i_mesh) { + throw NormalError("primal discrete functions are not defined on the same mesh"); + } + const std::shared_ptr i_dual_mesh = getCommonMesh({dual_lambdab, dual_mub, dual_U}); + if (not i_dual_mesh) { + throw NormalError("dual discrete functions are not defined on the same mesh"); + } + checkDiscretizationType({dual_lambdab, dual_mub, dual_U, U, source}, DiscreteFunctionType::P0); switch (i_mesh->dimension()) { case 1: { @@ -1964,8 +1970,8 @@ EnergyComputerHandler::EnergyComputerHandler( std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh); - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mub->mesh()) { - throw NormalError("mu_dual is not defined on the dual mesh"); + if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != i_dual_mesh) { + throw NormalError("dual variables are is not defined on the diamond dual of the primal mesh"); } m_energy_computer = @@ -1985,8 +1991,8 @@ EnergyComputerHandler::EnergyComputerHandler( std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh); - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mub->mesh()) { - throw NormalError("mu_dual is not defined on the dual mesh"); + if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != i_dual_mesh) { + throw NormalError("dual variables are is not defined on the diamond dual of the primal mesh"); } m_energy_computer = @@ -2006,8 +2012,8 @@ EnergyComputerHandler::EnergyComputerHandler( std::shared_ptr mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh); - if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != dual_mub->mesh()) { - throw NormalError("mu_dual is not defined on the dual mesh"); + if (DualMeshManager::instance().getDiamondDualMesh(*mesh) != i_dual_mesh) { + throw NormalError("dual variables are is not defined on the diamond dual of the primal mesh"); } m_energy_computer = -- GitLab