Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Staggered scheme for coupled processes with different orders of elements #2016

Merged
merged 33 commits into from Jan 20, 2018

Conversation

Projects
None yet
4 participants
@wenqing
Copy link
Member

wenqing commented Dec 8, 2017

This PR presents the framework for the staggered scheme for the case of coupled processes with different orders of elements.

As a first implementation within the framework, the staggered scheme is realized in HydroMechanicsProcess. The implementation is tested by solving an example with both of the staggered and monolithic schemes.

Besides, the calculation of stress and strain for output and also for the coupling is corrected.

Since the convergence criterion of the coupling loop is changed in this PR, the solutions of the existing benchmark about the staggered scheme have slightly difference from that by the code in the master branch. Therefore, the solutions of the benchmark in Parabolic/HT/StaggeredCoupling/ConstViscosity have been updated in this PR too.

@wenqing wenqing force-pushed the wenqing:staggeredHM branch 2 times, most recently from a0f7d38 to 508586c Dec 14, 2017

@wenqing wenqing force-pushed the wenqing:staggeredHM branch from 508586c to 048da58 Dec 15, 2017

@wenqing wenqing removed the WIP 🏗 label Dec 15, 2017

@wenqing wenqing force-pushed the wenqing:staggeredHM branch 2 times, most recently from f691ae3 to 1fa1aab Dec 18, 2017

@wenqing

This comment has been minimized.

Copy link
Member Author

wenqing commented Dec 18, 2017

@bilke
What is the error in CI and appveyor tests with the following message?

VTK/7.1.0@bilke/stable: WARN: Can't find a 'VTK/7.1.0@bilke/stable' package for the specified options and settings:

- Settings: arch=x86_64, build_type=Release, compiler=gcc, compiler.libcxx=libstdc++, compiler.version=4.9, os=Linux

- Options: fPIC=False, mpi=False, qt=False, shared=False, x11=True
@bilke

This comment has been minimized.

Copy link
Member

bilke commented Dec 18, 2017

@wenqing Should be fixed now!

@@ -0,0 +1,73 @@
---

This comment has been minimized.

@bilke

bilke Dec 18, 2017

Member

This file needs to have a so called frontmatter at the beginning. See other for reference.

This comment has been minimized.

@wenqing

wenqing Dec 19, 2017

Author Member

I added a formatter manually, please take a look.

@wenqing wenqing force-pushed the wenqing:staggeredHM branch 2 times, most recently from dce1eea to 8f49fc5 Dec 19, 2017

@wenqing

This comment has been minimized.

Copy link
Member Author

wenqing commented Dec 19, 2017

@bilke Thanks for fixing.

@bilke

This comment has been minimized.

Copy link
Member

bilke commented Dec 19, 2017

@wenqing There are still some error I have to fix left.

@wenqing wenqing force-pushed the wenqing:staggeredHM branch from 8f49fc5 to f977a6b Dec 20, 2017

@@ -0,0 +1,2 @@
Defines the convergence criteria of the global un-coupling loop of the staggered

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

"un-coupling loop" is not a common phrase. Do you mean weakly coupled processes:
"Defines the convergence criterion for each weakly coupled process." --- would be this correct?

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Changed accordingly.

@@ -76,12 +76,12 @@ class TimeDiscretizedODESystem<ODESystemTag::FirstOrderImplicitQuasilinear,
* \param ode the ODE to be wrapped.
* \param time_discretization the time discretization to be used.
*/
explicit TimeDiscretizedODESystem(ODE& ode, TimeDisc& time_discretization);
explicit TimeDiscretizedODESystem(const int equation_id, ODE& ode,
TimeDisc& time_discretization);

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

equation_id must be commented. There are corresponding doxygen warnings.

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Commented.

@@ -34,6 +37,14 @@ class AnalyticalJacobianAssembler final : public AbstractJacobianAssembler
const double dx_dx, std::vector<double>& local_M_data,
std::vector<double>& local_K_data, std::vector<double>& local_b_data,
std::vector<double>& local_Jac_data) override;

void assembleWithJacobianAndCoupling(

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

docs missing.
Maybe it would be better to call this assembleWithJacobianWeakCoupling because "Coupling" is not a product of this function like "Jacobian".

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Docs was added and the name was changed as well.

const std::vector<GlobalIndexType>& indices);
const std::vector<
std::reference_wrapper<const std::vector<GlobalIndexType>>>&
indices);

