diff --git a/src/mesh/ItemArrayUtils.hpp b/src/mesh/ItemArrayUtils.hpp index 6afa39c167c9d2fa95ab5ef53e502df27359dce7..abba5cfe2631bf6c472f3a4c60a5c0eec8703acb 100644 --- a/src/mesh/ItemArrayUtils.hpp +++ b/src/mesh/ItemArrayUtils.hpp @@ -308,7 +308,7 @@ sum(const ItemArray<DataType, item_type, ConnectivityPtr>& item_value) template <typename DataType, ItemType item_type, typename ConnectivityPtr> void -synchronize(ItemArray<DataType, item_type, ConnectivityPtr>& item_array) +synchronize(ItemArray<DataType, item_type, ConnectivityPtr> item_array) { static_assert(not std::is_const_v<DataType>, "cannot synchronize ItemArray of const data"); if (parallel::size() > 1) { diff --git a/src/mesh/ItemValueUtils.hpp b/src/mesh/ItemValueUtils.hpp index ef99f59e44f9f71a01267ca0c424220241573f60..3bff4868c30677ef8d11db3a94cb493c5bcb07a0 100644 --- a/src/mesh/ItemValueUtils.hpp +++ b/src/mesh/ItemValueUtils.hpp @@ -296,7 +296,7 @@ sum(const ItemValue<DataType, item_type, ConnectivityPtr>& item_value) template <typename DataType, ItemType item_type, typename ConnectivityPtr> void -synchronize(ItemValue<DataType, item_type, ConnectivityPtr>& item_value) +synchronize(ItemValue<DataType, item_type, ConnectivityPtr> item_value) { static_assert(not std::is_const_v<DataType>, "cannot synchronize ItemValue of const data"); if (parallel::size() > 1) { @@ -309,7 +309,7 @@ synchronize(ItemValue<DataType, item_type, ConnectivityPtr>& item_value) template <typename DataType, ItemType item_type, typename ConnectivityPtr> bool -isSynchronized(const ItemValue<DataType, item_type, ConnectivityPtr>& item_value) +isSynchronized(ItemValue<DataType, item_type, ConnectivityPtr> item_value) { bool is_synchronized = true; diff --git a/src/scheme/FluxingAdvectionSolver.cpp b/src/scheme/FluxingAdvectionSolver.cpp index 27f724ff3d5509b7d7bb6c242a0c7e8615bbcb0c..ab9a0b138ed6b34b5cb4cdc6fa07d7207b052f3a 100644 --- a/src/scheme/FluxingAdvectionSolver.cpp +++ b/src/scheme/FluxingAdvectionSolver.cpp @@ -3,6 +3,8 @@ #include <language/utils/EvaluateAtPoints.hpp> #include <mesh/Connectivity.hpp> #include <mesh/IMesh.hpp> +#include <mesh/ItemArrayUtils.hpp> +#include <mesh/ItemValueUtils.hpp> #include <mesh/Mesh.hpp> #include <mesh/MeshData.hpp> #include <mesh/MeshDataManager.hpp> @@ -150,16 +152,19 @@ template <> FaceValue<double> FluxingAdvectionSolver<1>::_computeAlgebraicFluxingVolume() { - Array<double> algebraic_fluxing_volumes{m_new_mesh->numberOfNodes()}; - NodeValue<double> fluxing_volume(m_new_mesh->connectivity(), algebraic_fluxing_volumes); + Array<double> fluxing_volumes{m_new_mesh->numberOfNodes()}; + NodeValue<double> nodal_fluxing_volume(m_new_mesh->connectivity(), fluxing_volumes); auto old_xr = m_old_mesh->xr(); auto new_xr = m_new_mesh->xr(); parallel_for( m_new_mesh->numberOfNodes(), - PUGS_LAMBDA(NodeId node_id) { fluxing_volume[node_id] = new_xr[node_id][0] - old_xr[node_id][0]; }); + PUGS_LAMBDA(NodeId node_id) { nodal_fluxing_volume[node_id] = new_xr[node_id][0] - old_xr[node_id][0]; }); - return FaceValue<double>{m_new_mesh->connectivity(), algebraic_fluxing_volumes}; + FaceValue<double> algebraic_fluxing_volumes(m_new_mesh->connectivity(), fluxing_volumes); + + synchronize(algebraic_fluxing_volumes); + return algebraic_fluxing_volumes; } template <> @@ -185,6 +190,7 @@ FluxingAdvectionSolver<2>::_computeAlgebraicFluxingVolume() algebraic_fluxing_volume[face_id] = 0.5 * det(M); }); + synchronize(algebraic_fluxing_volume); return algebraic_fluxing_volume; } @@ -290,11 +296,12 @@ FluxingAdvectionSolver<Dimension>::_computeCycleNumber(FaceValue<double> fluxing MeshData<Dimension>& mesh_data = MeshDataManager::instance().getMeshData(*m_old_mesh); const CellValue<const double> Vj = mesh_data.Vj(); - const CellValue<size_t> ratio(m_old_mesh->connectivity()); + CellValue<size_t> ratio(m_old_mesh->connectivity()); parallel_for( m_old_mesh->numberOfCells(), PUGS_LAMBDA(CellId cell_id) { ratio[cell_id] = std::ceil(total_negative_flux[cell_id] / Vj[cell_id]); }); + synchronize(ratio); size_t number_of_cycles = max(ratio); @@ -379,6 +386,7 @@ FluxingAdvectionSolver<Dimension>::_remapOne(const CellValue<const double>& step } }); + synchronize(new_q); old_q = new_q; } @@ -411,6 +419,8 @@ FluxingAdvectionSolver<Dimension>::_remapAllQuantities() } }); + synchronize(step_Vj); + CellValue<double> inv_Vj(m_old_mesh->connectivity()); parallel_for( m_new_mesh->numberOfCells(), PUGS_LAMBDA(CellId cell_id) { inv_Vj[cell_id] = 1 / step_Vj[cell_id]; });