[UnitTest] Added a test for WaterVapourDensity

Tests/MaterialLib/TestMPLWaterVapourDensity.cpp

+/**
+ * \file
+ * \copyright
+ * Copyright (c) 2012-2021, OpenGeoSys Community (http://www.opengeosys.org)
+ *            Distributed under a Modified BSD License.
+ *              See accompanying file LICENSE.txt or
+ *              http://www.opengeosys.org/project/license
+ *
+ * Created on March 4, 2021, 5:23 PM
+ */
+#include <gtest/gtest.h>
+
+#include <array>
+#include <cmath>
+#include <limits>
+
+#include "MaterialLib/MPL/Properties/Density/CreateWaterVapourDensity.h"
+#include "TestMPL.h"
+
+TEST(MaterialPropertyLib, WaterVapourDensity)
+{
+    char const xml[] =
+        "<property>"
+        "   <name>vapour_density</name>"
+        "   <type>WaterVapourDensity</type>"
+        "</property>";
+
+    std::unique_ptr<MaterialPropertyLib::Property> const property_ptr =
+        Tests::createTestProperty(
+            xml, MaterialPropertyLib::createWaterVapourDensity);
+    MaterialPropertyLib::Property const& property = *property_ptr;
+
+    double const T = 390.0;
+    double const p = -1.e+6;
+    double const rho_w = 1000.0;
+
+    MaterialPropertyLib::VariableArray variable_array;
+    ParameterLib::SpatialPosition const pos;
+    double const t = std::numeric_limits<double>::quiet_NaN();
+    double const dt = std::numeric_limits<double>::quiet_NaN();
+    variable_array[static_cast<int>(
+        MaterialPropertyLib::Variable::temperature)] = T;
+    variable_array[static_cast<int>(
+        MaterialPropertyLib::Variable::phase_pressure)] = p;
+    variable_array[static_cast<int>(MaterialPropertyLib::Variable::density)] =
+        rho_w;
+
+    // The derivative of the water vapour with respect of temperature
+    {
+        std::array const temperatures = {273.0, 293.0, 393.0, 420.0, 500.0};
+        std::array const rho_vw_expected = {4.875989e-03, 1.692871e-02,
+                                            1.276865, 2.882635, 1.920386e+01};
+
+        for (std::size_t i = 0; i < temperatures.size(); ++i)
+        {
+            double const T_i = temperatures[i];
+            variable_array[static_cast<int>(
+                MaterialPropertyLib::Variable::temperature)] = T_i;
+
+            double const rho_vw =
+                property.template value<double>(variable_array, pos, t, dt);
+
+            ASSERT_LE(std::fabs(rho_vw_expected[i] - rho_vw), 5e-6)
+                << "for expected water vapour density " << rho_vw_expected[i]
+                << " and for computed water vapour density " << rho_vw;
+
+            double const dT = 1.0e-4;
+            variable_array[static_cast<int>(
+                MaterialPropertyLib::Variable::temperature)] = T_i - dT;
+            double const rho_vw0 =
+                property.template value<double>(variable_array, pos, t, dt);
+
+            variable_array[static_cast<int>(
+                MaterialPropertyLib::Variable::temperature)] = T_i + dT;
+            double const rho_vw1 =
+                property.template value<double>(variable_array, pos, t, dt);
+
+            double const approximated_drho_wv_dT =
+                0.5 * (rho_vw1 - rho_vw0) / dT;
+
+            double const analytic_drho_wv_dT = property.template dValue<double>(
+                variable_array, MaterialPropertyLib::Variable::temperature, pos,
+                t, dt);
+
+            ASSERT_LE(std::fabs(approximated_drho_wv_dT - analytic_drho_wv_dT),
+                      1e-6)
+                << "for expected derivative of water vapour density with "
+                   "respect to temperature "
+                << approximated_drho_wv_dT
+                << " and for computed derivative of water vapour density "
+                   "with respect to temperature."
+                << analytic_drho_wv_dT;
+        }
+    }
+
+    // The derivative of the water vapour density with respect of pressure
+    {
+        std::array const pressures = {-1.e+7, -1.e+6, 1.e+5, 1.e+6,
+                                      2.e+6,  6.e+6,  1.e+7};
+        std::array const rho_vw_expected = {1.101824, 1.158320, 1.165421,
+                                            1.171263, 1.177789, 1.204257,
+                                            1.23132};
+
+        for (std::size_t i = 0; i < pressures.size(); ++i)
+        {
+            double const p_i = pressures[i];
+
+            variable_array[static_cast<int>(
+                MaterialPropertyLib::Variable::temperature)] = T;
+            variable_array[static_cast<int>(
+                MaterialPropertyLib::Variable::phase_pressure)] = p_i;
+
+            double const rho_vw =
+                property.template value<double>(variable_array, pos, t, dt);
+
+            ASSERT_LE(std::fabs(rho_vw_expected[i] - rho_vw), 5e-6)
+                << "for expected water vapour density " << rho_vw_expected[i]
+                << " and for computed water vapour density " << rho_vw;
+
+            double const dp = 1.0e-3;
+            variable_array[static_cast<int>(
+                MaterialPropertyLib::Variable::phase_pressure)] = p_i - dp;
+
+            double const rho_vw0 =
+                property.template value<double>(variable_array, pos, t, dt);
+
+            variable_array[static_cast<int>(
+                MaterialPropertyLib::Variable::phase_pressure)] = p_i + dp;
+            double const rho_vw1 =
+                property.template value<double>(variable_array, pos, t, dt);
+
+            double const approximated_drho_wv_dp =
+                0.5 * (rho_vw1 - rho_vw0) / dp;
+
+            double const analytic_drho_wv_dp = property.template dValue<double>(
+                variable_array, MaterialPropertyLib::Variable::phase_pressure,
+                pos, t, dt);
+
+            ASSERT_LE(std::fabs(approximated_drho_wv_dp - analytic_drho_wv_dp),
+                      1e-7)
+                << "for expected derivative of water vapour density with "
+                   "respect to pressure "
+                << approximated_drho_wv_dp
+                << " and for computed derivative of water vapour density "
+                   "with respect to pressure."
+                << analytic_drho_wv_dp;
+        }
+    }
+}