std::vector<std::vector<double>> getCurrentLocalSolutions(

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

Just for better understanding: why is Solutions used in plural form? Could this be somehow documented?

I understand that this is not introduced in this PR.

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Documented these two functions.

@@ -69,9 +70,12 @@ class GroundwaterFlowProcess final : public Process
MeshLib::addPropertyToMesh(*_balance_mesh, _balance_pv_name,
MeshLib::MeshItemType::Cell, 1,
init_values);
const int process_id = 0;

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

This I don't understand: on one side process id is 0 but the passed function's argument is not used.... I'd expect that the callee is deciding which process_id is being used...

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Added a comment.

// For the staggered scheme, both processes are assumed to use the same
// element order.
const int process_id = 0;
ProcessLib::ProcessVariable const& pv = getProcessVariables(process_id)[0];

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

How is the above comment related to the process_id? I'd expect the general comments on implementation to be in the header file and in doxygen comments....

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Added more a detailed description.

@@ -160,5 +165,54 @@ void HTProcess::setCoupledTermForTheStaggeredSchemeToLocalAssemblers()
_local_assemblers, _coupled_solutions);
}

NumLib::LocalToGlobalIndexMap* HTProcess::getDOFTableForExtrapolatorData(
bool& manage_storage) const

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

Because manage_storage is used only as output parameter, it would be better to return it too:

std::tuple<bool, DOF*> f()
{
   if (A) return {false, dof};
   return {true, other_dof};
}

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Changed the return type to std::tuple.


manage_storage = true;

return new NumLib::LocalToGlobalIndexMap(

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

Please use automatic storage management for on-heap allocated objects.

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Explained in the discussion.

for (std::size_t i = 0; i < number_of_coupled_solutions; i++)
{
const auto x_t0 = _xs_previous_timestep[process_id].get();
if (x_t0 == nullptr)

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

Minor:
if the type of _xs_previous_timestep[process_id] is std::unique_ptr or similar, there is no need to get the pointer for comparison or later access:

auto const& x_t0 = xs[id];
if (!x_t0) OGS_FATAL(...);

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

removed .get().

/// It only performs for the staggered scheme.
void setCoupledSolutionsOfPreviousTimeStep();

NumLib::LocalToGlobalIndexMap* getDOFTableForExtrapolatorData(

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

docs missing.

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Added a documentation.

@@ -566,6 +582,20 @@ class HydroMechanicsLocalAssembler : public LocalAssemblerInterface
return cache;
}

void assembleWithJacobianForDeformationEquations(

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

docs missing here and below for pressure eqs.

This comment has been minimized.

@wenqing

wenqing Jan 11, 2018

Author Member

Added.

@@ -595,3 +625,5 @@ class HydroMechanicsLocalAssembler : public LocalAssemblerInterface

} // namespace HydroMechanics
} // namespace ProcessLib

#include "HydroMechanicsFEM-impl.h"

This comment has been minimized.

@endJunction

endJunction Jan 5, 2018

Member

I'd either keep all implementations in the HMFEM.h or move all implementations to HMFEM-impl.h, but not halfway, with my preference being the former case.

This comment has been minimized.

@wenqing

wenqing Jan 10, 2018

Author Member

If it is the former case, we have a large header file mixed with declarations and definitions. I will move the definitions of some large members to -impl.h file.

@wenqing wenqing force-pushed the wenqing:staggeredHM branch 3 times, most recently from 30dabe6 to f697954 Jan 18, 2018

@wenqing

This comment has been minimized.

Copy link
Member Author

wenqing commented Jan 18, 2018

@chleh All your comments have been addressed (see the last commit).

@wenqing

This comment has been minimized.

Copy link
Member Author

wenqing commented Jan 18, 2018

@endJunction and @chleh Thanks your review and comments.

@@ -343,6 +360,11 @@ std::unique_ptr<UncoupledProcessesTimeLoop> createUncoupledProcessesTimeLoop(
//! \ogs_file_param{prj__time_loop__processes}
config.getConfigSubtree("processes"), processes, nonlinear_solvers);

if (coupling_config)
{
assert(global_coupling_conv_criteria.size() == per_process_data.size());

This comment has been minimized.

@chleh

chleh Jan 18, 2018

Collaborator

Thanks for adding the check. Actually I'd prefer using OGS_FATAL, because this is a cheap check which is performed only once and the outcome depends on user input. Therefore the check should also be performed in release mode.

This comment has been minimized.

@wenqing

wenqing Jan 18, 2018

Author Member

Changed to use OGS_FATAL.

@chleh

chleh approved these changes Jan 18, 2018

Copy link
Collaborator

chleh left a comment

Thanks for the changes. If you want maybe you might change the assert to OGS_FATAL.

@wenqing wenqing force-pushed the wenqing:staggeredHM branch 2 times, most recently from c601364 to 0401689 Jan 18, 2018

@wenqing wenqing force-pushed the wenqing:staggeredHM branch from 0401689 to 15aa457 Jan 18, 2018

{
if (global_coupling_conv_criteria.size() != per_process_data.size())
{
OGS_FATAL(

This comment has been minimized.

@wenqing

wenqing Jan 18, 2018

Author Member

Changed to OGS_FATAL

@endJunction endJunction merged commit 1876ffd into ufz:master Jan 20, 2018

2 checks passed

continuous-integration/appveyor/pr AppVeyor build succeeded
Details
continuous-integration/jenkins/pr-merge This commit looks good
Details
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.