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

[PL] Multi-component transport process #2304

Merged
merged 12 commits into from
Jan 23, 2019

Conversation

renchao-lu
Copy link
Contributor

This pull request aims at extending the functionality of component transport process which now can be applicable for multiple component migration. Note that, this work is done using the newly developed multi-phase/-component library (#2303).

@renchao-lu renchao-lu force-pushed the MultiComponentTransportProcess branch 3 times, most recently from 5d41928 to 8300206 Compare December 21, 2018 14:16
@renchao-lu renchao-lu force-pushed the MultiComponentTransportProcess branch 3 times, most recently from e39af45 to a8bcf7b Compare January 16, 2019 23:54
@renchao-lu renchao-lu force-pushed the MultiComponentTransportProcess branch 3 times, most recently from e12cc6c to 440b085 Compare January 17, 2019 19:06
static const int pressure_index = 0;
static const int pressure_size = ShapeFunction::NPOINTS;
// per component
static const int concentration_size = ShapeFunction::NPOINTS;
Copy link
Member

@endJunction endJunction Jan 18, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please decide how to call a thing: component or concentration. Stick to it.
Maybe single_component_size.
It might be, that I don't yet grasp the difference btw the component and concentration.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

rewritten the comment.

@@ -131,48 +152,119 @@ class LocalAssemblerData : public ComponentTransportLocalAssemblerInterface
auto local_b = MathLib::createZeroedVector<LocalVectorType>(
local_b_data, local_matrix_size);

// Nodal DOFs include pressure
auto const num_comp = num_nodal_dof - 1;
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
auto const num_comp = num_nodal_dof - 1;
auto const number_of_components = num_nodal_dof - 1;

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Fixed.

auto local_p = Eigen::Map<const NodalVectorType>(
&local_x[pressure_index], pressure_size);

for (int comp_id = 0; comp_id < num_comp; ++comp_id)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
for (int comp_id = 0; comp_id < num_comp; ++comp_id)
for (int component_id = 0; component_id < number_of_components; ++component_id)

Also I suggest to keep variable's usage and definition close: move definition of the number_of_components before the for-loop.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Fixed and stuck the condition variable 'number_of_component' with the follow-up for-loop.


// Use the fluid density model to compute the density
// TODO: concentration of which component as the argument for
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
// TODO: concentration of which component as the argument for
// TODO (renchao): concentration of which component as the argument for

Copy link
Contributor Author

@renchao-lu renchao-lu Jan 18, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Added.

auto const concentration_index = 1 * ShapeFunction::NPOINTS;
// assuming that fluid density always depends on concentration of the
// first component.
// get local_C and local_p
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Drop this line. The comment is useless above (where it is correct). Here it is wrong. Don't copy-paste comments.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Removed.


// get all process variables stored in a vector before allocation
// pressure first, concentration then
auto const unalloc_process_variables = findProcessVariables(
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What do you mean by unalloc?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Given a more accurate name for this variable.

namespace ProcessLib
{
namespace ComponentTransport
{
std::unique_ptr<Process> createComponentTransportProcess(
MeshLib::Mesh& mesh,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

clang-format

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

reformatted

auto var_names = pv_config.getConfigParameterList<std::string>(tag);

if (var_names.empty())
OGS_FATAL("Config tag <%s> is not found.", tag.c_str());
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Message is wrong. The tag was found, but the result is empty.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Rewritten the error message.

var_name.c_str(), tag.c_str());
}
DBUG("Found process variable \'%s\' for config tag <%s>.",
variable->getName().c_str(), tag.c_str());
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This part from find_if until here is same as in the
findProcessVariable(std::vector<ProcessVariable> const&, BaseLib::ConfigTree const&, std::string const&)
function. Duplications should be avoided.

Copy link
Contributor Author

@renchao-lu renchao-lu Jan 18, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Totally agree with you. Since the utility of the newly-added function has covered that of function findProcessVariable(std::vector const&, BaseLib::ConfigTree const&, std::string const&), we may consider to drop this function and use the new function in all the processes afterwards.

#include "ProcessLib/Utils/InitShapeMatrices.h"

namespace ProcessLib
{
class ProcessVariable;
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Since you include the file 'ProcessVariable.h' in line 28 I think the forward declaration isn't needed.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Dropped the unneeded forward declaration.

@@ -85,14 +97,20 @@ class LocalAssemblerData : public ComponentTransportLocalAssemblerInterface
std::size_t const local_matrix_size,
bool is_axially_symmetric,
unsigned const integration_order,
ComponentTransportProcessData const& process_data)
ComponentTransportProcessData const& process_data,
std::vector<std::vector<std::reference_wrapper<ProcessVariable>>> const&
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I think it is possible to pass
std::vector<std::reference_wrapper<ProcessVariable>> const& process_variables
instead of std::vector<std::vector<std::reference_wrapper<ProcessVariable>>>.
Then you could write later on (for instance in line 113):
_process_variables.size().

Copy link
Contributor Author

@renchao-lu renchao-lu Jan 18, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

_process_variables is derived from the base class of Process. If using std::vector<std::reference_wrapper> instead, we may encounter problems when adopting staggered scheme.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Dropped the two-dimensional vector.

* | c2p | 0 | c2c2| 0 |
* |-----|-----|-----|-----|
* | c3p | 0 | 0 | c3c3|
*/
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Nice documentation 🍰

auto local_C = Eigen::Map<const NodalVectorType>(
&local_x[concentration_index], concentration_size);

// Prevent duplicate computation in loops
Copy link
Member

@TomFischer TomFischer Jan 18, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What does this mean? As far as I can see the pure pressure parts, i.e., Mpp and Kpp are recomputed for each component. Is this avoidable?

Copy link
Contributor Author

@renchao-lu renchao-lu Jan 19, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Fixed. Block matrices Mpp, Kpp, and bp are computed in the first loop over components.

@renchao-lu renchao-lu force-pushed the MultiComponentTransportProcess branch from 440b085 to 9f92a7f Compare January 19, 2019 15:32
}
}

void assembleBlockMatrices(
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, clang format the complete function signature.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I rerun clang-format. There is no changes in function signature, though.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Seems to be a bug in clang format. If you put all the code of the function signature in one line and rerun clang format the result will change.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It works! I get the right clang-format outputs.


if (it != collected_process_variables.end())
OGS_FATAL(
"Number of components for process variable \"%s\" should be 1 "
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, use ' instead of \" for consistency with other message.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Fixed.


// Use the fluid density model to compute the density
// TODO (renchao): concentration of which component as the argument
// for calculation of fluid density
vars[static_cast<int>(
MaterialLib::Fluid::PropertyVariableType::C)] = C_int_pt;
vars[static_cast<int>(
MaterialLib::Fluid::PropertyVariableType::p)] = p_int_pt;
auto const density = _process_data.fluid_properties->getValue(
Copy link
Member

@TomFischer TomFischer Jan 21, 2019

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Here, the density calculation depends only on one of the concentrations. Is this correct? Sorry, I read the comment to late.

@@ -43,17 +44,66 @@ std::vector<std::reference_wrapper<ProcessVariable>> findProcessVariables(
std::vector<ProcessVariable> const& variables,
BaseLib::ConfigTree const& pv_config,
std::initializer_list<std::string>
tag_names)
tags)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, clang format.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I checked that too, it is the way how clang-format outputs it.

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Moving the function signature into one line and then running clang-format will produce what you expect. Once I move into next step to make a commit, the pre-commit hook will tell me that the staged content is not formatted correctly. I would suggest to keep it so far.

@@ -40,7 +40,12 @@ std::vector<std::reference_wrapper<ProcessVariable>> findProcessVariables(
std::vector<ProcessVariable> const& variables,
BaseLib::ConfigTree const& process_config,
std::initializer_list<std::string>
tag_names);
tags);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please, clang format.

@@ -201,16 +216,16 @@
<linear_solvers>
<linear_solver>
<name>general_linear_solver</name>
<lis>-i bicgstab -p ilut -tol 1e-8 -maxiter 20000</lis>
<lis>-i bicgstab -p ilut -tol 1e-16 -maxiter 20000</lis>
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Why did you choose such a strong tolerance criterion (1e-16)?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I have no idea on the setting of linear solver... @endJunction

@TomFischer
Copy link
Member

Please remove 'Tests/Data/Parabolic/ComponentTransport/SimpleSynthetics/Runtime_Record_Sheet.pdf' and 'Tests/Data/Parabolic/ComponentTransport/SimpleSynthetics/Runtime_Record_Sheet.pdf.x'.

vdbc_pcs_0_ts_9990_t_90000.000000_expected.vtu vdbc_pcs_0_ts_9990_t_90000.000000.vtu concentration Si 1e-5 1e-4
vdbc_pcs_0_ts_15990_t_150000.000000_expected.vtu vdbc_pcs_0_ts_15990_t_150000.000000.vtu concentration Si 1e-5 1e-4
vdbc_pcs_0_ts_21990_t_210000.000000_expected.vtu vdbc_pcs_0_ts_21990_t_210000.000000.vtu concentration Si 1e-5 1e-4
vdbc_pcs_0_ts_25990_t_250000.000000_expected.vtu vdbc_pcs_0_ts_25990_t_250000.000000.vtu concentration Si 1e-5 1e-4
)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Are the results of the new test cases (ConcentrationDiffusionOnly_10Comp.prj and ConcentrationDiffusionOnly_30Comp.prj) compared with expected results somewhere in Tests.cmake?

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Supplemented a more simple benchmark for testing which includes 3 components.

Copy link
Member

@TomFischer TomFischer left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Some things to discuss and fix. When tests are green and the fixes are done: ⏩

@renchao-lu
Copy link
Contributor Author

Two pdf files have been removed. @TomFischer

@renchao-lu renchao-lu force-pushed the MultiComponentTransportProcess branch 3 times, most recently from 1f236ca to 48a1eae Compare January 22, 2019 12:24
@@ -50,7 +50,8 @@ void ComponentTransportProcess::initializeConcreteProcess(
ProcessLib::createLocalAssemblers<LocalAssemblerData>(
mesh.getDimension(), mesh.getElements(), dof_table,
pv.getShapeFunctionOrder(), _local_assemblers,
mesh.isAxiallySymmetric(), integration_order, _process_data);
mesh.isAxiallySymmetric(), integration_order, _process_data,
_process_variables[0]);
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Suggested change
_process_variables[0]);
_process_variables[process_id]);

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Changed.

@renchao-lu renchao-lu force-pushed the MultiComponentTransportProcess branch from 656f3f1 to ad2460e Compare January 22, 2019 16:34
@endJunction
Copy link
Member

@renchao-lu When you finish your fixes, please tell, so we can finish the review. As I see now, there are still open questions.

@renchao-lu renchao-lu force-pushed the MultiComponentTransportProcess branch from 5421278 to 5f8250a Compare January 22, 2019 20:42
@renchao-lu renchao-lu force-pushed the MultiComponentTransportProcess branch from 5f8250a to d7a2f41 Compare January 22, 2019 21:53
@TomFischer TomFischer merged commit 5d7e608 into ufz:master Jan 23, 2019
@renchao-lu renchao-lu deleted the MultiComponentTransportProcess branch January 24, 2019 11:19
@ogsbot
Copy link
Member

ogsbot commented Jun 19, 2020

OpenGeoSys development has been moved to GitLab.

See this pull request on GitLab.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

4 participants