Balance for HT process #2132

TomFischer · 2018-05-28T05:39:46Z

PR implements getFlux in the in the class HTProcess and the local assembler for the HT process. The new feature is checked by two new test cases.

endJunction · 2018-05-29T14:39:00Z

ProcessLib/HT/HTFEM.h

+        fe.computeShapeFunctions(pnt_local_coords.getCoords(), shape_matrices,
+                                 GlobalDim, false);
+        std::vector<double> flux;
+        flux.resize(3);


One could move it down to near the return statement.
I'd say, that returning a Eigen::Vector3d is safer, then the std::vector...

endJunction · 2018-05-29T14:41:10Z

ProcessLib/HT/HTFEM.h

+        NumLib::shapeFunctionInterpolate(local_x, shape_matrices.N, T_int_pt, p_int_pt);
+
+        // set the values into array vars
+        vars[static_cast<int>(MaterialLib::Fluid::PropertyVariableType::T)] =


Would it make sense to pass the vars[T] directly to shapeFunctionInterpolate call?

NumLib::shapeFunctionInterpolate(local_x, shape_matrices.N, vars[static_cast<int>(MaterialLib::Fluid::PropertyVariableType::T)], vars[static_cast<int>(MaterialLib::Fluid::PropertyVariableType::p)]);

endJunction · 2018-05-29T14:43:33Z

ProcessLib/HT/HTFEM.h

+                    MaterialLib::Fluid::FluidPropertyType::Density, vars);
+            auto const b = this->_material_properties.specific_body_force;
+            // here it is assumed that the vector b is directed 'downwards'
+            q += K_over_mu * rho_w * b;


I'm not getting the intentions of the comment. What is 'downwards'? Maybe it's copied from some other context...

endJunction · 2018-05-29T14:45:44Z

ProcessLib/HT/HTProcess.cpp

+        DBUG("This is the thermal part of the staggered HTProcess.");
+        return;
+    }
+    if (_balance_mesh)  // computing the balance is optional


if (!_balance_mesh) { return; }

endJunction · 2018-05-29T14:47:29Z

ProcessLib/HT/HTProcess.cpp

+        auto* const balance_pv = MeshLib::getOrCreateMeshProperty<double>(
+            *_balance_mesh, _balance_pv_name, MeshLib::MeshItemType::Cell, 1);
+        // initialise the PropertyVector pv with zero values
+        std::fill(balance_pv->begin(), balance_pv->end(), 0.0);


Not sure here;
If the properties are set, then use NaN for initialization. If the properties are '+='-ed, then setting them 0 is of course necessary.

The properties are '+='-ed.

endJunction

Smaller comments, mostly for discussion.
Didn't look into the ctests.

Needed for the description of the balance data.

chleh

Two small comments and two bigger ones:

time dependence?
mass vs. heat flux?

chleh · 2018-05-30T11:09:55Z

ProcessLib/GroundwaterFlow/GroundwaterFlowFEM.h

@@ -105,7 +105,7 @@ class LocalAssemblerData : public GroundwaterFlowLocalAssemblerInterface
    /// Computes the flux in the point \c p_local_coords that is given in local
    /// coordinates using the values from \c local_x.
    // TODO add time dependency
-    std::vector<double> getFlux(


~~Isn't the flux 2D in two dimentsions?~~

chleh · 2018-05-30T11:14:28Z

ProcessLib/HT/HTFEM.h

@@ -88,6 +88,73 @@ class HTFEM : public HTLocalAssemblerInterface
        return Eigen::Map<const Eigen::RowVectorXd>(N.data(), N.size());
    }

+    /// Computes the flux in the point \c pnt_local_coords that is given in
+    /// local coordinates using the values from \c local_x.
+    // TODO add time dependency


Yep, that would be nice (for the sake of completeness). Could you please add it now?

chleh · 2018-05-30T11:17:46Z

ProcessLib/HT/HTFEM.h

+        auto const K =
+            this->_material_properties.porous_media_properties
+                .getIntrinsicPermeability(t, pos)
+                .getValue(t, pos, 0.0,


Is the 0.0 correct here?

Yes, I think so.

OK, class Permeability says that this can have an arbitrary value.

chleh · 2018-05-30T11:23:34Z

ProcessLib/HT/HTFEM.h

+        }
+
+        Eigen::Vector3d flux;
+        flux.head(GlobalDim) = q;


Maybe here the templated version head>GlobalDim>() would be a tiny bit better.
And I see: That's why flux is always 3D.

Regarding the flux being 3D always: I thought it would be easier in general to be it of constant size, like the node coordinates.

flux.head<GlobalDim>() would be better, that's true!

