From 72f1136ddff58a6aca117de2f3d6e9b9e7b33bf1 Mon Sep 17 00:00:00 2001 From: Christoph Lehmann Date: Mon, 24 Oct 2022 16:40:13 +0200 Subject: [PATCH 1/6] [PL] Pass global matrix pointers to assemble fcts --- ProcessLib/VectorMatrixAssembler.cpp | 32 ++++++++++++++-------------- ProcessLib/VectorMatrixAssembler.h | 4 ++-- 2 files changed, 18 insertions(+), 18 deletions(-) diff --git a/ProcessLib/VectorMatrixAssembler.cpp b/ProcessLib/VectorMatrixAssembler.cpp index 99f44c61453..422e900926f 100644 --- a/ProcessLib/VectorMatrixAssembler.cpp +++ b/ProcessLib/VectorMatrixAssembler.cpp @@ -41,7 +41,7 @@ void VectorMatrixAssembler::assemble( std::vector const& dof_tables, const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b) + GlobalMatrix* M, GlobalMatrix* K, GlobalVector* b) { std::vector> indices_of_processes; indices_of_processes.reserve(dof_tables.size()); @@ -83,20 +83,20 @@ void VectorMatrixAssembler::assemble( auto const r_c_indices = NumLib::LocalToGlobalIndexMap::RowColumnIndices(indices, indices); - if (!_local_M_data.empty()) + if (M && !_local_M_data.empty()) { auto const local_M = MathLib::toMatrix(_local_M_data, num_r_c, num_r_c); - M.add(r_c_indices, local_M); + M->add(r_c_indices, local_M); } - if (!_local_K_data.empty()) + if (K && !_local_K_data.empty()) { auto const local_K = MathLib::toMatrix(_local_K_data, num_r_c, num_r_c); - K.add(r_c_indices, local_K); + K->add(r_c_indices, local_K); } - if (!_local_b_data.empty()) + if (b && !_local_b_data.empty()) { assert(_local_b_data.size() == num_r_c); - b.add(indices, _local_b_data); + b->add(indices, _local_b_data); } _local_output(t, process_id, mesh_item_id, _local_M_data, _local_K_data, @@ -108,7 +108,7 @@ void VectorMatrixAssembler::assembleWithJacobian( std::vector const& dof_tables, const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalMatrix* M, GlobalMatrix* K, GlobalVector* b, GlobalMatrix* Jac) { std::vector> indices_of_processes; indices_of_processes.reserve(dof_tables.size()); @@ -153,26 +153,26 @@ void VectorMatrixAssembler::assembleWithJacobian( auto const r_c_indices = NumLib::LocalToGlobalIndexMap::RowColumnIndices(indices, indices); - if (!_local_M_data.empty()) + if (M && !_local_M_data.empty()) { auto const local_M = MathLib::toMatrix(_local_M_data, num_r_c, num_r_c); - M.add(r_c_indices, local_M); + M->add(r_c_indices, local_M); } - if (!_local_K_data.empty()) + if (K && !_local_K_data.empty()) { auto const local_K = MathLib::toMatrix(_local_K_data, num_r_c, num_r_c); - K.add(r_c_indices, local_K); + K->add(r_c_indices, local_K); } - if (!_local_b_data.empty()) + if (b && !_local_b_data.empty()) { assert(_local_b_data.size() == num_r_c); - b.add(indices, _local_b_data); + b->add(indices, _local_b_data); } - if (!_local_Jac_data.empty()) + if (Jac && !_local_Jac_data.empty()) { auto const local_Jac = MathLib::toMatrix(_local_Jac_data, num_r_c, num_r_c); - Jac.add(r_c_indices, local_Jac); + Jac->add(r_c_indices, local_Jac); } else { diff --git a/ProcessLib/VectorMatrixAssembler.h b/ProcessLib/VectorMatrixAssembler.h index b96fb2219af..2bd58de8a05 100644 --- a/ProcessLib/VectorMatrixAssembler.h +++ b/ProcessLib/VectorMatrixAssembler.h @@ -48,7 +48,7 @@ class VectorMatrixAssembler final std::vector const& dof_tables, double const t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b); + GlobalMatrix* M, GlobalMatrix* K, GlobalVector* b); //! Assembles \c M, \c K, \c b, and the Jacobian \c Jac of the residual. //! \note The Jacobian must be assembled. @@ -58,7 +58,7 @@ class VectorMatrixAssembler final std::vector const& dof_tables, const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac); + GlobalMatrix* M, GlobalMatrix* K, GlobalVector* b, GlobalMatrix* Jac); private: // temporary data only stored here in order to avoid frequent memory From 950934f5b3921d39d6d584c086b13b68cd0e1cd7 Mon Sep 17 00:00:00 2001 From: Dmitri Naumov Date: Mon, 24 Oct 2022 16:40:13 +0200 Subject: [PATCH 2/6] [PL] Update assemble* signatures to use pointers --- ProcessLib/Assembly/AssembledMatrixCache.h | 22 +++++++++---------- .../ComponentTransportProcess.cpp | 6 ++--- ProcessLib/HT/HTProcess.cpp | 6 ++--- .../HeatConduction/HeatConductionProcess.cpp | 4 ++-- .../HeatTransportBHEProcess.cpp | 6 ++--- .../HydroMechanics/HydroMechanicsProcess.cpp | 6 ++--- .../HydroMechanics/HydroMechanicsProcess.cpp | 4 ++-- .../SmallDeformationProcess.cpp | 6 ++--- .../LargeDeformationProcess.cpp | 6 ++--- ProcessLib/LiquidFlow/LiquidFlowProcess.cpp | 6 ++--- ProcessLib/PhaseField/PhaseFieldProcess.cpp | 6 ++--- .../RichardsComponentTransportProcess.cpp | 6 ++--- .../RichardsFlow/RichardsFlowProcess.cpp | 6 ++--- .../RichardsMechanicsProcess.cpp | 6 ++--- .../SmallDeformationProcess.cpp | 6 ++--- .../SmallDeformationNonlocalProcess.cpp | 6 ++--- .../SteadyStateDiffusion.cpp | 6 ++--- ProcessLib/StokesFlow/StokesFlowProcess.cpp | 4 ++-- ProcessLib/TES/TESProcess.cpp | 6 ++--- .../ThermalTwoPhaseFlowWithPPProcess.cpp | 6 ++--- .../ThermoHydroMechanicsProcess.cpp | 6 ++--- .../ThermoMechanicalPhaseFieldProcess.cpp | 6 ++--- .../ThermoMechanicsProcess.cpp | 6 ++--- .../ThermoRichardsFlowProcess.cpp | 4 ++-- .../TwoPhaseFlowWithPPProcess.cpp | 6 ++--- .../TwoPhaseFlowWithPrhoProcess.cpp | 6 ++--- 26 files changed, 82 insertions(+), 82 deletions(-) diff --git a/ProcessLib/Assembly/AssembledMatrixCache.h b/ProcessLib/Assembly/AssembledMatrixCache.h index fdc1912b422..55fa3bd3cd2 100644 --- a/ProcessLib/Assembly/AssembledMatrixCache.h +++ b/ProcessLib/Assembly/AssembledMatrixCache.h @@ -33,7 +33,7 @@ struct AssembledMatrixCache final void assemble( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, + GlobalMatrix* const M, GlobalMatrix* const K, GlobalVector* const b, std::vector const& dof_tables, VectorMatrixAssembler& global_assembler, VectorOfLocalAssemblers const& local_assemblers, @@ -53,7 +53,7 @@ template void AssembledMatrixCache::assemble( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, + GlobalMatrix* const M, GlobalMatrix* const K, GlobalVector* const b, std::vector const& dof_tables, VectorMatrixAssembler& global_assembler, VectorOfLocalAssemblers const& local_assemblers, @@ -70,9 +70,9 @@ void AssembledMatrixCache::assemble( local_assemblers, active_element_ids, dof_tables, t, dt, x, x_prev, process_id, M, K, b); - MathLib::finalizeMatrixAssembly(M); - MathLib::finalizeMatrixAssembly(K); - MathLib::finalizeVectorAssembly(b); + MathLib::finalizeMatrixAssembly(*M); + MathLib::finalizeMatrixAssembly(*K); + MathLib::finalizeVectorAssembly(*b); INFO("[time] Calling local assemblers took {:g} s", time_asm.elapsed()); @@ -83,9 +83,9 @@ void AssembledMatrixCache::assemble( BaseLib::RunTime time_save; time_save.start(); - K_ = MathLib::MatrixVectorTraits::newInstance(K); - M_ = MathLib::MatrixVectorTraits::newInstance(M); - b_ = MathLib::MatrixVectorTraits::newInstance(b); + K_ = MathLib::MatrixVectorTraits::newInstance(*K); + M_ = MathLib::MatrixVectorTraits::newInstance(*M); + b_ = MathLib::MatrixVectorTraits::newInstance(*b); INFO("[time] Saving global K, M, b took {:g} s", time_save.elapsed()); @@ -98,9 +98,9 @@ void AssembledMatrixCache::assemble( BaseLib::RunTime time_restore; time_restore.start(); - MathLib::LinAlg::copy(*K_, K); - MathLib::LinAlg::copy(*M_, M); - MathLib::LinAlg::copy(*b_, b); + MathLib::LinAlg::copy(*K_, *K); + MathLib::LinAlg::copy(*M_, *M); + MathLib::LinAlg::copy(*b_, *b); INFO("[time] Restoring global K, M, b took {:g} s", time_restore.elapsed()); diff --git a/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp b/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp index f2fab9455ff..55c2cdf83c1 100644 --- a/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp +++ b/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp @@ -256,8 +256,8 @@ void ComponentTransportProcess::assembleConcreteProcess( [&]() { return _local_to_global_index_map.get(); }); } - _asm_mat_cache.assemble(t, dt, x, x_prev, process_id, M, K, b, dof_tables, - _global_assembler, _local_assemblers, + _asm_mat_cache.assemble(t, dt, x, x_prev, process_id, &M, &K, &b, + dof_tables, _global_assembler, _local_assemblers, getActiveElementIDs()); // TODO (naumov) What about temperature? A test with FCT would reveal any @@ -296,7 +296,7 @@ void ComponentTransportProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); // b is the negated residumm used in the Newton's method. // Here negating b is to recover the primitive residuum. diff --git a/ProcessLib/HT/HTProcess.cpp b/ProcessLib/HT/HTProcess.cpp index e470cfc5e75..63d6ddfcc24 100644 --- a/ProcessLib/HT/HTProcess.cpp +++ b/ProcessLib/HT/HTProcess.cpp @@ -105,8 +105,8 @@ void HTProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_tables, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_tables, t, dt, x, x_prev, process_id, &M, &K, + &b); } void HTProcess::assembleWithJacobianConcreteProcess( @@ -131,7 +131,7 @@ void HTProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } std::tuple diff --git a/ProcessLib/HeatConduction/HeatConductionProcess.cpp b/ProcessLib/HeatConduction/HeatConductionProcess.cpp index 15af7514aed..a12e8a11fa0 100644 --- a/ProcessLib/HeatConduction/HeatConductionProcess.cpp +++ b/ProcessLib/HeatConduction/HeatConductionProcess.cpp @@ -103,7 +103,7 @@ void HeatConductionProcess::assembleConcreteProcess( std::vector dof_table = { _local_to_global_index_map.get()}; - _asm_mat_cache.assemble(t, dt, x, x_prev, process_id, M, K, b, dof_table, + _asm_mat_cache.assemble(t, dt, x, x_prev, process_id, &M, &K, &b, dof_table, _global_assembler, _local_assemblers, getActiveElementIDs()); } @@ -121,7 +121,7 @@ void HeatConductionProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); transformVariableFromGlobalVector(b, 0 /*variable id*/, *_local_to_global_index_map, *_heat_flux, diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp index 4e50103fd3e..03878612a98 100644 --- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp +++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp @@ -170,8 +170,8 @@ void HeatTransportBHEProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } void HeatTransportBHEProcess::assembleWithJacobianConcreteProcess( @@ -188,7 +188,7 @@ void HeatTransportBHEProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } void HeatTransportBHEProcess::computeSecondaryVariableConcrete( diff --git a/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp b/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp index 7219b8596b6..2addf268dd6 100644 --- a/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp +++ b/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp @@ -299,8 +299,8 @@ void HydroMechanicsProcess::assembleConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } template @@ -339,7 +339,7 @@ void HydroMechanicsProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); auto copyRhs = [&](int const variable_id, auto& output_vector) { diff --git a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp index adb2a187374..976676c913d 100644 --- a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp +++ b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp @@ -454,7 +454,7 @@ void HydroMechanicsProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - dof_table, t, dt, x, x_prev, process_id, M, K, b); + dof_table, t, dt, x, x_prev, process_id, &M, &K, &b); } template @@ -471,7 +471,7 @@ void HydroMechanicsProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); auto copyRhs = [&](int const variable_id, auto& output_vector) { diff --git a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp index b5284774aa1..6cf8ae89652 100644 --- a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp +++ b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp @@ -421,8 +421,8 @@ void SmallDeformationProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } template void SmallDeformationProcess:: @@ -439,7 +439,7 @@ void SmallDeformationProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } template void SmallDeformationProcess::preTimestepConcreteProcess( diff --git a/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp b/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp index 27757556bc6..169ad6ed5e8 100644 --- a/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp +++ b/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp @@ -142,8 +142,8 @@ void LargeDeformationProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); _global_output(t, process_id, M, K, b); } @@ -163,7 +163,7 @@ void LargeDeformationProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); transformVariableFromGlobalVector(b, 0, *_local_to_global_index_map, *_nodal_forces, std::negate()); diff --git a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp index 220050a8cd7..af379eab6fc 100644 --- a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp +++ b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp @@ -81,8 +81,8 @@ void LiquidFlowProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); MathLib::finalizeMatrixAssembly(M); MathLib::finalizeMatrixAssembly(K); @@ -108,7 +108,7 @@ void LiquidFlowProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } void LiquidFlowProcess::computeSecondaryVariableConcrete( diff --git a/ProcessLib/PhaseField/PhaseFieldProcess.cpp b/ProcessLib/PhaseField/PhaseFieldProcess.cpp index e4b3c29d583..1a35046a8fd 100644 --- a/ProcessLib/PhaseField/PhaseFieldProcess.cpp +++ b/ProcessLib/PhaseField/PhaseFieldProcess.cpp @@ -210,8 +210,8 @@ void PhaseFieldProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_tables, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_tables, t, dt, x, x_prev, process_id, &M, &K, + &b); } template @@ -242,7 +242,7 @@ void PhaseFieldProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); if (process_id == 0) { diff --git a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp index 34dc2483a02..f034a1c1360 100644 --- a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp +++ b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp @@ -77,8 +77,8 @@ void RichardsComponentTransportProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } void RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess( @@ -95,7 +95,7 @@ void RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } } // namespace RichardsComponentTransport diff --git a/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp b/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp index 01d28ff16c3..e00e2ac3da6 100644 --- a/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp +++ b/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp @@ -71,8 +71,8 @@ void RichardsFlowProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } void RichardsFlowProcess::assembleWithJacobianConcreteProcess( @@ -89,7 +89,7 @@ void RichardsFlowProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } } // namespace RichardsFlow diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp index 4ab85ce74c2..d04fe1f174a 100644 --- a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp +++ b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp @@ -283,8 +283,8 @@ void RichardsMechanicsProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } template @@ -322,7 +322,7 @@ void RichardsMechanicsProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); auto copyRhs = [&](int const variable_id, auto& output_vector) { diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp index fbeff5e241b..dc768172fa7 100644 --- a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp +++ b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp @@ -144,8 +144,8 @@ void SmallDeformationProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); _global_output(t, process_id, M, K, b); } @@ -165,7 +165,7 @@ void SmallDeformationProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); transformVariableFromGlobalVector(b, 0, *_local_to_global_index_map, *_nodal_forces, std::negate()); diff --git a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp index 36aacd014d2..cbf799fa300 100644 --- a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp +++ b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp @@ -189,8 +189,8 @@ void SmallDeformationNonlocalProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } template @@ -224,7 +224,7 @@ void SmallDeformationNonlocalProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); b.copyValues(*_nodal_forces); std::transform(_nodal_forces->begin(), _nodal_forces->end(), diff --git a/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.cpp b/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.cpp index 9c9071f5ff0..1fcf8b7de10 100644 --- a/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.cpp +++ b/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.cpp @@ -67,8 +67,8 @@ void SteadyStateDiffusion::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } void SteadyStateDiffusion::assembleWithJacobianConcreteProcess( @@ -84,7 +84,7 @@ void SteadyStateDiffusion::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } } // namespace SteadyStateDiffusion diff --git a/ProcessLib/StokesFlow/StokesFlowProcess.cpp b/ProcessLib/StokesFlow/StokesFlowProcess.cpp index 488a9719716..3fc37e549b2 100644 --- a/ProcessLib/StokesFlow/StokesFlowProcess.cpp +++ b/ProcessLib/StokesFlow/StokesFlowProcess.cpp @@ -143,8 +143,8 @@ void StokesFlowProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_tables, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_tables, t, dt, x, x_prev, process_id, &M, &K, + &b); } template diff --git a/ProcessLib/TES/TESProcess.cpp b/ProcessLib/TES/TESProcess.cpp index 00b1ab5dc87..466026afe15 100644 --- a/ProcessLib/TES/TESProcess.cpp +++ b/ProcessLib/TES/TESProcess.cpp @@ -216,8 +216,8 @@ void TESProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } void TESProcess::assembleWithJacobianConcreteProcess( @@ -232,7 +232,7 @@ void TESProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } void TESProcess::preTimestepConcreteProcess(std::vector const& x, diff --git a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp index 8716c977ec5..474a0426a91 100644 --- a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp +++ b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp @@ -93,8 +93,8 @@ void ThermalTwoPhaseFlowWithPPProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); _global_output(t, process_id, M, K, b); } @@ -113,7 +113,7 @@ void ThermalTwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); _global_output(t, process_id, M, K, b, &Jac); } diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp index 7087f1570aa..bd2f596277a 100644 --- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp +++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp @@ -346,8 +346,8 @@ void ThermoHydroMechanicsProcess::assembleConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } template @@ -392,7 +392,7 @@ void ThermoHydroMechanicsProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); auto copyRhs = [&](int const variable_id, auto& output_vector) { diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp index 3d3830ff152..72a414084ac 100644 --- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp +++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp @@ -204,8 +204,8 @@ void ThermoMechanicalPhaseFieldProcess:: // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } template @@ -247,7 +247,7 @@ void ThermoMechanicalPhaseFieldProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } template diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp index 680c2f460fc..21fa72010ca 100644 --- a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp +++ b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp @@ -216,8 +216,8 @@ void ThermoMechanicsProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } template @@ -273,7 +273,7 @@ void ThermoMechanicsProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); // TODO (naumov): Refactor the copy rhs part. This is copy from HM. auto copyRhs = [&](int const variable_id, auto& output_vector) diff --git a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp index 665dc126c5a..8c88b5ddc93 100644 --- a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp +++ b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp @@ -137,8 +137,8 @@ void ThermoRichardsFlowProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } void ThermoRichardsFlowProcess::assembleWithJacobianConcreteProcess( diff --git a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp index 29de8ba6841..2ea49f2db01 100644 --- a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp +++ b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp @@ -71,8 +71,8 @@ void TwoPhaseFlowWithPPProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } void TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess( @@ -89,7 +89,7 @@ void TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } } // namespace TwoPhaseFlowWithPP diff --git a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp index 05944f93e53..1bbd135aafc 100644 --- a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp +++ b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp @@ -71,8 +71,8 @@ void TwoPhaseFlowWithPrhoProcess::assembleConcreteProcess( // Call global assembler for each local assembly item. GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assemble, _local_assemblers, - getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, M, K, - b); + getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, + &b); } void TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess( @@ -89,7 +89,7 @@ void TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, M, K, b, Jac); + process_id, &M, &K, &b, &Jac); } } // namespace TwoPhaseFlowWithPrho From 3b1362c1880f5440182428dee0cc63b34f5638e0 Mon Sep 17 00:00:00 2001 From: Dmitri Naumov Date: Wed, 17 Jul 2024 11:23:33 +0200 Subject: [PATCH 3/6] [PL] Only asmJacobian is implemented for asmMixins The Picard implementation would be slightly different then the current assembleGeneric suggests anyway. Simplifying the code, small cleanups. --- ProcessLib/AssemblyMixin.h | 55 ++++++++++++++------------------------ 1 file changed, 20 insertions(+), 35 deletions(-) diff --git a/ProcessLib/AssemblyMixin.h b/ProcessLib/AssemblyMixin.h index 2d00a31f7df..b082a5873ff 100644 --- a/ProcessLib/AssemblyMixin.h +++ b/ProcessLib/AssemblyMixin.h @@ -137,23 +137,22 @@ class AssemblyMixin : private AssemblyMixinBase } // cppcheck-suppress functionStatic - void assemble(const double /*t*/, double const /*dt*/, + void assemble(double const t, double const dt, std::vector const& /*x*/, std::vector const& /*x_prev*/, - int const /*process_id*/, GlobalMatrix& /*M*/, + int const process_id, GlobalMatrix& /*M*/, GlobalMatrix& /*K*/, GlobalVector& /*b*/) { - /* DBUG("AssemblyMixin assemble(t={}, dt={}, process_id={}).", t, dt, process_id); - assembleGeneric(&Assembly::ParallelVectorMatrixAssembler::assemble, t, - dt, x, x_prev, process_id, M, K, b); - */ - OGS_FATAL("Not yet implemented."); + /// Implementation similar to assembleWithJacobian calling + /// Assembly::ParallelVectorMatrixAssembler::assemble function. + /// Residuum must be correctly computed. + OGS_FATAL("AssemblyMixin for Picard scheme is not yet implemented."); } - void assembleWithJacobian(const double t, double const dt, + void assembleWithJacobian(double const t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, GlobalMatrix& M, @@ -163,27 +162,6 @@ class AssemblyMixin : private AssemblyMixinBase DBUG("AssemblyMixin assembleWithJacobian(t={}, dt={}, process_id={}).", t, dt, process_id); - assembleGeneric( - &Assembly::ParallelVectorMatrixAssembler::assembleWithJacobian, t, - dt, x, x_prev, process_id, M, K, b, Jac); - } - -private: - Process& derived() { return static_cast(*this); } - Process const& derived() const - { - return static_cast(*this); - } - - /// Generic assembly routine covering both the case with and without - /// Jacobian assembly. - template - void assembleGeneric(Method global_assembler_method, const double t, - double const dt, std::vector const& x, - std::vector const& x_prev, - int const process_id, GlobalMatrix& M, GlobalMatrix& K, - GlobalVector& b, Jac&... jac_or_not_jac) - { // TODO why not getDOFTables(x.size()); ? std::vector const dof_tables{ derived()._local_to_global_index_map.get()}; @@ -199,9 +177,9 @@ class AssemblyMixin : private AssemblyMixinBase { b_submesh.setZero(); - (pvma_.*global_assembler_method)( - loc_asms, sad.active_element_ids, dof_tables, t, dt, x, - x_prev, process_id, M, K, b_submesh, jac_or_not_jac...); + pvma_.assembleWithJacobian(loc_asms, sad.active_element_ids, + dof_tables, t, dt, x, x_prev, + process_id, M, K, b_submesh, Jac); MathLib::LinAlg::axpy(b, 1.0, b_submesh); @@ -218,14 +196,21 @@ class AssemblyMixin : private AssemblyMixinBase ProcessLib::ProcessVariable const& pv = derived().getProcessVariables(process_id)[0]; - (pvma_.*global_assembler_method)( - loc_asms, pv.getActiveElementIDs(), dof_tables, t, dt, x, - x_prev, process_id, M, K, b, jac_or_not_jac...); + pvma_.assembleWithJacobian(loc_asms, pv.getActiveElementIDs(), + dof_tables, t, dt, x, x_prev, process_id, + M, K, b, Jac); } AssemblyMixinBase::copyResiduumVectorsToBulkMesh( b, *(dof_tables.front()), residuum_vectors_bulk_); } + +private: + Process& derived() { return static_cast(*this); } + Process const& derived() const + { + return static_cast(*this); + } }; } // namespace ProcessLib From 4466324b697b0687f0c8f87618513650d9128927 Mon Sep 17 00:00:00 2001 From: Dmitri Naumov Date: Mon, 15 Jul 2024 21:26:01 +0200 Subject: [PATCH 4/6] Remove M and K from Newton's assembleWithJacobian The matrices are not in use by any process. --- NumLib/ODESolver/ODESystem.h | 41 +----- NumLib/ODESolver/TimeDiscretizedODESystem.cpp | 32 ++--- NumLib/ODESolver/TimeDiscretizedODESystem.h | 4 - ProcessLib/AbstractJacobianAssembler.h | 4 - ProcessLib/AnalyticalJacobianAssembler.cpp | 8 +- ProcessLib/AnalyticalJacobianAssembler.h | 3 - ProcessLib/Assembly/MatrixAssemblyStats.h | 6 - ProcessLib/Assembly/MatrixElementCache.h | 16 +-- ProcessLib/Assembly/MatrixOutput.cpp | 22 ++-- ProcessLib/Assembly/MatrixOutput.h | 8 +- .../ParallelVectorMatrixAssembler.cpp | 38 ++---- .../Assembly/ParallelVectorMatrixAssembler.h | 2 +- ProcessLib/AssemblyMixin.h | 7 +- .../BoundaryCondition.h | 2 +- .../BoundaryConditionCollection.cpp | 2 +- .../BoundaryConditionCollection.h | 2 +- .../GenericNaturalBoundaryCondition-impl.h | 2 +- .../GenericNaturalBoundaryCondition.h | 2 +- ...icNaturalBoundaryConditionLocalAssembler.h | 2 +- ...onentFlowBoundaryConditionLocalAssembler.h | 2 +- .../NeumannBoundaryConditionLocalAssembler.h | 2 +- .../NormalTractionBoundaryCondition-impl.h | 2 +- .../NormalTractionBoundaryCondition.h | 2 +- ...lTractionBoundaryConditionLocalAssembler.h | 4 +- .../Python/PythonBoundaryCondition.cpp | 2 +- .../Python/PythonBoundaryCondition.h | 2 +- .../PythonBoundaryConditionLocalAssembler.h | 2 +- .../RobinBoundaryConditionLocalAssembler.h | 14 ++- ...ntNeumannBoundaryConditionLocalAssembler.h | 2 +- .../CentralDifferencesJacobianAssembler.cpp | 4 +- .../CentralDifferencesJacobianAssembler.h | 2 - .../CompareJacobiansJacobianAssembler.cpp | 62 +-------- .../CompareJacobiansJacobianAssembler.h | 2 - .../ForwardDifferencesJacobianAssembler.cpp | 4 +- .../ForwardDifferencesJacobianAssembler.h | 2 - ProcessLib/LocalAssemblerInterface.cpp | 4 - ProcessLib/LocalAssemblerInterface.h | 3 - ProcessLib/Process.cpp | 12 +- ProcessLib/Process.h | 6 +- .../SteadyStateDiffusion.cpp | 4 +- .../SteadyStateDiffusion.h | 3 +- ProcessLib/VectorMatrixAssembler.cpp | 26 ++-- ProcessLib/VectorMatrixAssembler.h | 4 +- Tests/NumLib/ODEs.h | 118 +++++++++++------- Tests/ProcessLib/TestJacobianAssembler.cpp | 69 ++-------- 45 files changed, 192 insertions(+), 370 deletions(-) diff --git a/NumLib/ODESolver/ODESystem.h b/NumLib/ODESolver/ODESystem.h index 62eaf4bf448..6625978d341 100644 --- a/NumLib/ODESolver/ODESystem.h +++ b/NumLib/ODESolver/ODESystem.h @@ -103,55 +103,18 @@ class ODESystem const& x, std::vector const& x_prev, - int const process_id, GlobalMatrix& M, - GlobalMatrix& K, GlobalVector& b, + int const process_id, GlobalVector& b, GlobalMatrix& Jac) = 0; }; diff --git a/NumLib/ODESolver/TimeDiscretizedODESystem.cpp b/NumLib/ODESolver/TimeDiscretizedODESystem.cpp index f45e91ffbef..a536e2ff034 100644 --- a/NumLib/ODESolver/TimeDiscretizedODESystem.cpp +++ b/NumLib/ODESolver/TimeDiscretizedODESystem.cpp @@ -55,10 +55,6 @@ TimeDiscretizedODESystem::~TimeDiscretizedODESystem() { NumLib::GlobalMatrixProvider::provider.releaseMatrix(*_Jac); - NumLib::GlobalMatrixProvider::provider.releaseMatrix(*_M); - NumLib::GlobalMatrixProvider::provider.releaseMatrix(*_K); NumLib::GlobalVectorProvider::provider.releaseVector(*_b); } @@ -79,35 +73,29 @@ void TimeDiscretizedODESystem const& x_prev, int const process_id) { - namespace LinAlg = MathLib::LinAlg; - auto const t = _time_disc.getCurrentTime(); auto const dt = _time_disc.getCurrentTimeIncrement(); auto const& x_curr = *x_new_timestep[process_id]; - _M->setZero(); - _K->setZero(); _b->setZero(); _Jac->setZero(); _ode.preAssemble(t, dt, x_curr); - _ode.assembleWithJacobian(t, dt, x_new_timestep, x_prev, process_id, *_M, - *_K, *_b, *_Jac); + _ode.assembleWithJacobian(t, dt, x_new_timestep, x_prev, process_id, *_b, + *_Jac); - LinAlg::finalizeAssembly(*_M); - LinAlg::finalizeAssembly(*_K); - LinAlg::finalizeAssembly(*_b); + MathLib::LinAlg::finalizeAssembly(*_b); MathLib::LinAlg::finalizeAssembly(*_Jac); } -void TimeDiscretizedODESystem< - ODESystemTag::FirstOrderImplicitQuasilinear, - NonlinearSolverTag::Newton>::getResidual(GlobalVector const& x_new_timestep, - GlobalVector const& x_prev, - GlobalVector& res) const +void TimeDiscretizedODESystem:: + getResidual(GlobalVector const& /*x_new_timestep*/, + GlobalVector const& /*x_prev*/, + GlobalVector& res) const { - double const dt = _time_disc.getCurrentTimeIncrement(); - _mat_trans->computeResidual(*_M, *_K, *_b, dt, x_new_timestep, x_prev, res); + MathLib::LinAlg::copy(*_b, res); // res = b + MathLib::LinAlg::scale(res, -1.); // res = -b } void TimeDiscretizedODESystem< diff --git a/NumLib/ODESolver/TimeDiscretizedODESystem.h b/NumLib/ODESolver/TimeDiscretizedODESystem.h index 66dfc810002..4fffd07cbcb 100644 --- a/NumLib/ODESolver/TimeDiscretizedODESystem.h +++ b/NumLib/ODESolver/TimeDiscretizedODESystem.h @@ -149,13 +149,9 @@ class TimeDiscretizedODESystem const& local_x, std::vector const& local_x_prev, - std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) = 0; @@ -40,8 +38,6 @@ class AbstractJacobianAssembler LocalAssemblerInterface& /*local_assembler*/, double const /*t*/, double const /*dt*/, Eigen::VectorXd const& /*local_x*/, Eigen::VectorXd const& /*local_x_prev*/, int const /*process_id*/, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& /*local_b_data*/, std::vector& /*local_Jac_data*/) { diff --git a/ProcessLib/AnalyticalJacobianAssembler.cpp b/ProcessLib/AnalyticalJacobianAssembler.cpp index 8ee648ffc0e..d908f3beb1d 100644 --- a/ProcessLib/AnalyticalJacobianAssembler.cpp +++ b/ProcessLib/AnalyticalJacobianAssembler.cpp @@ -17,24 +17,20 @@ namespace ProcessLib void AnalyticalJacobianAssembler::assembleWithJacobian( LocalAssemblerInterface& local_assembler, double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) { local_assembler.assembleWithJacobian(t, dt, local_x, local_x_prev, - local_M_data, local_K_data, local_b_data, local_Jac_data); } void AnalyticalJacobianAssembler::assembleWithJacobianForStaggeredScheme( LocalAssemblerInterface& local_assembler, double const t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& local_x_prev, - int const process_id, std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, + int const process_id, std::vector& local_b_data, std::vector& local_Jac_data) { local_assembler.assembleWithJacobianForStaggeredScheme( - t, dt, local_x, local_x_prev, process_id, local_M_data, local_K_data, - local_b_data, local_Jac_data); + t, dt, local_x, local_x_prev, process_id, local_b_data, local_Jac_data); } std::unique_ptr AnalyticalJacobianAssembler::copy() diff --git a/ProcessLib/AnalyticalJacobianAssembler.h b/ProcessLib/AnalyticalJacobianAssembler.h index a189521ee44..95d51699199 100644 --- a/ProcessLib/AnalyticalJacobianAssembler.h +++ b/ProcessLib/AnalyticalJacobianAssembler.h @@ -33,8 +33,6 @@ class AnalyticalJacobianAssembler final : public AbstractJacobianAssembler double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; @@ -42,7 +40,6 @@ class AnalyticalJacobianAssembler final : public AbstractJacobianAssembler LocalAssemblerInterface& local_assembler, double const t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/Assembly/MatrixAssemblyStats.h b/ProcessLib/Assembly/MatrixAssemblyStats.h index f2087855469..5be93732ea6 100644 --- a/ProcessLib/Assembly/MatrixAssemblyStats.h +++ b/ProcessLib/Assembly/MatrixAssemblyStats.h @@ -47,15 +47,11 @@ struct Stats struct MultiStats { - Stats M; - Stats K; Stats b; Stats Jac; MultiStats& operator+=(MultiStats const& other) { - M += other.M; - K += other.K; b += other.b; Jac += other.Jac; @@ -64,8 +60,6 @@ struct MultiStats void print() const { - M.print("M"); - K.print("K"); b.print("b"); Jac.print("J"); } diff --git a/ProcessLib/Assembly/MatrixElementCache.h b/ProcessLib/Assembly/MatrixElementCache.h index e7505aafedb..355ce5c1381 100644 --- a/ProcessLib/Assembly/MatrixElementCache.h +++ b/ProcessLib/Assembly/MatrixElementCache.h @@ -257,31 +257,21 @@ class MultiMatrixElementCache final using GlobalVectorView = ConcurrentMatrixView<1>; public: - MultiMatrixElementCache(GlobalMatrixView& M, GlobalMatrixView& K, - GlobalVectorView& b, GlobalMatrixView& Jac, + MultiMatrixElementCache(GlobalVectorView& b, GlobalMatrixView& Jac, MultiStats& stats) - : cache_M_(M, stats.M), - cache_K_(K, stats.K), - cache_b_(b, stats.b), - cache_Jac_(Jac, stats.Jac) + : cache_b_(b, stats.b), cache_Jac_(Jac, stats.Jac) { } - void add(std::vector const& local_M_data, - std::vector const& local_K_data, - std::vector const& local_b_data, + void add(std::vector const& local_b_data, std::vector const& local_Jac_data, std::vector const& indices) { - cache_M_.add(local_M_data, indices); - cache_K_.add(local_K_data, indices); cache_b_.add(local_b_data, indices); cache_Jac_.add(local_Jac_data, indices); } private: - MatrixElementCache<2> cache_M_; - MatrixElementCache<2> cache_K_; MatrixElementCache<1> cache_b_; MatrixElementCache<2> cache_Jac_; }; diff --git a/ProcessLib/Assembly/MatrixOutput.cpp b/ProcessLib/Assembly/MatrixOutput.cpp index 46d17cd0f95..6e01624cd83 100644 --- a/ProcessLib/Assembly/MatrixOutput.cpp +++ b/ProcessLib/Assembly/MatrixOutput.cpp @@ -175,8 +175,8 @@ GlobalMatrixOutput::GlobalMatrixOutput() } void GlobalMatrixOutput::operator()(double const t, int const process_id, - GlobalMatrix const& M, - GlobalMatrix const& K, + GlobalMatrix const* M, + GlobalMatrix const* K, GlobalVector const& b, GlobalMatrix const* const Jac) { @@ -188,20 +188,22 @@ void GlobalMatrixOutput::operator()(double const t, int const process_id, #ifndef USE_PETSC ++counter_; + if (M) { auto fh = openGlobalMatrixOutputFile(filenamePrefix_, counter_, t, process_id, "M", "mat"); fh << "M "; - outputGlobalMatrix(M, fh); + outputGlobalMatrix(*M, fh); } + if (K) { auto fh = openGlobalMatrixOutputFile(filenamePrefix_, counter_, t, process_id, "K", "mat"); fh << "K "; - outputGlobalMatrix(K, fh); + outputGlobalMatrix(*K, fh); } { @@ -286,8 +288,8 @@ LocalMatrixOutput::LocalMatrixOutput() void LocalMatrixOutput::operator()( double const t, int const process_id, std::size_t const element_id, - std::vector const& local_M_data, - std::vector const& local_K_data, + std::vector const* local_M_data, + std::vector const* local_K_data, std::vector const& local_b_data, std::vector const* const local_Jac_data) { @@ -305,16 +307,16 @@ void LocalMatrixOutput::operator()( fmt::print(fh, "## t = {:.15g}, process id = {}, element id = {}\n\n", t, process_id, element_id); - if (!local_M_data.empty()) + if (local_M_data && !local_M_data->empty()) { DBUG("... M"); - fmt::print(fh, "# M\n{}\n\n", toSquareMatrixRowMajor(local_M_data)); + fmt::print(fh, "# M\n{}\n\n", toSquareMatrixRowMajor(*local_M_data)); } - if (!local_K_data.empty()) + if (local_K_data && !local_K_data->empty()) { DBUG("... K"); - fmt::print(fh, "# K\n{}\n\n", toSquareMatrixRowMajor(local_K_data)); + fmt::print(fh, "# K\n{}\n\n", toSquareMatrixRowMajor(*local_K_data)); } if (!local_b_data.empty()) diff --git a/ProcessLib/Assembly/MatrixOutput.h b/ProcessLib/Assembly/MatrixOutput.h index 6ade067497a..aee4e411932 100644 --- a/ProcessLib/Assembly/MatrixOutput.h +++ b/ProcessLib/Assembly/MatrixOutput.h @@ -24,8 +24,8 @@ struct GlobalMatrixOutput { GlobalMatrixOutput(); - void operator()(double const t, int const process_id, GlobalMatrix const& M, - GlobalMatrix const& K, GlobalVector const& b, + void operator()(double const t, int const process_id, GlobalMatrix const* M, + GlobalMatrix const* K, GlobalVector const& b, GlobalMatrix const* const Jac = nullptr); private: @@ -45,8 +45,8 @@ struct LocalMatrixOutput void operator()(double const t, int const process_id, std::size_t const element_id, - std::vector const& local_M_data, - std::vector const& local_K_data, + std::vector const* local_M_data, + std::vector const* local_K_data, std::vector const& local_b_data, std::vector const* const local_Jac_data = nullptr); diff --git a/ProcessLib/Assembly/ParallelVectorMatrixAssembler.cpp b/ProcessLib/Assembly/ParallelVectorMatrixAssembler.cpp index 18d0fa334ee..ee4d4f94b8f 100644 --- a/ProcessLib/Assembly/ParallelVectorMatrixAssembler.cpp +++ b/ProcessLib/Assembly/ParallelVectorMatrixAssembler.cpp @@ -27,7 +27,6 @@ void assembleWithJacobianOneElement( ProcessLib::LocalAssemblerInterface& local_assembler, const NumLib::LocalToGlobalIndexMap& dof_table, const double t, const double dt, const GlobalVector& x, const GlobalVector& x_prev, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data, std::vector& indices, ProcessLib::AbstractJacobianAssembler& jacobian_assembler, @@ -35,16 +34,14 @@ void assembleWithJacobianOneElement( { indices = NumLib::getIndices(mesh_item_id, dof_table); - local_M_data.clear(); - local_K_data.clear(); local_b_data.clear(); local_Jac_data.clear(); auto const local_x = x.get(indices); auto const local_x_prev = x_prev.get(indices); - jacobian_assembler.assembleWithJacobian( - local_assembler, t, dt, local_x, local_x_prev, local_M_data, - local_K_data, local_b_data, local_Jac_data); + jacobian_assembler.assembleWithJacobian(local_assembler, t, dt, local_x, + local_x_prev, local_b_data, + local_Jac_data); if (local_Jac_data.empty()) { @@ -54,8 +51,7 @@ void assembleWithJacobianOneElement( "current process."); } - cache.add(local_M_data, local_K_data, local_b_data, local_Jac_data, - indices); + cache.add(local_b_data, local_Jac_data, indices); } int getNumberOfThreads() @@ -121,7 +117,7 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( std::vector const& dof_tables, const double t, double const dt, std::vector const& xs, std::vector const& x_prevs, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { // checks ////////////////////////////////////////////////////////////////// if (process_id != 0) @@ -150,8 +146,6 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( // algorithm /////////////////////////////////////////////////////////////// auto stats = CumulativeStats::create(); - ConcurrentMatrixView M_view(M); - ConcurrentMatrixView K_view(K); ConcurrentMatrixView b_view(b); ConcurrentMatrixView Jac_view(Jac); @@ -168,8 +162,6 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( // temporary data only stored here in order to avoid frequent memory // reallocations. - std::vector local_M_data; - std::vector local_K_data; std::vector local_b_data; std::vector local_Jac_data; std::vector indices; @@ -178,7 +170,7 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( auto const jac_asm = jacobian_assembler_.copy(); auto stats_this_thread = stats->clone(); - MultiMatrixElementCache cache{M_view, K_view, b_view, Jac_view, + MultiMatrixElementCache cache{b_view, Jac_view, stats_this_thread->data}; // TODO corner case: what if all elements on a submesh are deactivated? @@ -204,8 +196,7 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( { assembleWithJacobianOneElement( element_id, loc_asm, dof_table, t, dt, x, x_prev, - local_M_data, local_K_data, local_b_data, - local_Jac_data, indices, *jac_asm, cache); + local_b_data, local_Jac_data, indices, *jac_asm, cache); } catch (...) { @@ -214,9 +205,8 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( continue; } - local_matrix_output_(t, process_id, element_id, local_M_data, - local_K_data, local_b_data, - &local_Jac_data); + local_matrix_output_(t, process_id, element_id, nullptr, + nullptr, local_b_data, &local_Jac_data); } } else @@ -242,8 +232,7 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( { assembleWithJacobianOneElement( element_id, loc_asm, dof_table, t, dt, x, x_prev, - local_M_data, local_K_data, local_b_data, - local_Jac_data, indices, *jac_asm, cache); + local_b_data, local_Jac_data, indices, *jac_asm, cache); } catch (...) { @@ -252,16 +241,15 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( continue; } - local_matrix_output_(t, process_id, element_id, local_M_data, - local_K_data, local_b_data, - &local_Jac_data); + local_matrix_output_(t, process_id, element_id, nullptr, + nullptr, local_b_data, &local_Jac_data); } } } // OpenMP parallel section stats->print(); - global_matrix_output_(t, process_id, M, K, b, &Jac); + global_matrix_output_(t, process_id, nullptr, nullptr, b, &Jac); exception.rethrow(); } } // namespace ProcessLib::Assembly diff --git a/ProcessLib/Assembly/ParallelVectorMatrixAssembler.h b/ProcessLib/Assembly/ParallelVectorMatrixAssembler.h index c00ae80eebc..c7601f3be85 100644 --- a/ProcessLib/Assembly/ParallelVectorMatrixAssembler.h +++ b/ProcessLib/Assembly/ParallelVectorMatrixAssembler.h @@ -30,7 +30,7 @@ class ParallelVectorMatrixAssembler std::vector const& dof_tables, const double t, double const dt, std::vector const& xs, std::vector const& x_prevs, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac); + GlobalVector& b, GlobalMatrix& Jac); private: AbstractJacobianAssembler& jacobian_assembler_; diff --git a/ProcessLib/AssemblyMixin.h b/ProcessLib/AssemblyMixin.h index b082a5873ff..72e4f1a397f 100644 --- a/ProcessLib/AssemblyMixin.h +++ b/ProcessLib/AssemblyMixin.h @@ -155,8 +155,7 @@ class AssemblyMixin : private AssemblyMixinBase void assembleWithJacobian(double const t, double const dt, std::vector const& x, std::vector const& x_prev, - int const process_id, GlobalMatrix& M, - GlobalMatrix& K, GlobalVector& b, + int const process_id, GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssemblyMixin assembleWithJacobian(t={}, dt={}, process_id={}).", @@ -179,7 +178,7 @@ class AssemblyMixin : private AssemblyMixinBase pvma_.assembleWithJacobian(loc_asms, sad.active_element_ids, dof_tables, t, dt, x, x_prev, - process_id, M, K, b_submesh, Jac); + process_id, b_submesh, Jac); MathLib::LinAlg::axpy(b, 1.0, b_submesh); @@ -198,7 +197,7 @@ class AssemblyMixin : private AssemblyMixinBase pvma_.assembleWithJacobian(loc_asms, pv.getActiveElementIDs(), dof_tables, t, dt, x, x_prev, process_id, - M, K, b, Jac); + b, Jac); } AssemblyMixinBase::copyResiduumVectorsToBulkMesh( diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryCondition.h b/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryCondition.h index b5ddae8ff6f..f4fce0481d2 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryCondition.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryCondition.h @@ -39,7 +39,7 @@ class BoundaryCondition //! \c K and the vector \c b. virtual void applyNaturalBC(const double /*t*/, std::vector const& /*x*/, - int const /*process_id*/, GlobalMatrix& /*K*/, + int const /*process_id*/, GlobalMatrix* /*K*/, GlobalVector& /*b*/, GlobalMatrix* /*Jac*/) { // By default it is assumed that the BC is not a natural BC. Therefore diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryConditionCollection.cpp b/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryConditionCollection.cpp index 73808f20cbd..78b4d3e3ce4 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryConditionCollection.cpp +++ b/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryConditionCollection.cpp @@ -14,7 +14,7 @@ namespace ProcessLib { void BoundaryConditionCollection::applyNaturalBC( const double t, std::vector const& x, int const process_id, - GlobalMatrix& K, GlobalVector& b, GlobalMatrix* Jac) const + GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) const { for (auto const& bc : _boundary_conditions) { diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryConditionCollection.h b/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryConditionCollection.h index 1d6eddca2fd..e955f045b33 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryConditionCollection.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/BoundaryConditionCollection.h @@ -28,7 +28,7 @@ class BoundaryConditionCollection final } void applyNaturalBC(const double t, std::vector const& x, - int const process_id, GlobalMatrix& K, GlobalVector& b, + int const process_id, GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) const; std::vector> const* diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryCondition-impl.h b/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryCondition-impl.h index d10bbc98412..74b3c8856c5 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryCondition-impl.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryCondition-impl.h @@ -82,7 +82,7 @@ void GenericNaturalBoundaryCondition const& x, int const process_id, - GlobalMatrix& K, + GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) { diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryCondition.h b/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryCondition.h index 6cc2d85de29..a3bb00c54cf 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryCondition.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryCondition.h @@ -36,7 +36,7 @@ class GenericNaturalBoundaryCondition final : public BoundaryCondition /// Calls local assemblers which calculate their contributions to the global /// matrix and the right-hand-side. void applyNaturalBC(const double t, std::vector const& x, - int const process_id, GlobalMatrix& K, GlobalVector& b, + int const process_id, GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) override; private: diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryConditionLocalAssembler.h index 73fc9ca0ff0..d9d64f4bdcb 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryConditionLocalAssembler.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/GenericNaturalBoundaryConditionLocalAssembler.h @@ -27,7 +27,7 @@ class GenericNaturalBoundaryConditionLocalAssemblerInterface std::size_t const id, NumLib::LocalToGlobalIndexMap const& dof_table_boundary, double const t, std::vector const& x, int const process_id, - GlobalMatrix& K, GlobalVector& b, GlobalMatrix* Jac) = 0; + GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) = 0; }; template diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/HCNonAdvectiveFreeComponentFlowBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryConditionAndSourceTerm/HCNonAdvectiveFreeComponentFlowBoundaryConditionLocalAssembler.h index 72ae581b919..2c9add1f9eb 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/HCNonAdvectiveFreeComponentFlowBoundaryConditionLocalAssembler.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/HCNonAdvectiveFreeComponentFlowBoundaryConditionLocalAssembler.h @@ -57,7 +57,7 @@ class HCNonAdvectiveFreeComponentFlowBoundaryConditionLocalAssembler final void assemble(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const& dof_table_boundary, double const t, std::vector const& x, - int const process_id, GlobalMatrix& /*K*/, GlobalVector& b, + int const process_id, GlobalMatrix* /*K*/, GlobalVector& b, GlobalMatrix* /*Jac*/) override { NodalVectorType _local_rhs = NodalVectorType::Zero(_local_matrix_size); diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/NeumannBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryConditionAndSourceTerm/NeumannBoundaryConditionLocalAssembler.h index 30b7628f92a..0b4c693083e 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/NeumannBoundaryConditionLocalAssembler.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/NeumannBoundaryConditionLocalAssembler.h @@ -55,7 +55,7 @@ class NeumannBoundaryConditionLocalAssembler final void assemble(std::size_t const id, NumLib::LocalToGlobalIndexMap const& dof_table_boundary, double const t, std::vector const& /*x*/, - int const /*process_id*/, GlobalMatrix& /*K*/, + int const /*process_id*/, GlobalMatrix* /*K*/, GlobalVector& b, GlobalMatrix* /*Jac*/) override { _local_rhs.setZero(); diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryCondition-impl.h b/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryCondition-impl.h index b82e164f17d..39110b2b7e7 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryCondition-impl.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryCondition-impl.h @@ -62,7 +62,7 @@ template class LocalAssemblerImplementation> void NormalTractionBoundaryCondition:: applyNaturalBC(const double t, std::vector const& x, - int const /*process_id*/, GlobalMatrix& K, GlobalVector& b, + int const /*process_id*/, GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) { GlobalExecutor::executeMemberOnDereferenced( diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryCondition.h b/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryCondition.h index 21d3450b779..e5c646c07d0 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryCondition.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryCondition.h @@ -44,7 +44,7 @@ class NormalTractionBoundaryCondition final : public BoundaryCondition /// Calls local assemblers which calculate their contributions to the global /// matrix and the right-hand-side. void applyNaturalBC(const double t, std::vector const& x, - int const process_id, GlobalMatrix& K, GlobalVector& b, + int const process_id, GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) override; private: diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryConditionLocalAssembler.h index 578e06e4537..8cf7db87f67 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryConditionLocalAssembler.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/NormalTractionBoundaryConditionLocalAssembler.h @@ -42,7 +42,7 @@ class NormalTractionBoundaryConditionLocalAssemblerInterface virtual void assemble( std::size_t const id, NumLib::LocalToGlobalIndexMap const& dof_table_boundary, double const t, - std::vector const& /*x*/, GlobalMatrix& /*K*/, + std::vector const& /*x*/, GlobalMatrix* /*K*/, GlobalVector& b, GlobalMatrix* /*Jac*/) = 0; virtual ~NormalTractionBoundaryConditionLocalAssemblerInterface() = default; }; @@ -115,7 +115,7 @@ class NormalTractionBoundaryConditionLocalAssembler final void assemble(std::size_t const id, NumLib::LocalToGlobalIndexMap const& dof_table_boundary, double const t, std::vector const& /*x*/, - GlobalMatrix& /*K*/, GlobalVector& local_rhs, + GlobalMatrix* /*K*/, GlobalVector& local_rhs, GlobalMatrix* /*Jac*/) override { _local_rhs.setZero(); diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryCondition.cpp b/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryCondition.cpp index 0bb86a848ef..00059cd9bf6 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryCondition.cpp +++ b/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryCondition.cpp @@ -245,7 +245,7 @@ double PythonBoundaryCondition::interpolateToHigherOrderNode( void PythonBoundaryCondition::applyNaturalBC( const double t, std::vector const& x, int const process_id, - GlobalMatrix& K, GlobalVector& b, GlobalMatrix* Jac) + GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) { FlushStdoutGuard guard(_flush_stdout); diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryCondition.h b/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryCondition.h index 0f553108a78..d664e2b1c7f 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryCondition.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryCondition.h @@ -69,7 +69,7 @@ class PythonBoundaryCondition final : public BoundaryCondition NumLib::IndexValueVector& bc_values) const override; void applyNaturalBC(const double t, std::vector const& x, - int const process_id, GlobalMatrix& K, GlobalVector& b, + int const process_id, GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) override; private: diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryConditionLocalAssembler.h index 9993da0a62f..e08c5031e14 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryConditionLocalAssembler.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/Python/PythonBoundaryConditionLocalAssembler.h @@ -40,7 +40,7 @@ class PythonBoundaryConditionLocalAssembler final void assemble(std::size_t const boundary_element_id, NumLib::LocalToGlobalIndexMap const& dof_table_boundary, double const t, std::vector const& xs, - int const process_id, GlobalMatrix& /*K*/, GlobalVector& b, + int const process_id, GlobalMatrix* /*K*/, GlobalVector& b, GlobalMatrix* const Jac) override { impl_.assemble(boundary_element_id, dof_table_boundary, t, diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/RobinBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryConditionAndSourceTerm/RobinBoundaryConditionLocalAssembler.h index 045656f2129..4bc2e0e4be2 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/RobinBoundaryConditionLocalAssembler.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/RobinBoundaryConditionLocalAssembler.h @@ -50,7 +50,7 @@ class RobinBoundaryConditionLocalAssembler final void assemble(std::size_t const id, NumLib::LocalToGlobalIndexMap const& dof_table_boundary, double const t, std::vector const& xs, - int const process_id, GlobalMatrix& K, GlobalVector& b, + int const process_id, GlobalMatrix* K, GlobalVector& b, GlobalMatrix* Jac) override { _local_K.setZero(); @@ -116,11 +116,17 @@ class RobinBoundaryConditionLocalAssembler final indices), _local_K); } + else if (K != nullptr) + { + K->add(NumLib::LocalToGlobalIndexMap::RowColumnIndices(indices, + indices), + _local_K); + } else { - K.add(NumLib::LocalToGlobalIndexMap::RowColumnIndices(indices, - indices), - _local_K); + OGS_FATAL( + "In Robin boundary condition assembler, both, the Jacobian and " + "the K-matrices are null, but one matrix must be provided."); } } diff --git a/ProcessLib/BoundaryConditionAndSourceTerm/VariableDependentNeumannBoundaryConditionLocalAssembler.h b/ProcessLib/BoundaryConditionAndSourceTerm/VariableDependentNeumannBoundaryConditionLocalAssembler.h index 2c9ee4ed74f..b6e23676970 100644 --- a/ProcessLib/BoundaryConditionAndSourceTerm/VariableDependentNeumannBoundaryConditionLocalAssembler.h +++ b/ProcessLib/BoundaryConditionAndSourceTerm/VariableDependentNeumannBoundaryConditionLocalAssembler.h @@ -57,7 +57,7 @@ class VariableDependentNeumannBoundaryConditionLocalAssembler final void assemble(std::size_t const mesh_item_id, NumLib::LocalToGlobalIndexMap const& dof_table_boundary, double const t, std::vector const& x, - int const process_id, GlobalMatrix& /*K*/, GlobalVector& b, + int const process_id, GlobalMatrix* /*K*/, GlobalVector& b, GlobalMatrix* /*Jac*/) override { NodalVectorType _local_rhs(_local_matrix_size); diff --git a/ProcessLib/CentralDifferencesJacobianAssembler.cpp b/ProcessLib/CentralDifferencesJacobianAssembler.cpp index 8c49297c015..f6d1bb8cc8f 100644 --- a/ProcessLib/CentralDifferencesJacobianAssembler.cpp +++ b/ProcessLib/CentralDifferencesJacobianAssembler.cpp @@ -31,9 +31,11 @@ void CentralDifferencesJacobianAssembler::assembleWithJacobian( LocalAssemblerInterface& local_assembler, const double t, double const dt, const std::vector& local_x_data, const std::vector& local_x_prev_data, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) { + std::vector local_M_data(local_Jac_data.size()); + std::vector local_K_data(local_Jac_data.size()); + // TODO do not check in every call. if (local_x_data.size() % _absolute_epsilons.size() != 0) { diff --git a/ProcessLib/CentralDifferencesJacobianAssembler.h b/ProcessLib/CentralDifferencesJacobianAssembler.h index 78c020ffdcd..78691abd416 100644 --- a/ProcessLib/CentralDifferencesJacobianAssembler.h +++ b/ProcessLib/CentralDifferencesJacobianAssembler.h @@ -55,8 +55,6 @@ class CentralDifferencesJacobianAssembler final double const t, double const dt, std::vector const& local_x_data, std::vector const& local_x_prev_data, - std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/CompareJacobiansJacobianAssembler.cpp b/ProcessLib/CompareJacobiansJacobianAssembler.cpp index e06e472a49d..af682af6631 100644 --- a/ProcessLib/CompareJacobiansJacobianAssembler.cpp +++ b/ProcessLib/CompareJacobiansJacobianAssembler.cpp @@ -127,46 +127,26 @@ namespace ProcessLib void CompareJacobiansJacobianAssembler::assembleWithJacobian( LocalAssemblerInterface& local_assembler, double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) { ++_counter; auto const num_dof = local_x.size(); - auto to_mat = [num_dof](std::vector const& data) - { - if (data.empty()) - { - return Eigen::Map const>( - nullptr, 0, 0); - } - - return MathLib::toMatrix(data, num_dof, num_dof); - }; // First assembly -- the one whose results will be added to the global // equation system finally. _asm1->assembleWithJacobian(local_assembler, t, dt, local_x, local_x_prev, - local_M_data, local_K_data, local_b_data, - local_Jac_data); + local_b_data, local_Jac_data); - auto const local_M1 = to_mat(local_M_data); - auto const local_K1 = to_mat(local_K_data); auto const local_b1 = MathLib::toVector(local_b_data); - std::vector local_M_data2; - std::vector local_K_data2; std::vector local_b_data2; std::vector local_Jac_data2; // Second assembly -- used for checking only. _asm2->assembleWithJacobian(local_assembler, t, dt, local_x, local_x_prev, - local_M_data2, local_K_data2, local_b_data2, - local_Jac_data2); + local_b_data2, local_Jac_data2); - auto const local_M2 = to_mat(local_M_data2); - auto const local_K2 = to_mat(local_K_data2); auto const local_b2 = MathLib::toVector(local_b_data2); auto const local_Jac1 = MathLib::toMatrix(local_Jac_data, num_dof, num_dof); @@ -242,35 +222,17 @@ void CompareJacobiansJacobianAssembler::assembleWithJacobian( } }; - check_equality(local_M1, local_M2); - check_equality(local_K1, local_K2); check_equality(local_b1, local_b2); Eigen::VectorXd res1 = Eigen::VectorXd::Zero(num_dof); auto const x = MathLib::toVector(local_x); auto const x_dot = ((x - MathLib::toVector(local_x_prev)) / dt).eval(); - if (local_M1.size() != 0) - { - res1.noalias() += local_M1 * x_dot; - } - if (local_K1.size() != 0) - { - res1.noalias() += local_K1 * x; - } if (local_b1.size() != 0) { res1.noalias() -= local_b1; } Eigen::VectorXd res2 = Eigen::VectorXd::Zero(num_dof); - if (local_M2.size() != 0) - { - res2.noalias() += local_M2 * x_dot; - } - if (local_K2.size() != 0) - { - res2.noalias() += local_K2 * x; - } if (local_b2.size() != 0) { res2.noalias() -= local_b2; @@ -302,8 +264,8 @@ void CompareJacobiansJacobianAssembler::assembleWithJacobian( "might be\n" << "# (a) that the assembly routine has side " "effects or\n" - << "# (b) that the assembly routines for M, K " - "and b themselves differ.\n" + << "# (b) that the assembly routines for b " + "themselves differ.\n" << "#######################################################\n" << '\n'; } @@ -347,22 +309,6 @@ void CompareJacobiansJacobianAssembler::assembleWithJacobian( _log_file << '\n'; - dump_py(_log_file, "M_1", local_M1); - dump_py(_log_file, "M_2", local_M2); - if (fatal_error && local_M1.size() == local_M2.size()) - { - dump_py(_log_file, "delta_M", local_M2 - local_M1); - _log_file << '\n'; - } - - dump_py(_log_file, "K_1", local_K1); - dump_py(_log_file, "K_2", local_K2); - if (fatal_error && local_K1.size() == local_K2.size()) - { - dump_py(_log_file, "delta_K", local_K2 - local_K1); - _log_file << '\n'; - } - dump_py(_log_file, "b_1", local_b_data); dump_py(_log_file, "b_2", local_b_data2); if (fatal_error && local_b1.size() == local_b2.size()) diff --git a/ProcessLib/CompareJacobiansJacobianAssembler.h b/ProcessLib/CompareJacobiansJacobianAssembler.h index 77039ada136..8b7fcd31328 100644 --- a/ProcessLib/CompareJacobiansJacobianAssembler.h +++ b/ProcessLib/CompareJacobiansJacobianAssembler.h @@ -56,8 +56,6 @@ class CompareJacobiansJacobianAssembler final : public AbstractJacobianAssembler double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/ForwardDifferencesJacobianAssembler.cpp b/ProcessLib/ForwardDifferencesJacobianAssembler.cpp index 2419fba457f..53fd7298bb8 100644 --- a/ProcessLib/ForwardDifferencesJacobianAssembler.cpp +++ b/ProcessLib/ForwardDifferencesJacobianAssembler.cpp @@ -29,9 +29,11 @@ void ForwardDifferencesJacobianAssembler::assembleWithJacobian( LocalAssemblerInterface& local_assembler, const double t, double const dt, const std::vector& local_x_data, const std::vector& local_x_prev_data, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) { + std::vector local_M_data(local_Jac_data.size()); + std::vector local_K_data(local_Jac_data.size()); + // TODO do not check in every call. if (local_x_data.size() % _absolute_epsilons.size() != 0) { diff --git a/ProcessLib/ForwardDifferencesJacobianAssembler.h b/ProcessLib/ForwardDifferencesJacobianAssembler.h index 60b9688e552..78254cdf7ab 100644 --- a/ProcessLib/ForwardDifferencesJacobianAssembler.h +++ b/ProcessLib/ForwardDifferencesJacobianAssembler.h @@ -47,8 +47,6 @@ class ForwardDifferencesJacobianAssembler final double const t, double const dt, std::vector const& local_x_data, std::vector const& local_x_prev_data, - std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/LocalAssemblerInterface.cpp b/ProcessLib/LocalAssemblerInterface.cpp index fdd5d1a9107..202d306c1fd 100644 --- a/ProcessLib/LocalAssemblerInterface.cpp +++ b/ProcessLib/LocalAssemblerInterface.cpp @@ -46,8 +46,6 @@ void LocalAssemblerInterface::assembleWithJacobian( double const /*t*/, double const /*dt*/, std::vector const& /*local_x*/, std::vector const& /*local_x_prev*/, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& /*local_b_data*/, std::vector& /*local_Jac_data*/) { @@ -59,8 +57,6 @@ void LocalAssemblerInterface::assembleWithJacobian( void LocalAssemblerInterface::assembleWithJacobianForStaggeredScheme( double const /*t*/, double const /*dt*/, Eigen::VectorXd const& /*local_x*/, Eigen::VectorXd const& /*local_x_prev*/, int const /*process_id*/, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& /*local_b_data*/, std::vector& /*local_Jac_data*/) { diff --git a/ProcessLib/LocalAssemblerInterface.h b/ProcessLib/LocalAssemblerInterface.h index ba852d125fd..7832b323790 100644 --- a/ProcessLib/LocalAssemblerInterface.h +++ b/ProcessLib/LocalAssemblerInterface.h @@ -62,15 +62,12 @@ class LocalAssemblerInterface virtual void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data); virtual void assembleWithJacobianForStaggeredScheme( double const t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data); virtual void computeSecondaryVariable( diff --git a/ProcessLib/Process.cpp b/ProcessLib/Process.cpp index 9e9277e0169..12afbb610ad 100644 --- a/ProcessLib/Process.cpp +++ b/ProcessLib/Process.cpp @@ -276,7 +276,7 @@ void Process::assemble(const double t, double const dt, assembleConcreteProcess(t, dt, x, x_prev, process_id, M, K, b); // the last argument is for the jacobian, nullptr is for a unused jacobian - _boundary_conditions[process_id].applyNaturalBC(t, x, process_id, K, b, + _boundary_conditions[process_id].applyNaturalBC(t, x, process_id, &K, b, nullptr); // the last argument is for the jacobian, nullptr is for a unused jacobian @@ -287,20 +287,18 @@ void Process::assemble(const double t, double const dt, void Process::assembleWithJacobian(const double t, double const dt, std::vector const& x, std::vector const& x_prev, - int const process_id, GlobalMatrix& M, - GlobalMatrix& K, GlobalVector& b, + int const process_id, GlobalVector& b, GlobalMatrix& Jac) { assert(x.size() == x_prev.size()); setLocalAccessibleVectors(x, x_prev); - assembleWithJacobianConcreteProcess(t, dt, x, x_prev, process_id, M, K, b, - Jac); + assembleWithJacobianConcreteProcess(t, dt, x, x_prev, process_id, b, Jac); // TODO: apply BCs to Jacobian. - _boundary_conditions[process_id].applyNaturalBC(t, x, process_id, K, b, - &Jac); + _boundary_conditions[process_id].applyNaturalBC(t, x, process_id, nullptr, + b, &Jac); // the last argument is for the jacobian, nullptr is for a unused jacobian _source_term_collections[process_id].integrate(t, *x[process_id], b, &Jac); diff --git a/ProcessLib/Process.h b/ProcessLib/Process.h index be39d2bcb13..3b2f076105d 100644 --- a/ProcessLib/Process.h +++ b/ProcessLib/Process.h @@ -123,8 +123,7 @@ class Process void assembleWithJacobian(const double t, double const dt, std::vector const& x, std::vector const& x_prev, - int const process_id, GlobalMatrix& M, - GlobalMatrix& K, GlobalVector& b, + int const process_id, GlobalVector& b, GlobalMatrix& Jac) final; /* Computes data necessary for output. @@ -251,8 +250,7 @@ class Process virtual void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) = 0; + GlobalVector& b, GlobalMatrix& Jac) = 0; virtual void preTimestepConcreteProcess( std::vector const& /*x*/, diff --git a/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.cpp b/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.cpp index 1fcf8b7de10..a72a0c2c380 100644 --- a/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.cpp +++ b/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.cpp @@ -74,7 +74,7 @@ void SteadyStateDiffusion::assembleConcreteProcess( void SteadyStateDiffusion::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian SteadyStateDiffusion."); @@ -84,7 +84,7 @@ void SteadyStateDiffusion::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } } // namespace SteadyStateDiffusion diff --git a/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.h b/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.h index a1ee63186f7..2320aca42af 100644 --- a/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.h +++ b/ProcessLib/SteadyStateDiffusion/SteadyStateDiffusion.h @@ -96,8 +96,7 @@ class SteadyStateDiffusion final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; SteadyStateDiffusionData _process_data; diff --git a/ProcessLib/VectorMatrixAssembler.cpp b/ProcessLib/VectorMatrixAssembler.cpp index 422e900926f..35e00a85b33 100644 --- a/ProcessLib/VectorMatrixAssembler.cpp +++ b/ProcessLib/VectorMatrixAssembler.cpp @@ -99,7 +99,7 @@ void VectorMatrixAssembler::assemble( b->add(indices, _local_b_data); } - _local_output(t, process_id, mesh_item_id, _local_M_data, _local_K_data, + _local_output(t, process_id, mesh_item_id, &_local_M_data, &_local_K_data, _local_b_data); } @@ -108,7 +108,7 @@ void VectorMatrixAssembler::assembleWithJacobian( std::vector const& dof_tables, const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix* M, GlobalMatrix* K, GlobalVector* b, GlobalMatrix* Jac) + GlobalVector* b, GlobalMatrix* Jac) { std::vector> indices_of_processes; indices_of_processes.reserve(dof_tables.size()); @@ -119,8 +119,6 @@ void VectorMatrixAssembler::assembleWithJacobian( auto const& indices = indices_of_processes[process_id]; - _local_M_data.clear(); - _local_K_data.clear(); _local_b_data.clear(); _local_Jac_data.clear(); @@ -131,8 +129,8 @@ void VectorMatrixAssembler::assembleWithJacobian( auto const local_x = x[process_id]->get(indices); auto const local_x_prev = x_prev[process_id]->get(indices); _jacobian_assembler.assembleWithJacobian( - local_assembler, t, dt, local_x, local_x_prev, _local_M_data, - _local_K_data, _local_b_data, _local_Jac_data); + local_assembler, t, dt, local_x, local_x_prev, _local_b_data, + _local_Jac_data); } else // Staggered scheme { @@ -146,23 +144,13 @@ void VectorMatrixAssembler::assembleWithJacobian( _jacobian_assembler.assembleWithJacobianForStaggeredScheme( local_assembler, t, dt, local_x, local_x_prev, process_id, - _local_M_data, _local_K_data, _local_b_data, _local_Jac_data); + _local_b_data, _local_Jac_data); } auto const num_r_c = indices.size(); auto const r_c_indices = NumLib::LocalToGlobalIndexMap::RowColumnIndices(indices, indices); - if (M && !_local_M_data.empty()) - { - auto const local_M = MathLib::toMatrix(_local_M_data, num_r_c, num_r_c); - M->add(r_c_indices, local_M); - } - if (K && !_local_K_data.empty()) - { - auto const local_K = MathLib::toMatrix(_local_K_data, num_r_c, num_r_c); - K->add(r_c_indices, local_K); - } if (b && !_local_b_data.empty()) { assert(_local_b_data.size() == num_r_c); @@ -181,8 +169,8 @@ void VectorMatrixAssembler::assembleWithJacobian( "errors in the local assembler of the current process."); } - _local_output(t, process_id, mesh_item_id, _local_M_data, _local_K_data, - _local_b_data, &_local_Jac_data); + _local_output(t, process_id, mesh_item_id, nullptr, nullptr, _local_b_data, + &_local_Jac_data); } } // namespace ProcessLib diff --git a/ProcessLib/VectorMatrixAssembler.h b/ProcessLib/VectorMatrixAssembler.h index 2bd58de8a05..75c1c04b539 100644 --- a/ProcessLib/VectorMatrixAssembler.h +++ b/ProcessLib/VectorMatrixAssembler.h @@ -50,7 +50,7 @@ class VectorMatrixAssembler final std::vector const& x_prev, int const process_id, GlobalMatrix* M, GlobalMatrix* K, GlobalVector* b); - //! Assembles \c M, \c K, \c b, and the Jacobian \c Jac of the residual. + //! Assembles \c b, and the Jacobian \c Jac of the residual. //! \note The Jacobian must be assembled. void assembleWithJacobian( std::size_t const mesh_item_id, @@ -58,7 +58,7 @@ class VectorMatrixAssembler final std::vector const& dof_tables, const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix* M, GlobalMatrix* K, GlobalVector* b, GlobalMatrix* Jac); + GlobalVector* b, GlobalMatrix* Jac); private: // temporary data only stored here in order to avoid frequent memory diff --git a/Tests/NumLib/ODEs.h b/Tests/NumLib/ODEs.h index 5e4c7a4ce79..2453a4a66a7 100644 --- a/Tests/NumLib/ODEs.h +++ b/Tests/NumLib/ODEs.h @@ -19,6 +19,43 @@ template class ODETraits; +inline GlobalVector computeResidumFromMKb(double const dt, + GlobalVector const& x_curr, + GlobalVector const& x_prev, + GlobalMatrix const& M, + GlobalMatrix const& K, + GlobalVector const& b) +{ + using namespace MathLib::LinAlg; + GlobalVector res(b); + + // res = -M * x_dot - K * x_curr + b + GlobalVector x_dot; + copy(x_curr, x_dot); // x_dot = x + axpy(x_dot, -1., x_prev); // x_dot = x - x_prev + scale(x_dot, 1. / dt); // x_dot = (x - x_prev)/dt + matMult(M, x_dot, res); // res = M*x_dot + matMultAdd(K, x_curr, res, res); // res = M*x_dot + K*x + axpy(res, -1., b); // res = M*x_dot + K*x - b + scale(res, -1.); // res = -M*x_dot - K*x + b + return res; +} + +inline GlobalMatrix computeJacobianFromMK(double const dt, + GlobalMatrix const& M, + GlobalMatrix const& K) +{ + using namespace MathLib::LinAlg; + GlobalMatrix J(M); + // compute J = M*1/dt + K + copy(M, J); + scale(J, 1. / dt); + axpy(J, 1., K); + finalizeAssembly(J); + + return J; +} + // ODE 1 ////////////////////////////////////////////////////////// class ODE1 final : public NumLib::ODESystem< NumLib::ODESystemTag::FirstOrderImplicitQuasilinear, @@ -48,22 +85,23 @@ class ODE1 final : public NumLib::ODESystem< } void assembleWithJacobian(const double /*t*/, double const dt, - std::vector const& /*x_curr*/, - std::vector const& /*xdot*/, - int const /*process_id*/, GlobalMatrix& M, - GlobalMatrix& K, GlobalVector& b, + std::vector const& x_curr, + std::vector const& x_prev, + int const process_id, GlobalVector& b, GlobalMatrix& Jac) override { - namespace LinAlg = MathLib::LinAlg; - + MathLib::LinAlg::setLocalAccessibleVector(*x_curr[process_id]); + GlobalMatrix M(Jac); + GlobalMatrix K(Jac); setMKbValues(M, K, b); - // compute Jac = M*1/dt + K - LinAlg::finalizeAssembly(M); - LinAlg::copy(M, Jac); - LinAlg::scale(Jac, 1. / dt); - LinAlg::finalizeAssembly(K); - LinAlg::axpy(Jac, 1., K); + MathLib::LinAlg::finalizeAssembly(M); + MathLib::LinAlg::finalizeAssembly(K); + MathLib::LinAlg::copy(computeJacobianFromMK(dt, M, K), Jac); + MathLib::LinAlg::copy( + computeResidumFromMKb(dt, *x_curr[process_id], *x_prev[process_id], + M, K, b), + b); } MathLib::MatrixSpecifications getMatrixSpecifications( @@ -138,27 +176,23 @@ class ODE2 final : public NumLib::ODESystem< } void assembleWithJacobian(const double /*t*/, double const dt, - std::vector const& x, - std::vector const& /*xdot*/, - int const process_id, GlobalMatrix& M, - GlobalMatrix& K, GlobalVector& b, + std::vector const& x_curr, + std::vector const& x_prev, + int const process_id, GlobalVector& b, GlobalMatrix& Jac) override { - MathLib::LinAlg::setLocalAccessibleVector(*x[process_id]); - setMKbValues(*x[process_id], M, K, b); - - namespace LinAlg = MathLib::LinAlg; - - LinAlg::finalizeAssembly(M); - // compute Jac = M*1/dt + dK_dx + K - LinAlg::copy(M, Jac); - LinAlg::scale(Jac, 1. / dt); - - MathLib::addToMatrix(Jac, {(*x[process_id])[0]}); // add dK_dx + MathLib::LinAlg::setLocalAccessibleVector(*x_curr[process_id]); + GlobalMatrix M(Jac); + GlobalMatrix K(Jac); + setMKbValues(*x_curr[process_id], M, K, b); - LinAlg::finalizeAssembly(K); - LinAlg::finalizeAssembly(Jac); - LinAlg::axpy(Jac, 1., K); + MathLib::LinAlg::finalizeAssembly(M); + MathLib::LinAlg::finalizeAssembly(K); + MathLib::LinAlg::copy(computeJacobianFromMK(dt, M, K), Jac); + MathLib::LinAlg::copy( + computeResidumFromMKb(dt, *x_curr[process_id], *x_prev[process_id], + M, K, b), + b); } MathLib::MatrixSpecifications getMatrixSpecifications( @@ -242,21 +276,22 @@ class ODE3 final : public NumLib::ODESystem< void assembleWithJacobian(const double /*t*/, double const dt, std::vector const& x_curr, - std::vector const& xdot, - int const process_id, GlobalMatrix& M, - GlobalMatrix& K, GlobalVector& b, + std::vector const& x_prev, + int const process_id, GlobalVector& b, GlobalMatrix& Jac) override { MathLib::LinAlg::setLocalAccessibleVector(*x_curr[process_id]); - MathLib::LinAlg::setLocalAccessibleVector(*xdot[process_id]); + MathLib::LinAlg::setLocalAccessibleVector(*x_prev[process_id]); + GlobalMatrix M(Jac); + GlobalMatrix K(Jac); setMKbValues(*x_curr[process_id], M, K, b); auto const u = (*x_curr[process_id])[0]; auto const v = (*x_curr[process_id])[1]; - auto const du = (*xdot[process_id])[0]; - auto const dv = (*xdot[process_id])[1]; + auto const du = (*x_prev[process_id])[0]; + auto const dv = (*x_prev[process_id])[1]; namespace LinAlg = MathLib::LinAlg; @@ -264,16 +299,15 @@ class ODE3 final : public NumLib::ODESystem< // with dxdot/dx = 1/dt and dx/dx = 1 LinAlg::finalizeAssembly(M); - LinAlg::copy(M, Jac); - LinAlg::scale(Jac, 1. / dt); // Jac = M * 1/dt + LinAlg::finalizeAssembly(K); + LinAlg::copy(computeJacobianFromMK(dt, M, K), Jac); + LinAlg::copy(computeResidumFromMKb(dt, *x_curr[process_id], + *x_prev[process_id], M, K, b), + b); // add dM/dx \cdot \dot x MathLib::addToMatrix(Jac, {du - dv, 0.0, 0.0, dv - du}); - LinAlg::finalizeAssembly(K); - LinAlg::finalizeAssembly(Jac); - LinAlg::axpy(Jac, 1., K); // add K - // add dK/dx \cdot \dot x MathLib::addToMatrix(Jac, {-du - dv, -du, 0.0, du + 2.0 * dv}); diff --git a/Tests/ProcessLib/TestJacobianAssembler.cpp b/Tests/ProcessLib/TestJacobianAssembler.cpp index 2e764e99b5e..2fcb24a4c01 100644 --- a/Tests/ProcessLib/TestJacobianAssembler.cpp +++ b/Tests/ProcessLib/TestJacobianAssembler.cpp @@ -331,11 +331,12 @@ class LocalAssemblerB final : public ProcessLib::LocalAssemblerInterface void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_xdot, - std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override { + std::vector local_M_data(local_Jac_data.size()); + std::vector local_K_data(local_Jac_data.size()); + assemble(t, dt, local_x, local_xdot, local_M_data, local_K_data, local_b_data); @@ -383,41 +384,17 @@ struct ProcessLibForwardDifferencesJacobianAssembler : public ::testing::Test double const eps = std::numeric_limits::epsilon(); - std::vector M_data_cd; - std::vector K_data_cd; std::vector b_data_cd; std::vector Jac_data_cd; - std::vector M_data_ana; - std::vector K_data_ana; std::vector b_data_ana; std::vector Jac_data_ana; double const t = 0.0; - jac_asm_cd.assembleWithJacobian(loc_asm, t, dt, x, xdot, M_data_cd, - K_data_cd, b_data_cd, Jac_data_cd); - - jac_asm_ana.assembleWithJacobian(loc_asm, t, dt, x, xdot, M_data_ana, - K_data_ana, b_data_ana, Jac_data_ana); + jac_asm_cd.assembleWithJacobian(loc_asm, t, dt, x, xdot, b_data_cd, + Jac_data_cd); - if (LocAsm::asmM) - { - ASSERT_EQ(x.size() * x.size(), M_data_cd.size()); - ASSERT_EQ(x.size() * x.size(), M_data_ana.size()); - for (std::size_t i = 0; i < x.size() * x.size(); ++i) - { - EXPECT_NEAR(M_data_ana[i], M_data_cd[i], eps); - } - } - - if (LocAsm::asmK) - { - ASSERT_EQ(x.size() * x.size(), K_data_cd.size()); - ASSERT_EQ(x.size() * x.size(), K_data_ana.size()); - for (std::size_t i = 0; i < x.size() * x.size(); ++i) - { - EXPECT_NEAR(K_data_ana[i], K_data_cd[i], eps); - } - } + jac_asm_ana.assembleWithJacobian(loc_asm, t, dt, x, xdot, b_data_ana, + Jac_data_ana); if (LocAsm::asmb) { @@ -467,41 +444,17 @@ struct ProcessLibCentralDifferencesJacobianAssembler : public ::testing::Test double const eps = std::numeric_limits::epsilon(); - std::vector M_data_cd; - std::vector K_data_cd; std::vector b_data_cd; std::vector Jac_data_cd; - std::vector M_data_ana; - std::vector K_data_ana; std::vector b_data_ana; std::vector Jac_data_ana; double const t = 0.0; - jac_asm_cd.assembleWithJacobian(loc_asm, t, dt, x, xdot, M_data_cd, - K_data_cd, b_data_cd, Jac_data_cd); + jac_asm_cd.assembleWithJacobian(loc_asm, t, dt, x, xdot, b_data_cd, + Jac_data_cd); - jac_asm_ana.assembleWithJacobian(loc_asm, t, dt, x, xdot, M_data_ana, - K_data_ana, b_data_ana, Jac_data_ana); - - if (LocAsm::asmM) - { - ASSERT_EQ(x.size() * x.size(), M_data_cd.size()); - ASSERT_EQ(x.size() * x.size(), M_data_ana.size()); - for (std::size_t i = 0; i < x.size() * x.size(); ++i) - { - EXPECT_NEAR(M_data_ana[i], M_data_cd[i], eps); - } - } - - if (LocAsm::asmK) - { - ASSERT_EQ(x.size() * x.size(), K_data_cd.size()); - ASSERT_EQ(x.size() * x.size(), K_data_ana.size()); - for (std::size_t i = 0; i < x.size() * x.size(); ++i) - { - EXPECT_NEAR(K_data_ana[i], K_data_cd[i], eps); - } - } + jac_asm_ana.assembleWithJacobian(loc_asm, t, dt, x, xdot, b_data_ana, + Jac_data_ana); if (LocAsm::asmb) { From b0cfd47ed518e87c5a74f2dcad2c630bec52336e Mon Sep 17 00:00:00 2001 From: Dmitri Naumov Date: Tue, 16 Jul 2024 14:43:33 +0200 Subject: [PATCH 5/6] [PL] Remove M and K from asmWithJacobian functions --- .../ComponentTransportFEM.h | 2 -- .../ComponentTransportProcess.cpp | 4 ++-- .../ComponentTransportProcess.h | 3 +-- ProcessLib/HT/HTProcess.cpp | 4 ++-- ProcessLib/HT/HTProcess.h | 3 +-- ProcessLib/HeatConduction/HeatConductionFEM.h | 2 -- .../HeatConduction/HeatConductionProcess.cpp | 4 ++-- .../HeatConduction/HeatConductionProcess.h | 3 +-- .../HeatTransportBHEProcess.cpp | 4 ++-- .../HeatTransportBHEProcess.h | 3 +-- .../HeatTransportBHELocalAssemblerBHE-impl.h | 4 ++-- .../HeatTransportBHELocalAssemblerBHE.h | 2 -- .../HeatTransportBHELocalAssemblerSoil-impl.h | 7 ++++--- .../HeatTransportBHELocalAssemblerSoil.h | 2 -- .../HydroMechanics/HydroMechanicsFEM-impl.h | 14 ++++++-------- ProcessLib/HydroMechanics/HydroMechanicsFEM.h | 3 --- .../HydroMechanics/HydroMechanicsProcess.cpp | 4 ++-- .../HydroMechanics/HydroMechanicsProcess.h | 3 +-- .../HydroMechanics/HydroMechanicsProcess.cpp | 4 ++-- .../HydroMechanics/HydroMechanicsProcess.h | 3 +-- .../HydroMechanicsLocalAssemblerInterface.h | 2 -- .../SmallDeformationLocalAssemblerInterface.h | 2 -- ...mallDeformationLocalAssemblerMatrix-impl.h | 2 -- .../SmallDeformationLocalAssemblerMatrix.h | 2 -- .../SmallDeformationProcess.cpp | 4 ++-- .../SmallDeformationProcess.h | 3 +-- .../LargeDeformation/LargeDeformationFEM.h | 2 -- .../LargeDeformationProcess.cpp | 8 ++++---- .../LargeDeformationProcess.h | 3 +-- ProcessLib/LiquidFlow/LiquidFlowProcess.cpp | 4 ++-- ProcessLib/LiquidFlow/LiquidFlowProcess.h | 3 +-- ProcessLib/PhaseField/PhaseFieldFEM-impl.h | 2 -- ProcessLib/PhaseField/PhaseFieldFEM.h | 1 - ProcessLib/PhaseField/PhaseFieldProcess.cpp | 4 ++-- ProcessLib/PhaseField/PhaseFieldProcess.h | 3 +-- .../RichardsComponentTransportProcess.cpp | 4 ++-- .../RichardsComponentTransportProcess.h | 3 +-- .../RichardsFlow/RichardsFlowProcess.cpp | 4 ++-- ProcessLib/RichardsFlow/RichardsFlowProcess.h | 3 +-- .../RichardsMechanicsFEM-impl.h | 19 ++++++------------- .../RichardsMechanics/RichardsMechanicsFEM.h | 19 ++----------------- .../RichardsMechanicsProcess.cpp | 4 ++-- .../RichardsMechanicsProcess.h | 3 +-- .../SmallDeformation/SmallDeformationFEM.h | 2 -- .../SmallDeformationProcess.cpp | 8 ++++---- .../SmallDeformationProcess.h | 3 +-- .../SmallDeformationNonlocalFEM.h | 2 -- .../SmallDeformationNonlocalProcess.cpp | 4 ++-- .../SmallDeformationNonlocalProcess.h | 3 +-- ProcessLib/StokesFlow/StokesFlowProcess.cpp | 3 +-- ProcessLib/StokesFlow/StokesFlowProcess.h | 3 +-- ProcessLib/TES/TESProcess.cpp | 4 ++-- ProcessLib/TES/TESProcess.h | 3 +-- ProcessLib/TH2M/TH2MFEM-impl.h | 2 -- ProcessLib/TH2M/TH2MFEM.h | 2 -- ProcessLib/TH2M/TH2MProcess.cpp | 4 ++-- ProcessLib/TH2M/TH2MProcess.h | 3 +-- .../ThermalTwoPhaseFlowWithPPProcess.cpp | 8 ++++---- .../ThermalTwoPhaseFlowWithPPProcess.h | 3 +-- .../ThermoHydroMechanicsFEM-impl.h | 15 +++++++-------- .../ThermoHydroMechanicsFEM.h | 2 -- .../ThermoHydroMechanicsProcess.cpp | 4 ++-- .../ThermoHydroMechanicsProcess.h | 3 +-- .../ThermoMechanicalPhaseFieldFEM-impl.h | 12 ++++++------ .../ThermoMechanicalPhaseFieldFEM.h | 1 - .../ThermoMechanicalPhaseFieldProcess.cpp | 4 ++-- .../ThermoMechanicalPhaseFieldProcess.h | 3 +-- .../ThermoMechanics/ThermoMechanicsFEM-impl.h | 14 ++++++-------- .../ThermoMechanics/ThermoMechanicsFEM.h | 3 --- .../ThermoMechanicsProcess.cpp | 4 ++-- .../ThermoMechanics/ThermoMechanicsProcess.h | 3 +-- .../ThermoRichardsFlowFEM-impl.h | 2 -- .../ThermoRichardsFlowFEM.h | 2 -- .../ThermoRichardsFlowProcess.cpp | 6 +++--- .../ThermoRichardsFlowProcess.h | 3 +-- .../ThermoRichardsMechanicsFEM-impl.h | 2 -- .../ThermoRichardsMechanicsFEM.h | 2 -- .../ThermoRichardsMechanicsProcess.cpp | 5 ++--- .../ThermoRichardsMechanicsProcess.h | 3 +-- .../TwoPhaseFlowWithPPProcess.cpp | 4 ++-- .../TwoPhaseFlowWithPPProcess.h | 3 +-- .../TwoPhaseFlowWithPrhoProcess.cpp | 4 ++-- .../TwoPhaseFlowWithPrhoProcess.h | 3 +-- 83 files changed, 123 insertions(+), 223 deletions(-) diff --git a/ProcessLib/ComponentTransport/ComponentTransportFEM.h b/ProcessLib/ComponentTransport/ComponentTransportFEM.h index d0a90d31864..974ca4be42e 100644 --- a/ProcessLib/ComponentTransport/ComponentTransportFEM.h +++ b/ProcessLib/ComponentTransport/ComponentTransportFEM.h @@ -1285,8 +1285,6 @@ class LocalAssemblerData : public ComponentTransportLocalAssemblerInterface void assembleWithJacobianForStaggeredScheme( double const t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_b_data, std::vector& local_Jac_data) override { diff --git a/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp b/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp index 55c2cdf83c1..b635727bfa3 100644 --- a/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp +++ b/ProcessLib/ComponentTransport/ComponentTransportProcess.cpp @@ -276,7 +276,7 @@ void ComponentTransportProcess::assembleConcreteProcess( void ComponentTransportProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian ComponentTransportProcess."); @@ -296,7 +296,7 @@ void ComponentTransportProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); // b is the negated residumm used in the Newton's method. // Here negating b is to recover the primitive residuum. diff --git a/ProcessLib/ComponentTransport/ComponentTransportProcess.h b/ProcessLib/ComponentTransport/ComponentTransportProcess.h index 224933c61b9..6ed3b112afd 100644 --- a/ProcessLib/ComponentTransport/ComponentTransportProcess.h +++ b/ProcessLib/ComponentTransport/ComponentTransportProcess.h @@ -166,8 +166,7 @@ class ComponentTransportProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preOutputConcreteProcess(const double t, double const dt, std::vector const& x, diff --git a/ProcessLib/HT/HTProcess.cpp b/ProcessLib/HT/HTProcess.cpp index 63d6ddfcc24..2f41777e6ce 100644 --- a/ProcessLib/HT/HTProcess.cpp +++ b/ProcessLib/HT/HTProcess.cpp @@ -112,7 +112,7 @@ void HTProcess::assembleConcreteProcess( void HTProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian HTProcess."); @@ -131,7 +131,7 @@ void HTProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } std::tuple diff --git a/ProcessLib/HT/HTProcess.h b/ProcessLib/HT/HTProcess.h index 950793e0411..2e7c65b9fa3 100644 --- a/ProcessLib/HT/HTProcess.h +++ b/ProcessLib/HT/HTProcess.h @@ -98,8 +98,7 @@ class HTProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; /** * @copydoc ProcessLib::Process::getDOFTableForExtrapolatorData() diff --git a/ProcessLib/HeatConduction/HeatConductionFEM.h b/ProcessLib/HeatConduction/HeatConductionFEM.h index bb93cc45209..4096af46a59 100644 --- a/ProcessLib/HeatConduction/HeatConductionFEM.h +++ b/ProcessLib/HeatConduction/HeatConductionFEM.h @@ -152,8 +152,6 @@ class LocalAssemblerData : public HeatConductionLocalAssemblerInterface void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) override { diff --git a/ProcessLib/HeatConduction/HeatConductionProcess.cpp b/ProcessLib/HeatConduction/HeatConductionProcess.cpp index a12e8a11fa0..727736df4ee 100644 --- a/ProcessLib/HeatConduction/HeatConductionProcess.cpp +++ b/ProcessLib/HeatConduction/HeatConductionProcess.cpp @@ -111,7 +111,7 @@ void HeatConductionProcess::assembleConcreteProcess( void HeatConductionProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian HeatConductionProcess."); @@ -121,7 +121,7 @@ void HeatConductionProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); transformVariableFromGlobalVector(b, 0 /*variable id*/, *_local_to_global_index_map, *_heat_flux, diff --git a/ProcessLib/HeatConduction/HeatConductionProcess.h b/ProcessLib/HeatConduction/HeatConductionProcess.h index c0a461d4377..dc60006530c 100644 --- a/ProcessLib/HeatConduction/HeatConductionProcess.h +++ b/ProcessLib/HeatConduction/HeatConductionProcess.h @@ -70,8 +70,7 @@ class HeatConductionProcess final : public Process const double t, double const /*dt*/, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preOutputConcreteProcess(const double t, double const dt, std::vector const& x, diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp index 03878612a98..5b2a54b1110 100644 --- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp +++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.cpp @@ -177,7 +177,7 @@ void HeatTransportBHEProcess::assembleConcreteProcess( void HeatTransportBHEProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian HeatTransportBHE process."); @@ -188,7 +188,7 @@ void HeatTransportBHEProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } void HeatTransportBHEProcess::computeSecondaryVariableConcrete( diff --git a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h index 328e27a5995..e59d3af65fc 100644 --- a/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h +++ b/ProcessLib/HeatTransportBHE/HeatTransportBHEProcess.h @@ -63,8 +63,7 @@ class HeatTransportBHEProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void createBHEBoundaryConditionTopBottom( std::vector> const& all_bhe_nodes); diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h index f26a5e5d950..e4fb2d485fe 100644 --- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h +++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE-impl.h @@ -219,8 +219,6 @@ void HeatTransportBHELocalAssemblerBHE:: assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_rhs_data, std::vector& local_Jac_data) { @@ -236,6 +234,8 @@ void HeatTransportBHELocalAssemblerBHE:: auto local_rhs = MathLib::createZeroedVector( local_rhs_data, local_matrix_size); + std::vector local_M_data(local_Jac_data.size()); + std::vector local_K_data(local_Jac_data.size()); assemble(t, dt, local_x, local_x_prev, local_M_data, local_K_data, local_rhs_data /*not going to be used*/); diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE.h index bef12205b7e..2664222a5f6 100644 --- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE.h +++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerBHE.h @@ -73,8 +73,6 @@ class HeatTransportBHELocalAssemblerBHE void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h index 56351127972..fc1443e76f9 100644 --- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h +++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil-impl.h @@ -195,9 +195,8 @@ void HeatTransportBHELocalAssemblerSoil::assemble( template void HeatTransportBHELocalAssemblerSoil::assembleWithJacobian( double const t, double const dt, std::vector const& local_x, - std::vector const& local_x_prev, std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_rhs_data, - std::vector& local_Jac_data) + std::vector const& local_x_prev, + std::vector& local_rhs_data, std::vector& local_Jac_data) { assert(local_x.size() == ShapeFunction::NPOINTS); auto const local_matrix_size = local_x.size(); @@ -212,6 +211,8 @@ void HeatTransportBHELocalAssemblerSoil::assembleWithJacobian( auto local_rhs = MathLib::createZeroedVector( local_rhs_data, local_matrix_size); + std::vector local_M_data(local_Jac_data.size()); + std::vector local_K_data(local_Jac_data.size()); assemble(t, dt, local_x, local_x_prev, local_M_data, local_K_data, local_rhs_data /*not going to be used*/); diff --git a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h index 82488e738cd..8eed921cf2e 100644 --- a/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h +++ b/ProcessLib/HeatTransportBHE/LocalAssemblers/HeatTransportBHELocalAssemblerSoil.h @@ -62,8 +62,6 @@ class HeatTransportBHELocalAssemblerSoil void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h index 07094dc53be..867359cb049 100644 --- a/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h +++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM-impl.h @@ -103,8 +103,6 @@ void HydroMechanicsLocalAssembler const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) { @@ -791,12 +789,12 @@ template void HydroMechanicsLocalAssembler:: - assembleWithJacobianForStaggeredScheme( - const double t, double const dt, Eigen::VectorXd const& local_x, - Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, - std::vector& local_b_data, std::vector& local_Jac_data) + assembleWithJacobianForStaggeredScheme(const double t, double const dt, + Eigen::VectorXd const& local_x, + Eigen::VectorXd const& local_x_prev, + int const process_id, + std::vector& local_b_data, + std::vector& local_Jac_data) { // For the equations with pressure if (process_id == _process_data.hydraulic_process_id) diff --git a/ProcessLib/HydroMechanics/HydroMechanicsFEM.h b/ProcessLib/HydroMechanics/HydroMechanicsFEM.h index 45dcdc49d94..fdb3010d97a 100644 --- a/ProcessLib/HydroMechanics/HydroMechanicsFEM.h +++ b/ProcessLib/HydroMechanics/HydroMechanicsFEM.h @@ -215,15 +215,12 @@ class HydroMechanicsLocalAssembler void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) override; void assembleWithJacobianForStaggeredScheme( const double t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp b/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp index 2addf268dd6..d5a123a71ab 100644 --- a/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp +++ b/ProcessLib/HydroMechanics/HydroMechanicsProcess.cpp @@ -308,7 +308,7 @@ void HydroMechanicsProcess:: assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { // For the monolithic scheme bool const use_monolithic_scheme = _process_data.isMonolithicSchemeUsed(); @@ -339,7 +339,7 @@ void HydroMechanicsProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); auto copyRhs = [&](int const variable_id, auto& output_vector) { diff --git a/ProcessLib/HydroMechanics/HydroMechanicsProcess.h b/ProcessLib/HydroMechanics/HydroMechanicsProcess.h index 70df7d70b7e..775fcd7f49e 100644 --- a/ProcessLib/HydroMechanics/HydroMechanicsProcess.h +++ b/ProcessLib/HydroMechanics/HydroMechanicsProcess.h @@ -85,8 +85,7 @@ class HydroMechanicsProcess final : public Process const double t, double const /*dt*/, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& x, double const t, double const dt, diff --git a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp index 976676c913d..dda3ef2d93f 100644 --- a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp +++ b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.cpp @@ -461,7 +461,7 @@ template void HydroMechanicsProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian HydroMechanicsProcess."); @@ -471,7 +471,7 @@ void HydroMechanicsProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); auto copyRhs = [&](int const variable_id, auto& output_vector) { diff --git a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.h b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.h index fe8b621d6ee..884f224705a 100644 --- a/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.h +++ b/ProcessLib/LIE/HydroMechanics/HydroMechanicsProcess.h @@ -73,8 +73,7 @@ class HydroMechanicsProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& x, double const t, double const dt, int const process_id) override; diff --git a/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerInterface.h b/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerInterface.h index 34f37086cec..7bb3281a40a 100644 --- a/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerInterface.h +++ b/ProcessLib/LIE/HydroMechanics/LocalAssembler/HydroMechanicsLocalAssemblerInterface.h @@ -67,8 +67,6 @@ class HydroMechanicsLocalAssemblerInterface void assembleWithJacobian(double const t, double const dt, std::vector const& local_x_, std::vector const& local_x_prev_, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_b_data, std::vector& local_Jac_data) override { diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h index a9c7999bcbc..4c1c98d4539 100644 --- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h +++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerInterface.h @@ -41,8 +41,6 @@ class SmallDeformationLocalAssemblerInterface void assembleWithJacobian(double const t, double const dt, std::vector const& local_x_, std::vector const& /*local_x_prev*/, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_b_data, std::vector& local_Jac_data) override { diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix-impl.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix-impl.h index fc19e7c026c..750193bc5ce 100644 --- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix-impl.h +++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix-impl.h @@ -94,8 +94,6 @@ void SmallDeformationLocalAssemblerMatrix:: assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& /*local_x_prev*/, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_b_data, std::vector& local_Jac_data) { diff --git a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h index b624622d5ea..9ea94d6e315 100644 --- a/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h +++ b/ProcessLib/LIE/SmallDeformation/LocalAssembler/SmallDeformationLocalAssemblerMatrix.h @@ -74,8 +74,6 @@ class SmallDeformationLocalAssemblerMatrix void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& /*local_x_prev*/, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_b_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp index 6cf8ae89652..ce7ae11c8e7 100644 --- a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp +++ b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.cpp @@ -429,7 +429,7 @@ void SmallDeformationProcess:: assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian SmallDeformationProcess."); @@ -439,7 +439,7 @@ void SmallDeformationProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } template void SmallDeformationProcess::preTimestepConcreteProcess( diff --git a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.h b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.h index 61f76db9a4d..50834590970 100644 --- a/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.h +++ b/ProcessLib/LIE/SmallDeformation/SmallDeformationProcess.h @@ -70,8 +70,7 @@ class SmallDeformationProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& x, double const t, double const dt, diff --git a/ProcessLib/LargeDeformation/LargeDeformationFEM.h b/ProcessLib/LargeDeformation/LargeDeformationFEM.h index 31c20f977d0..bde4331a96d 100644 --- a/ProcessLib/LargeDeformation/LargeDeformationFEM.h +++ b/ProcessLib/LargeDeformation/LargeDeformationFEM.h @@ -242,8 +242,6 @@ class LargeDeformationLocalAssembler void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_b_data, std::vector& local_Jac_data) override { diff --git a/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp b/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp index 169ad6ed5e8..bf5a13101ff 100644 --- a/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp +++ b/ProcessLib/LargeDeformation/LargeDeformationProcess.cpp @@ -145,7 +145,7 @@ void LargeDeformationProcess::assembleConcreteProcess( getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, &b); - _global_output(t, process_id, M, K, b); + _global_output(t, process_id, &M, &K, b); } template @@ -153,7 +153,7 @@ void LargeDeformationProcess:: assembleWithJacobianConcreteProcess( double const t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian LargeDeformationProcess."); @@ -163,12 +163,12 @@ void LargeDeformationProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); transformVariableFromGlobalVector(b, 0, *_local_to_global_index_map, *_nodal_forces, std::negate()); - _global_output(t, process_id, M, K, b, &Jac); + _global_output(t, process_id, nullptr, nullptr, b, &Jac); } template diff --git a/ProcessLib/LargeDeformation/LargeDeformationProcess.h b/ProcessLib/LargeDeformation/LargeDeformationProcess.h index c30f829d583..27e41b7e1a8 100644 --- a/ProcessLib/LargeDeformation/LargeDeformationProcess.h +++ b/ProcessLib/LargeDeformation/LargeDeformationProcess.h @@ -58,8 +58,7 @@ class LargeDeformationProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void postTimestepConcreteProcess(std::vector const& x, std::vector const& x_prev, diff --git a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp index af379eab6fc..8a7ed8a18b8 100644 --- a/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp +++ b/ProcessLib/LiquidFlow/LiquidFlowProcess.cpp @@ -97,7 +97,7 @@ void LiquidFlowProcess::assembleConcreteProcess( void LiquidFlowProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian LiquidFlowProcess."); @@ -108,7 +108,7 @@ void LiquidFlowProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } void LiquidFlowProcess::computeSecondaryVariableConcrete( diff --git a/ProcessLib/LiquidFlow/LiquidFlowProcess.h b/ProcessLib/LiquidFlow/LiquidFlowProcess.h index 8ac5910005e..b61365e5b6b 100644 --- a/ProcessLib/LiquidFlow/LiquidFlowProcess.h +++ b/ProcessLib/LiquidFlow/LiquidFlowProcess.h @@ -100,8 +100,7 @@ class LiquidFlowProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; LiquidFlowData _process_data; diff --git a/ProcessLib/PhaseField/PhaseFieldFEM-impl.h b/ProcessLib/PhaseField/PhaseFieldFEM-impl.h index cccd13e0d99..56b62e26e96 100644 --- a/ProcessLib/PhaseField/PhaseFieldFEM-impl.h +++ b/ProcessLib/PhaseField/PhaseFieldFEM-impl.h @@ -22,8 +22,6 @@ void PhaseFieldLocalAssembler:: assembleWithJacobianForStaggeredScheme( double const t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& /*local_x_prev*/, int const process_id, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_b_data, std::vector& local_Jac_data) { // For the equations with phase field. diff --git a/ProcessLib/PhaseField/PhaseFieldFEM.h b/ProcessLib/PhaseField/PhaseFieldFEM.h index a3d1e3170c7..17bfe2b9f1e 100644 --- a/ProcessLib/PhaseField/PhaseFieldFEM.h +++ b/ProcessLib/PhaseField/PhaseFieldFEM.h @@ -277,7 +277,6 @@ class PhaseFieldLocalAssembler : public PhaseFieldLocalAssemblerInterface void assembleWithJacobianForStaggeredScheme( double const t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/PhaseField/PhaseFieldProcess.cpp b/ProcessLib/PhaseField/PhaseFieldProcess.cpp index 1a35046a8fd..38233cf3db2 100644 --- a/ProcessLib/PhaseField/PhaseFieldProcess.cpp +++ b/ProcessLib/PhaseField/PhaseFieldProcess.cpp @@ -218,7 +218,7 @@ template void PhaseFieldProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { std::vector dof_tables; @@ -242,7 +242,7 @@ void PhaseFieldProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); if (process_id == 0) { diff --git a/ProcessLib/PhaseField/PhaseFieldProcess.h b/ProcessLib/PhaseField/PhaseFieldProcess.h index 1d354872871..9f2a2db958e 100644 --- a/ProcessLib/PhaseField/PhaseFieldProcess.h +++ b/ProcessLib/PhaseField/PhaseFieldProcess.h @@ -70,8 +70,7 @@ class PhaseFieldProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& x, double const t, double const dt, diff --git a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp index f034a1c1360..69c43e6d5b1 100644 --- a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp +++ b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.cpp @@ -84,7 +84,7 @@ void RichardsComponentTransportProcess::assembleConcreteProcess( void RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian RichardsComponentTransportProcess."); @@ -95,7 +95,7 @@ void RichardsComponentTransportProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } } // namespace RichardsComponentTransport diff --git a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.h b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.h index d84d50f6318..9e9143fa126 100644 --- a/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.h +++ b/ProcessLib/RichardsComponentTransport/RichardsComponentTransportProcess.h @@ -137,8 +137,7 @@ class RichardsComponentTransportProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; RichardsComponentTransportProcessData _process_data; diff --git a/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp b/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp index e00e2ac3da6..6b41d1f8394 100644 --- a/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp +++ b/ProcessLib/RichardsFlow/RichardsFlowProcess.cpp @@ -78,7 +78,7 @@ void RichardsFlowProcess::assembleConcreteProcess( void RichardsFlowProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian RichardsFlowProcess."); @@ -89,7 +89,7 @@ void RichardsFlowProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } } // namespace RichardsFlow diff --git a/ProcessLib/RichardsFlow/RichardsFlowProcess.h b/ProcessLib/RichardsFlow/RichardsFlowProcess.h index 70bee188df1..dae1ee4bdb1 100644 --- a/ProcessLib/RichardsFlow/RichardsFlowProcess.h +++ b/ProcessLib/RichardsFlow/RichardsFlowProcess.h @@ -56,8 +56,7 @@ class RichardsFlowProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; RichardsFlowProcessData _process_data; diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h index b225b114f76..84be48989be 100644 --- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h +++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM-impl.h @@ -1036,8 +1036,6 @@ void RichardsMechanicsLocalAssembler const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) { @@ -1451,8 +1449,6 @@ void RichardsMechanicsLocalAssembler& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& /*local_b_data*/, std::vector& /*local_Jac_data*/) { @@ -1467,8 +1463,6 @@ void RichardsMechanicsLocalAssembler& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& /*local_b_data*/, std::vector& /*local_Jac_data*/) { @@ -1479,24 +1473,23 @@ template void RichardsMechanicsLocalAssembler:: - assembleWithJacobianForStaggeredScheme( - double const t, double const dt, Eigen::VectorXd const& local_x, - Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& local_M_data, std::vector& local_K_data, - std::vector& local_b_data, std::vector& local_Jac_data) + assembleWithJacobianForStaggeredScheme(double const t, double const dt, + Eigen::VectorXd const& local_x, + Eigen::VectorXd const& local_x_prev, + int const process_id, + std::vector& local_b_data, + std::vector& local_Jac_data) { // For the equations with pressure if (process_id == 0) { assembleWithJacobianForPressureEquations(t, dt, local_x, local_x_prev, - local_M_data, local_K_data, local_b_data, local_Jac_data); return; } // For the equations with deformation assembleWithJacobianForDeformationEquations(t, dt, local_x, local_x_prev, - local_M_data, local_K_data, local_b_data, local_Jac_data); } diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h index f11da1d504d..4cec7a14c8e 100644 --- a/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h +++ b/ProcessLib/RichardsMechanics/RichardsMechanicsFEM.h @@ -104,15 +104,12 @@ class RichardsMechanicsLocalAssembler void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) override; void assembleWithJacobianForStaggeredScheme( double const t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; @@ -186,19 +183,13 @@ class RichardsMechanicsLocalAssembler * @param dt Time increment * @param local_x Nodal values of \f$x\f$ of an element. * @param local_x_prev Nodal values of \f$x_{prev}\f$ of an element. - * @param local_M_data Mass matrix of an element, which takes the form of - * \f$ \int N^T N\mathrm{d}\Omega\f$. Not used. - * @param local_K_data Laplacian matrix of an element, which takes the - * form of \f$ \int (\nabla N)^T K \nabla N\mathrm{d}\Omega\f$. - * Not used. * @param local_b_data Right hand side vector of an element. * @param local_Jac_data Element Jacobian matrix for the Newton-Raphson * method. */ void assembleWithJacobianForDeformationEquations( double const t, double const dt, Eigen::VectorXd const& local_x, - Eigen::VectorXd const& local_x_prev, std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, + Eigen::VectorXd const& local_x_prev, std::vector& local_b_data, std::vector& local_Jac_data); /** @@ -218,19 +209,13 @@ class RichardsMechanicsLocalAssembler * @param dt Time increment * @param local_x Nodal values of \f$x\f$ of an element. * @param local_x_prev Nodal values of \f$x_{prev}\f$ of an element. - * @param local_M_data Mass matrix of an element, which takes the form of - * \f$ \int N^T N\mathrm{d}\Omega\f$. Not used. - * @param local_K_data Laplacian matrix of an element, which takes the - * form of \f$ \int (\nabla N)^T K \nabla N\mathrm{d}\Omega\f$. - * Not used. * @param local_b_data Right hand side vector of an element. * @param local_Jac_data Element Jacobian matrix for the Newton-Raphson * method. */ void assembleWithJacobianForPressureEquations( double const t, double const dt, Eigen::VectorXd const& local_x, - Eigen::VectorXd const& local_x_prev, std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_b_data, + Eigen::VectorXd const& local_x_prev, std::vector& local_b_data, std::vector& local_Jac_data); private: diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp index d04fe1f174a..7d7877c3015 100644 --- a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp +++ b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.cpp @@ -292,7 +292,7 @@ void RichardsMechanicsProcess:: assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { // For the monolithic scheme if (_use_monolithic_scheme) @@ -322,7 +322,7 @@ void RichardsMechanicsProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); auto copyRhs = [&](int const variable_id, auto& output_vector) { diff --git a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.h b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.h index 4bde5cffffd..11212f11c61 100644 --- a/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.h +++ b/ProcessLib/RichardsMechanics/RichardsMechanicsProcess.h @@ -88,8 +88,7 @@ class RichardsMechanicsProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void postTimestepConcreteProcess(std::vector const& x, std::vector const& x_prev, diff --git a/ProcessLib/SmallDeformation/SmallDeformationFEM.h b/ProcessLib/SmallDeformation/SmallDeformationFEM.h index 4a436dfd6da..646d18e319e 100644 --- a/ProcessLib/SmallDeformation/SmallDeformationFEM.h +++ b/ProcessLib/SmallDeformation/SmallDeformationFEM.h @@ -229,8 +229,6 @@ class SmallDeformationLocalAssembler void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_b_data, std::vector& local_Jac_data) override { diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp index dc768172fa7..b588aef7a42 100644 --- a/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp +++ b/ProcessLib/SmallDeformation/SmallDeformationProcess.cpp @@ -147,7 +147,7 @@ void SmallDeformationProcess::assembleConcreteProcess( getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, &b); - _global_output(t, process_id, M, K, b); + _global_output(t, process_id, &M, &K, b); } template @@ -155,7 +155,7 @@ void SmallDeformationProcess:: assembleWithJacobianConcreteProcess( double const t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian SmallDeformationProcess."); @@ -165,12 +165,12 @@ void SmallDeformationProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); transformVariableFromGlobalVector(b, 0, *_local_to_global_index_map, *_nodal_forces, std::negate()); - _global_output(t, process_id, M, K, b, &Jac); + _global_output(t, process_id, nullptr, nullptr, b, &Jac); } template diff --git a/ProcessLib/SmallDeformation/SmallDeformationProcess.h b/ProcessLib/SmallDeformation/SmallDeformationProcess.h index ec5a51ca5a1..3693283872a 100644 --- a/ProcessLib/SmallDeformation/SmallDeformationProcess.h +++ b/ProcessLib/SmallDeformation/SmallDeformationProcess.h @@ -58,8 +58,7 @@ class SmallDeformationProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void postTimestepConcreteProcess(std::vector const& x, std::vector const& x_prev, diff --git a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h index 5c0b788c637..6661b54d76d 100644 --- a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h +++ b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalFEM.h @@ -440,8 +440,6 @@ class SmallDeformationNonlocalLocalAssembler void assembleWithJacobian(double const t, double const /*dt*/, std::vector const& local_x, std::vector const& /*local_x_prev*/, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_b_data, std::vector& local_Jac_data) override { diff --git a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp index cbf799fa300..83751f61293 100644 --- a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp +++ b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.cpp @@ -213,7 +213,7 @@ void SmallDeformationNonlocalProcess:: assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian SmallDeformationNonlocalProcess."); @@ -224,7 +224,7 @@ void SmallDeformationNonlocalProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); b.copyValues(*_nodal_forces); std::transform(_nodal_forces->begin(), _nodal_forces->end(), diff --git a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.h b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.h index f505faffa9e..545feeab101 100644 --- a/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.h +++ b/ProcessLib/SmallDeformationNonlocal/SmallDeformationNonlocalProcess.h @@ -62,8 +62,7 @@ class SmallDeformationNonlocalProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void postTimestepConcreteProcess(std::vector const& x, std::vector const& x_prev, diff --git a/ProcessLib/StokesFlow/StokesFlowProcess.cpp b/ProcessLib/StokesFlow/StokesFlowProcess.cpp index 3fc37e549b2..b34d904f47c 100644 --- a/ProcessLib/StokesFlow/StokesFlowProcess.cpp +++ b/ProcessLib/StokesFlow/StokesFlowProcess.cpp @@ -152,8 +152,7 @@ void StokesFlowProcess::assembleWithJacobianConcreteProcess( const double /*t*/, double const /*dt*/, std::vector const& /*x*/, std::vector const& /*x_prev*/, int const /*process_id*/, - GlobalMatrix& /*M*/, GlobalMatrix& /*K*/, GlobalVector& /*b*/, - GlobalMatrix& /*Jac*/) + GlobalVector& /*b*/, GlobalMatrix& /*Jac*/) { OGS_FATAL( "Assembly of Jacobian matrix has not yet been implemented for " diff --git a/ProcessLib/StokesFlow/StokesFlowProcess.h b/ProcessLib/StokesFlow/StokesFlowProcess.h index e6a3126c32d..fae43b5f7f9 100644 --- a/ProcessLib/StokesFlow/StokesFlowProcess.h +++ b/ProcessLib/StokesFlow/StokesFlowProcess.h @@ -80,8 +80,7 @@ class StokesFlowProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; std::vector _base_nodes; std::unique_ptr _mesh_subset_base_nodes; diff --git a/ProcessLib/TES/TESProcess.cpp b/ProcessLib/TES/TESProcess.cpp index 466026afe15..c2f59694243 100644 --- a/ProcessLib/TES/TESProcess.cpp +++ b/ProcessLib/TES/TESProcess.cpp @@ -223,7 +223,7 @@ void TESProcess::assembleConcreteProcess( void TESProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { std::vector dof_table = { _local_to_global_index_map.get()}; @@ -232,7 +232,7 @@ void TESProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } void TESProcess::preTimestepConcreteProcess(std::vector const& x, diff --git a/ProcessLib/TES/TESProcess.h b/ProcessLib/TES/TESProcess.h index 249021ab363..ab05be02beb 100644 --- a/ProcessLib/TES/TESProcess.h +++ b/ProcessLib/TES/TESProcess.h @@ -68,8 +68,7 @@ class TESProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; GlobalVector const& computeVapourPartialPressure( const double t, diff --git a/ProcessLib/TH2M/TH2MFEM-impl.h b/ProcessLib/TH2M/TH2MFEM-impl.h index 42d127276f9..a39ae5a5af6 100644 --- a/ProcessLib/TH2M/TH2MFEM-impl.h +++ b/ProcessLib/TH2M/TH2MFEM-impl.h @@ -1210,8 +1210,6 @@ void TH2MLocalAssembler const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) { diff --git a/ProcessLib/TH2M/TH2MFEM.h b/ProcessLib/TH2M/TH2MFEM.h index bad9ef32fd5..ac63333df31 100644 --- a/ProcessLib/TH2M/TH2MFEM.h +++ b/ProcessLib/TH2M/TH2MFEM.h @@ -109,8 +109,6 @@ class TH2MLocalAssembler : public LocalAssemblerInterface void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/TH2M/TH2MProcess.cpp b/ProcessLib/TH2M/TH2MProcess.cpp index 3162a1876a4..c7f66f32522 100644 --- a/ProcessLib/TH2M/TH2MProcess.cpp +++ b/ProcessLib/TH2M/TH2MProcess.cpp @@ -252,7 +252,7 @@ template void TH2MProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { if (!_use_monolithic_scheme) { @@ -260,7 +260,7 @@ void TH2MProcess::assembleWithJacobianConcreteProcess( } AssemblyMixin>::assembleWithJacobian( - t, dt, x, x_prev, process_id, M, K, b, Jac); + t, dt, x, x_prev, process_id, b, Jac); } template diff --git a/ProcessLib/TH2M/TH2MProcess.h b/ProcessLib/TH2M/TH2MProcess.h index c1831c70564..6b8942efd44 100644 --- a/ProcessLib/TH2M/TH2MProcess.h +++ b/ProcessLib/TH2M/TH2MProcess.h @@ -91,8 +91,7 @@ class TH2MProcess final : public Process, void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& x, double const t, double const dt, diff --git a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp index 474a0426a91..20165029fdd 100644 --- a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp +++ b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.cpp @@ -96,13 +96,13 @@ void ThermalTwoPhaseFlowWithPPProcess::assembleConcreteProcess( getActiveElementIDs(), dof_table, t, dt, x, x_prev, process_id, &M, &K, &b); - _global_output(t, process_id, M, K, b); + _global_output(t, process_id, &M, &K, b); } void ThermalTwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian ThermalTwoPhaseFlowWithPPProcess."); @@ -113,9 +113,9 @@ void ThermalTwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); - _global_output(t, process_id, M, K, b, &Jac); + _global_output(t, process_id, nullptr, nullptr, b, &Jac); } void ThermalTwoPhaseFlowWithPPProcess::preTimestepConcreteProcess( diff --git a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.h b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.h index d58e6fc1405..9b41114ee3a 100644 --- a/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.h +++ b/ProcessLib/ThermalTwoPhaseFlowWithPP/ThermalTwoPhaseFlowWithPPProcess.h @@ -71,8 +71,7 @@ class ThermalTwoPhaseFlowWithPPProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& x, const double t, const double delta_t, diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h index ae3e8597c75..1475f9c878f 100644 --- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h +++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM-impl.h @@ -549,14 +549,13 @@ ConstitutiveRelationsValues ThermoHydroMechanicsLocalAssembler< template void ThermoHydroMechanicsLocalAssembler< - ShapeFunctionDisplacement, ShapeFunctionPressure, DisplacementDim>:: - assembleWithJacobian(double const t, double const dt, - std::vector const& local_x, - std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, - std::vector& local_rhs_data, - std::vector& local_Jac_data) + ShapeFunctionDisplacement, ShapeFunctionPressure, + DisplacementDim>::assembleWithJacobian(double const t, double const dt, + std::vector const& local_x, + std::vector const& + local_x_prev, + std::vector& local_rhs_data, + std::vector& local_Jac_data) { assert(local_x.size() == pressure_size + displacement_size + temperature_size); diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h index 2fe17a9a636..335757843ea 100644 --- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h +++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsFEM.h @@ -109,8 +109,6 @@ class ThermoHydroMechanicsLocalAssembler void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp index bd2f596277a..2fe1b4fd4cb 100644 --- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp +++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.cpp @@ -355,7 +355,7 @@ void ThermoHydroMechanicsProcess:: assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { // For the monolithic scheme if (_use_monolithic_scheme) @@ -392,7 +392,7 @@ void ThermoHydroMechanicsProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); auto copyRhs = [&](int const variable_id, auto& output_vector) { diff --git a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.h b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.h index afc726448f2..c71b18e4ff8 100644 --- a/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.h +++ b/ProcessLib/ThermoHydroMechanics/ThermoHydroMechanicsProcess.h @@ -87,8 +87,7 @@ class ThermoHydroMechanicsProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& x, double const t, double const dt, diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM-impl.h b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM-impl.h index d2f503df1cf..4f676c7bb73 100644 --- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM-impl.h +++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM-impl.h @@ -19,12 +19,12 @@ namespace ThermoMechanicalPhaseField { template void ThermoMechanicalPhaseFieldLocalAssembler:: - assembleWithJacobianForStaggeredScheme( - double const t, double const dt, Eigen::VectorXd const& local_x, - Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, - std::vector& local_b_data, std::vector& local_Jac_data) + assembleWithJacobianForStaggeredScheme(double const t, double const dt, + Eigen::VectorXd const& local_x, + Eigen::VectorXd const& local_x_prev, + int const process_id, + std::vector& local_b_data, + std::vector& local_Jac_data) { if (process_id == _phase_field_process_id) { diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h index ec0c5fe6fcc..d40b691dee4 100644 --- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h +++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldFEM.h @@ -237,7 +237,6 @@ class ThermoMechanicalPhaseFieldLocalAssembler void assembleWithJacobianForStaggeredScheme( double const t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp index 72a414084ac..fcd0a3d5f94 100644 --- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp +++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.cpp @@ -213,7 +213,7 @@ void ThermoMechanicalPhaseFieldProcess:: assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { // For the staggered scheme if (process_id == _mechanics_related_process_id) @@ -247,7 +247,7 @@ void ThermoMechanicalPhaseFieldProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } template diff --git a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.h b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.h index 62f54a89c1d..ec88268d9b1 100644 --- a/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.h +++ b/ProcessLib/ThermoMechanicalPhaseField/ThermoMechanicalPhaseFieldProcess.h @@ -101,8 +101,7 @@ class ThermoMechanicalPhaseFieldProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& x, double const t, double const dt, diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h index f9e4f067363..125ac03d23c 100644 --- a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h +++ b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM-impl.h @@ -116,8 +116,6 @@ void ThermoMechanicsLocalAssembler:: assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) { @@ -329,12 +327,12 @@ void ThermoMechanicsLocalAssembler:: template void ThermoMechanicsLocalAssembler:: - assembleWithJacobianForStaggeredScheme( - const double t, double const dt, Eigen::VectorXd const& local_x, - Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, - std::vector& local_b_data, std::vector& local_Jac_data) + assembleWithJacobianForStaggeredScheme(const double t, double const dt, + Eigen::VectorXd const& local_x, + Eigen::VectorXd const& local_x_prev, + int const process_id, + std::vector& local_b_data, + std::vector& local_Jac_data) { // For the equations with pressure if (process_id == _process_data.heat_conduction_process_id) diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h index a8a45231d3a..681b6043ec5 100644 --- a/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h +++ b/ProcessLib/ThermoMechanics/ThermoMechanicsFEM.h @@ -182,15 +182,12 @@ class ThermoMechanicsLocalAssembler void assembleWithJacobianForStaggeredScheme( double const t, double const dt, Eigen::VectorXd const& local_x, Eigen::VectorXd const& local_x_prev, int const process_id, - std::vector& local_M_data, std::vector& local_K_data, std::vector& local_b_data, std::vector& local_Jac_data) override; void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& local_M_data, - std::vector& local_K_data, std::vector& local_rhs_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp index 21fa72010ca..3160fca6646 100644 --- a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp +++ b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.cpp @@ -225,7 +225,7 @@ void ThermoMechanicsProcess:: assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleJacobian ThermoMechanicsProcess."); @@ -273,7 +273,7 @@ void ThermoMechanicsProcess:: GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_tables, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); // TODO (naumov): Refactor the copy rhs part. This is copy from HM. auto copyRhs = [&](int const variable_id, auto& output_vector) diff --git a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h index ff35d33c61a..bd5930b7ac7 100644 --- a/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h +++ b/ProcessLib/ThermoMechanics/ThermoMechanicsProcess.h @@ -77,8 +77,7 @@ class ThermoMechanicsProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& x, double const t, double const dt, diff --git a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM-impl.h b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM-impl.h index 58c8b2724fa..3ccc832d9db 100644 --- a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM-impl.h +++ b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM-impl.h @@ -155,8 +155,6 @@ void ThermoRichardsFlowLocalAssembler:: assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) { diff --git a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM.h b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM.h index 8adfb0e0cb8..9bc94a65e9d 100644 --- a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM.h +++ b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowFEM.h @@ -70,8 +70,6 @@ class ThermoRichardsFlowLocalAssembler : public LocalAssemblerInterface void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp index 8c88b5ddc93..762fae97635 100644 --- a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp +++ b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.cpp @@ -144,7 +144,7 @@ void ThermoRichardsFlowProcess::assembleConcreteProcess( void ThermoRichardsFlowProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { std::vector dof_tables; @@ -154,8 +154,8 @@ void ThermoRichardsFlowProcess::assembleWithJacobianConcreteProcess( dof_tables.emplace_back(_local_to_global_index_map.get()); _pvma.assembleWithJacobian(_local_assemblers, getActiveElementIDs(), - dof_tables, t, dt, x, x_prev, process_id, M, K, - b, Jac); + dof_tables, t, dt, x, x_prev, process_id, b, + Jac); auto copyRhs = [&](int const variable_id, auto& output_vector) { diff --git a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.h b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.h index 12cbffd2d99..c4c3cb94266 100644 --- a/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.h +++ b/ProcessLib/ThermoRichardsFlow/ThermoRichardsFlowProcess.h @@ -122,8 +122,7 @@ class ThermoRichardsFlowProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void postTimestepConcreteProcess(std::vector const& x, std::vector const& x_prev, diff --git a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM-impl.h b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM-impl.h index 3a9c9776320..421998d8795 100644 --- a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM-impl.h +++ b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM-impl.h @@ -202,8 +202,6 @@ void ThermoRichardsMechanicsLocalAssembler const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) { diff --git a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM.h b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM.h index 410cad16ad4..0bb2215d1ff 100644 --- a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM.h +++ b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsFEM.h @@ -244,8 +244,6 @@ class ThermoRichardsMechanicsLocalAssembler void assembleWithJacobian(double const t, double const dt, std::vector const& local_x, std::vector const& local_x_prev, - std::vector& /*local_M_data*/, - std::vector& /*local_K_data*/, std::vector& local_rhs_data, std::vector& local_Jac_data) override; diff --git a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp index 3232236c077..3813235dd1b 100644 --- a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp +++ b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.cpp @@ -218,14 +218,13 @@ void ThermoRichardsMechanicsProcess:: assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { AssemblyMixin>::assembleWithJacobian(t, dt, x, x_prev, process_id, - M, K, b, - Jac); + b, Jac); } template diff --git a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.h b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.h index 45a99603e9c..a13c2b0d338 100644 --- a/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.h +++ b/ProcessLib/ThermoRichardsMechanics/ThermoRichardsMechanicsProcess.h @@ -184,8 +184,7 @@ class ThermoRichardsMechanicsProcess final void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; void preTimestepConcreteProcess(std::vector const& /*x*/, const double /*t*/, diff --git a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp index 2ea49f2db01..268193e4644 100644 --- a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp +++ b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.cpp @@ -78,7 +78,7 @@ void TwoPhaseFlowWithPPProcess::assembleConcreteProcess( void TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian TwoPhaseFlowWithPPProcess."); @@ -89,7 +89,7 @@ void TwoPhaseFlowWithPPProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } } // namespace TwoPhaseFlowWithPP diff --git a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.h b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.h index bcbf185c829..d3e85ca1f42 100644 --- a/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.h +++ b/ProcessLib/TwoPhaseFlowWithPP/TwoPhaseFlowWithPPProcess.h @@ -67,8 +67,7 @@ class TwoPhaseFlowWithPPProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; TwoPhaseFlowWithPPProcessData _process_data; diff --git a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp index 1bbd135aafc..6026da1846e 100644 --- a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp +++ b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.cpp @@ -78,7 +78,7 @@ void TwoPhaseFlowWithPrhoProcess::assembleConcreteProcess( void TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, GlobalMatrix& Jac) + GlobalVector& b, GlobalMatrix& Jac) { DBUG("AssembleWithJacobian TwoPhaseFlowWithPrhoProcess."); @@ -89,7 +89,7 @@ void TwoPhaseFlowWithPrhoProcess::assembleWithJacobianConcreteProcess( GlobalExecutor::executeSelectedMemberDereferenced( _global_assembler, &VectorMatrixAssembler::assembleWithJacobian, _local_assemblers, getActiveElementIDs(), dof_table, t, dt, x, x_prev, - process_id, &M, &K, &b, &Jac); + process_id, &b, &Jac); } } // namespace TwoPhaseFlowWithPrho diff --git a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.h b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.h index 14b0b8d8f70..111a0f5c8a7 100644 --- a/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.h +++ b/ProcessLib/TwoPhaseFlowWithPrho/TwoPhaseFlowWithPrhoProcess.h @@ -65,8 +65,7 @@ class TwoPhaseFlowWithPrhoProcess final : public Process void assembleWithJacobianConcreteProcess( const double t, double const dt, std::vector const& x, std::vector const& x_prev, int const process_id, - GlobalMatrix& M, GlobalMatrix& K, GlobalVector& b, - GlobalMatrix& Jac) override; + GlobalVector& b, GlobalMatrix& Jac) override; TwoPhaseFlowWithPrhoProcessData _process_data; From 3e9a70cd397af307a97011096814bcc3e4deae2a Mon Sep 17 00:00:00 2001 From: Dmitri Naumov Date: Wed, 17 Jul 2024 23:20:13 +0200 Subject: [PATCH 6/6] [PL] Separate matrix output for M,K,b and J,b A little code duplication in favor not to deal with pointers. --- ProcessLib/Assembly/MatrixOutput.cpp | 49 ++++++++++++++----- ProcessLib/Assembly/MatrixOutput.h | 10 ++-- .../ParallelVectorMatrixAssembler.cpp | 8 +-- ProcessLib/VectorMatrixAssembler.cpp | 5 +- 4 files changed, 50 insertions(+), 22 deletions(-) diff --git a/ProcessLib/Assembly/MatrixOutput.cpp b/ProcessLib/Assembly/MatrixOutput.cpp index 6e01624cd83..56ed505c69f 100644 --- a/ProcessLib/Assembly/MatrixOutput.cpp +++ b/ProcessLib/Assembly/MatrixOutput.cpp @@ -286,12 +286,11 @@ LocalMatrixOutput::LocalMatrixOutput() outputFilename); } -void LocalMatrixOutput::operator()( - double const t, int const process_id, std::size_t const element_id, - std::vector const* local_M_data, - std::vector const* local_K_data, - std::vector const& local_b_data, - std::vector const* const local_Jac_data) +void LocalMatrixOutput::operator()(double const t, int const process_id, + std::size_t const element_id, + std::vector const& local_M_data, + std::vector const& local_K_data, + std::vector const& local_b_data) { [[likely]] if (!isOutputRequested(element_id)) { @@ -307,16 +306,16 @@ void LocalMatrixOutput::operator()( fmt::print(fh, "## t = {:.15g}, process id = {}, element id = {}\n\n", t, process_id, element_id); - if (local_M_data && !local_M_data->empty()) + if (!local_M_data.empty()) { DBUG("... M"); - fmt::print(fh, "# M\n{}\n\n", toSquareMatrixRowMajor(*local_M_data)); + fmt::print(fh, "# M\n{}\n\n", toSquareMatrixRowMajor(local_M_data)); } - if (local_K_data && !local_K_data->empty()) + if (!local_K_data.empty()) { DBUG("... K"); - fmt::print(fh, "# K\n{}\n\n", toSquareMatrixRowMajor(*local_K_data)); + fmt::print(fh, "# K\n{}\n\n", toSquareMatrixRowMajor(local_K_data)); } if (!local_b_data.empty()) @@ -324,12 +323,38 @@ void LocalMatrixOutput::operator()( DBUG("... b"); fmt::print(fh, "# b\n{}\n\n", MathLib::toVector(local_b_data)); } +} + +void LocalMatrixOutput::operator()(double const t, int const process_id, + std::size_t const element_id, + std::vector const& local_b_data, + std::vector const& local_Jac_data) +{ + [[likely]] if (!isOutputRequested(element_id)) + { + return; + } + + std::lock_guard lock_guard{mutex_}; + + auto& fh = outputFile_; + + DBUG("Writing to local matrix debug output file..."); + + fmt::print(fh, "## t = {:.15g}, process id = {}, element id = {}\n\n", t, + process_id, element_id); + + if (!local_b_data.empty()) + { + DBUG("... b"); + fmt::print(fh, "# b\n{}\n\n", MathLib::toVector(local_b_data)); + } - if (local_Jac_data && !local_Jac_data->empty()) + if (!local_Jac_data.empty()) { DBUG("... Jac"); fmt::print(fh, "# Jac\n{}\n\n\n", - toSquareMatrixRowMajor(*local_Jac_data)); + toSquareMatrixRowMajor(local_Jac_data)); } } diff --git a/ProcessLib/Assembly/MatrixOutput.h b/ProcessLib/Assembly/MatrixOutput.h index aee4e411932..2b3d285f7b1 100644 --- a/ProcessLib/Assembly/MatrixOutput.h +++ b/ProcessLib/Assembly/MatrixOutput.h @@ -45,10 +45,14 @@ struct LocalMatrixOutput void operator()(double const t, int const process_id, std::size_t const element_id, - std::vector const* local_M_data, - std::vector const* local_K_data, + std::vector const& local_M_data, + std::vector const& local_K_data, + std::vector const& local_b_data); + + void operator()(double const t, int const process_id, + std::size_t const element_id, std::vector const& local_b_data, - std::vector const* const local_Jac_data = nullptr); + std::vector const& local_Jac_data); ~LocalMatrixOutput(); diff --git a/ProcessLib/Assembly/ParallelVectorMatrixAssembler.cpp b/ProcessLib/Assembly/ParallelVectorMatrixAssembler.cpp index ee4d4f94b8f..e799f5153ce 100644 --- a/ProcessLib/Assembly/ParallelVectorMatrixAssembler.cpp +++ b/ProcessLib/Assembly/ParallelVectorMatrixAssembler.cpp @@ -205,8 +205,8 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( continue; } - local_matrix_output_(t, process_id, element_id, nullptr, - nullptr, local_b_data, &local_Jac_data); + local_matrix_output_(t, process_id, element_id, local_b_data, + local_Jac_data); } } else @@ -241,8 +241,8 @@ void ParallelVectorMatrixAssembler::assembleWithJacobian( continue; } - local_matrix_output_(t, process_id, element_id, nullptr, - nullptr, local_b_data, &local_Jac_data); + local_matrix_output_(t, process_id, element_id, local_b_data, + local_Jac_data); } } } // OpenMP parallel section diff --git a/ProcessLib/VectorMatrixAssembler.cpp b/ProcessLib/VectorMatrixAssembler.cpp index 35e00a85b33..815a32159b2 100644 --- a/ProcessLib/VectorMatrixAssembler.cpp +++ b/ProcessLib/VectorMatrixAssembler.cpp @@ -99,7 +99,7 @@ void VectorMatrixAssembler::assemble( b->add(indices, _local_b_data); } - _local_output(t, process_id, mesh_item_id, &_local_M_data, &_local_K_data, + _local_output(t, process_id, mesh_item_id, _local_M_data, _local_K_data, _local_b_data); } @@ -169,8 +169,7 @@ void VectorMatrixAssembler::assembleWithJacobian( "errors in the local assembler of the current process."); } - _local_output(t, process_id, mesh_item_id, nullptr, nullptr, _local_b_data, - &_local_Jac_data); + _local_output(t, process_id, mesh_item_id, _local_b_data, _local_Jac_data); } } // namespace ProcessLib