| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | // -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- | ||
| 2 | // vi: set et ts=4 sw=4 sts=4: | ||
| 3 | // | ||
| 4 | // SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder | ||
| 5 | // SPDX-License-Identifier: GPL-3.0-or-later | ||
| 6 | // | ||
| 7 | |||
| 8 | /*! | ||
| 9 | * \file | ||
| 10 | * \ingroup Fluidmatrixinteractions | ||
| 11 | * \ingroup PoreNetworkModels | ||
| 12 | * \brief Specification of threshold capillary pressures for the PNM. | ||
| 13 | */ | ||
| 14 | #ifndef DUMUX_PNM_THRESHOLD_CAPILLARY_PRESSURES_HH | ||
| 15 | #define DUMUX_PNM_THRESHOLD_CAPILLARY_PRESSURES_HH | ||
| 16 | |||
| 17 | #include <cmath> | ||
| 18 | |||
| 19 | #include <dune/common/exceptions.hh> | ||
| 20 | #include <dumux/porenetwork/common/throatproperties.hh> | ||
| 21 | |||
| 22 | namespace Dumux::PoreNetwork::ThresholdCapillaryPressures { | ||
| 23 | |||
| 24 | /*! \brief The snap-off capillary pressure of a pore throat with regular cross section shape | ||
| 25 | * (with same corner angles and side length) | ||
| 26 | * For details, see Eq. 4.8 in Multiphase Flow in | ||
| 27 | * Blunt, M. J. (2017). Multiphase flow in permeable media: A pore-scale perspective. Cambridge university press. | ||
| 28 | * https://doi.org/10.1017/9781316145098 | ||
| 29 | */ | ||
| 30 | template <class Scalar> | ||
| 31 | 959946 | Scalar pcSnapoffRegularShape(const Scalar surfaceTension, | |
| 32 | const Scalar contactAngle, | ||
| 33 | const Scalar inscribedRadius, | ||
| 34 | const Throat::Shape shape) noexcept | ||
| 35 | { | ||
| 36 | using std::cos; | ||
| 37 | using std::sin; | ||
| 38 | |||
| 39 |
1/2✓ Branch 1 taken 959946 times.
✗ Branch 2 not taken.
|
959946 | const Scalar cornerHalfAngle = Throat::cornerHalfAngles<Scalar>(shape)[0]; |
| 40 | |||
| 41 | // check if snap-off is possible with such contact angle and corner half-angle | ||
| 42 | // if the criteriun is not fulfilled, return the lowest negative value to make sure that snap-off will not happen | ||
| 43 |
1/2✓ Branch 0 taken 959946 times.
✗ Branch 1 not taken.
|
959946 | if (cornerHalfAngle + contactAngle >= 0.5*M_PI) |
| 44 | return std::numeric_limits<Scalar>::lowest(); | ||
| 45 | |||
| 46 | 959946 | const Scalar cosContactAngle = cos(contactAngle); | |
| 47 | 959946 | const Scalar sinContactAngle = sin(contactAngle); | |
| 48 | |||
| 49 | // Get tangent of corner half-angle | ||
| 50 | 959946 | const Scalar tanCornerHalfAngle = [&]{ switch (shape){ // optimized tan(corner half-angle) for certain shapes | |
| 51 | case Throat::Shape::equilateralTriangle: return 0.577; | ||
| 52 | 959946 | case Throat::Shape::square: return 1.0; | |
| 53 | ✗ | default: using std::tan; return tan(cornerHalfAngle); | |
| 54 |
1/6✓ Branch 0 taken 959946 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
959946 | }}(); |
| 55 | |||
| 56 | 959946 | return surfaceTension / inscribedRadius * (cosContactAngle - sinContactAngle * tanCornerHalfAngle); | |
| 57 | } | ||
| 58 | |||
| 59 | |||
| 60 | /*! \brief The snap-off capillary pressure of a pore throat | ||
| 61 | * It checks if the cross section shape of the throat is a regular or irregular shape and call | ||
| 62 | * the proper pc snap-off accordingly | ||
| 63 | */ | ||
| 64 | template <class Scalar> | ||
| 65 | 959946 | Scalar pcSnapoff(const Scalar surfaceTension, | |
| 66 | const Scalar contactAngle, | ||
| 67 | const Scalar inscribedRadius, | ||
| 68 | const Throat::Shape shape) | ||
| 69 | { | ||
| 70 |
1/2✓ Branch 1 taken 959946 times.
✗ Branch 2 not taken.
|
959946 | if (Throat::isRegularShape(shape)) // call the pc snap-off calculated for regular shapes(same corner angles and side length) |
| 71 | 959946 | return pcSnapoffRegularShape(surfaceTension, contactAngle, inscribedRadius, shape); | |
| 72 | else | ||
| 73 | ✗ | DUNE_THROW(Dune::NotImplemented, "Pc snap-off is not implemented for this irregular shape"); | |
| 74 | } | ||
| 75 | |||
| 76 | /*! \brief The entry capillary pressure of a pore throat. | ||
| 77 | * | ||
| 78 | * For details, see Eq. 11 in Rabbani et al., 2016: https://doi.org/10.1016/j.jcis.2016.03.053 | ||
| 79 | * or Eq A-7 in Oren et al., 1998: https://doi.org/10.2118/52052-PA | ||
| 80 | */ | ||
| 81 | template<class Scalar> | ||
| 82 | 960833 | Scalar pcEntry(const Scalar surfaceTension, | |
| 83 | const Scalar contactAngle, | ||
| 84 | const Scalar inscribedRadius, | ||
| 85 | const Scalar shapeFactor) noexcept | ||
| 86 | { | ||
| 87 | using std::sin; | ||
| 88 | using std::cos; | ||
| 89 | using std::sqrt; | ||
| 90 | |||
| 91 | 960833 | const Scalar cosContactAngle = cos(contactAngle); | |
| 92 | 960833 | const Scalar sinContactAngle = sin(contactAngle); | |
| 93 | |||
| 94 | 960833 | const Scalar D = M_PI - 3*contactAngle + 3*sinContactAngle*cosContactAngle - (cosContactAngle*cosContactAngle) / (4*shapeFactor); | |
| 95 | 960833 | const Scalar F = (1 + sqrt(1 + 4*shapeFactor*D / (cosContactAngle*cosContactAngle))) / (1 + 2*sqrt(M_PI*shapeFactor)); | |
| 96 | 960833 | return surfaceTension / inscribedRadius * cosContactAngle * (1 + 2*sqrt(M_PI*shapeFactor)) * F; | |
| 97 | } | ||
| 98 | |||
| 99 | } // namespace Dumux::PoreNetwork::ThresholdCapillaryPressures | ||
| 100 | |||
| 101 | #endif | ||
| 102 |