chleh · 2018-05-30T11:31:35Z

ProcessLib/HT/HTFEM.h

+    /// Computes the flux in the point \c pnt_local_coords that is given in
+    /// local coordinates using the values from \c local_x.
+    // TODO add time dependency
+    Eigen::Vector3d getFlux(


For the HT process the problem is that there are potentially two different interesting fluxes:
The mass/volume flux and the heat flux.
Maybe there should be different names for them?

Or one could pass the variable_id to getFlux?

You are right, there are potentially two different interesting fluxes. Maybe other processes have also interesting quantities fluxes. For this reason we would need rather general names. I opened an issue for discussion and for the time being I would keep the implementation of this PR. Is this okay for you?

TomFischer · 2018-05-30T13:14:37Z

Added the time to all getFlux and the balance functions, see commit a823565. For the second bigger comment I added the issue https://github.com/ufz/ogs/issues/2134 for discussion.

chleh

👍 thanks for the changes.

ogsbot · 2020-06-19T12:52:51Z

OpenGeoSys development has been moved to GitLab.

See this pull request on GitLab.

TomFischer force-pushed the BalanceForHTProcess branch from 829473b to 100c15f Compare May 28, 2018 05:50

TomFischer added the WIP 🏗️ label May 28, 2018

TomFischer force-pushed the BalanceForHTProcess branch 2 times, most recently from 3ab0a3c to d64bad6 Compare May 28, 2018 09:42

TomFischer removed the WIP 🏗️ label May 28, 2018

endJunction reviewed May 29, 2018

View reviewed changes

endJunction approved these changes May 29, 2018

View reviewed changes

TomFischer added 11 commits May 30, 2018 12:23

[PL] HT: Implement HTProcess::getFlux().

81d275e

[PL] Clang format parts of HTProcess.

d8b55d5

[PL] Clang format parts of HTLocalAssemblerInterface.

744ddd4

[PL] Implement getFlux() in HT process local asm.

38413ee

[PL/HT] Enable reading information for the balance.

7546ae2

[PL/HT] Enable the balance calculation in HT.

f608621

[PL] Add t and dt to postTimestep{ConcreteProcess}.

2689f06

Needed for the description of the balance data.

[PL] Include time in balance output file name.

1d9aa9a

[PL/HT] Add tests for the balance/flux calculation.

0be9c2b

[PL] Early exit if calc. of balance isn't needed.

75fa95f

[PL/CalculateSurfaceFlux] Correct comment.

fdeafde

TomFischer force-pushed the BalanceForHTProcess branch from 59a9966 to 3901500 Compare May 30, 2018 10:23

chleh reviewed May 30, 2018

View reviewed changes

[PL] Change return type of getFlux().

3f11cb8

TomFischer force-pushed the BalanceForHTProcess branch from 3901500 to a823565 Compare May 30, 2018 13:09

endJunction approved these changes May 30, 2018

View reviewed changes

chleh approved these changes May 30, 2018

View reviewed changes

[PL] Time dependency for the flux calculation.

dd8c5ca

TomFischer force-pushed the BalanceForHTProcess branch from a823565 to dd8c5ca Compare May 31, 2018 03:54

TomFischer merged commit f2c767a into ufz:master May 31, 2018

TomFischer deleted the BalanceForHTProcess branch May 31, 2018 04:36

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Balance for HT process #2132

Balance for HT process #2132

TomFischer commented May 28, 2018

endJunction May 29, 2018

endJunction May 29, 2018 •

edited

Loading

endJunction May 29, 2018

endJunction May 29, 2018

endJunction May 29, 2018

TomFischer May 30, 2018

endJunction left a comment

chleh left a comment

chleh May 30, 2018

chleh May 30, 2018

chleh May 30, 2018

TomFischer May 30, 2018

chleh May 30, 2018 •

edited

Loading

chleh May 30, 2018

endJunction May 30, 2018

chleh May 30, 2018

endJunction May 30, 2018

TomFischer May 30, 2018

chleh May 30, 2018

TomFischer commented May 30, 2018

chleh left a comment

ogsbot commented Jun 19, 2020

Balance for HT process #2132

Balance for HT process #2132

Conversation

TomFischer commented May 28, 2018

Choose a reason for hiding this comment

endJunction May 29, 2018 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

endJunction left a comment

Choose a reason for hiding this comment

chleh left a comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

chleh May 30, 2018 • edited Loading

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

Choose a reason for hiding this comment

TomFischer commented May 30, 2018

chleh left a comment

Choose a reason for hiding this comment

ogsbot commented Jun 19, 2020

endJunction May 29, 2018 •

edited

Loading

chleh May 30, 2018 •

edited

Loading