From dd95be9f858c1c856d59ccbab837d0dc709608ea Mon Sep 17 00:00:00 2001
From: Wenqing Wang <wenqing.wang@ufz.de>
Date: Tue, 13 Dec 2016 13:26:14 +0100
Subject: [PATCH] [FL] Removed a static array buffer that can break multithread
environment.
---
.../Fluid/Viscosity/WaterViscosityIAPWS.cpp | 47 ++++++++++++++-----
1 file changed, 36 insertions(+), 11 deletions(-)
diff --git a/MaterialLib/Fluid/Viscosity/WaterViscosityIAPWS.cpp b/MaterialLib/Fluid/Viscosity/WaterViscosityIAPWS.cpp
index b8181b51b63..c6e1e885e03 100644
--- a/MaterialLib/Fluid/Viscosity/WaterViscosityIAPWS.cpp
+++ b/MaterialLib/Fluid/Viscosity/WaterViscosityIAPWS.cpp
@@ -12,6 +12,7 @@
#include "WaterViscosityIAPWS.h"
+#include <array>
#include <cmath>
namespace MaterialLib
@@ -29,10 +30,13 @@ static const double Hij[6][7] = {
static double computeBarMu0Factor(const double barT);
-static double series_factorT[6];
-static double series_factorRho[7];
-static void computeSeriesFactorForMu1(const double barT, const double bar_rho);
-static double computeBarMu1Factor();
+static std::array<double, 6> computeSeriesFactorTForMu1(const double barT,
+ const double bar_rho);
+static std::array<double, 7> computeSeriesFactorRhoForMu1(const double barT,
+ const double bar_rho);
+static double computeBarMu1Factor(
+ const std::array<double, 6>& series_factorT,
+ const std::array<double, 7>& series_factorRho);
static double computedBarMu_dbarT(const double barT, double bar_rho);
static double computedBarMu_dbarRho(const double barT, double bar_rho);
@@ -46,8 +50,10 @@ double WaterViscosityIAPWS::getValue(const ArrayType& var_vals) const
const double mu0 = 100. * std::sqrt(bar_T) / computeBarMu0Factor(bar_T);
- computeSeriesFactorForMu1(bar_T, bar_rho);
- const double mu1 = std::exp(bar_rho * computeBarMu1Factor());
+ const auto& series_factorT = computeSeriesFactorTForMu1(bar_T, bar_rho);
+ const auto& series_factorRho = computeSeriesFactorRhoForMu1(bar_T, bar_rho);
+ const double mu1 = std::exp(
+ bar_rho * computeBarMu1Factor(series_factorT, series_factorRho));
return mu0 * mu1 * _ref_mu;
}
@@ -65,7 +71,6 @@ double WaterViscosityIAPWS::getdValue(const ArrayType& var_vals,
case PropertyVariableType::T:
return _ref_mu * computedBarMu_dbarT(bar_T, bar_rho) / _ref_T;
case PropertyVariableType::rho:
- computeSeriesFactorForMu1(bar_T, bar_rho);
return _ref_mu * computedBarMu_dbarRho(bar_T, bar_rho) / _ref_rho;
default:
return 0.;
@@ -84,8 +89,10 @@ double computeBarMu0Factor(const double barT)
return sum_val;
}
-void computeSeriesFactorForMu1(const double barT, const double bar_rho)
+std::array<double, 6> computeSeriesFactorTForMu1(const double barT,
+ const double bar_rho)
{
+ std::array<double, 6> series_factorT;
series_factorT[0] = 1.;
const double barT_fac = 1 / barT - 1.0;
for (int i = 1; i < 6; i++)
@@ -93,14 +100,24 @@ void computeSeriesFactorForMu1(const double barT, const double bar_rho)
series_factorT[i] = series_factorT[i - 1] * barT_fac;
}
+ return series_factorT;
+}
+
+std::array<double, 7> computeSeriesFactorRhoForMu1(const double barT,
+ const double bar_rho)
+{
+ std::array<double, 7> series_factorRho;
series_factorRho[0] = 1.;
for (int i = 1; i < 7; i++)
{
series_factorRho[i] = series_factorRho[i - 1] * (bar_rho - 1.0);
}
+
+ return series_factorRho;
}
-double computeBarMu1Factor()
+double computeBarMu1Factor(const std::array<double, 6>& series_factorT,
+ const std::array<double, 7>& series_factorRho)
{
double sum_val = 0.;
for (int i = 0; i < 6; i++)
@@ -133,6 +150,9 @@ double computedBarMu_dbarT(const double barT, double bar_rho)
50. / (mu0_factor * sqrt_barT) -
100. * sqrt_barT * dmu0_factor_dbarT / (mu0_factor * mu0_factor);
+ const auto& series_factorT = computeSeriesFactorTForMu1(barT, bar_rho);
+ const auto& series_factorRho = computeSeriesFactorRhoForMu1(barT, bar_rho);
+
double dmu1_factor_dbarT = 0.0;
for (int i = 1; i < 6; i++)
{
@@ -145,7 +165,8 @@ double computedBarMu_dbarT(const double barT, double bar_rho)
sum_val_j / (barT * barT);
}
- const double mu1_factor = computeBarMu1Factor();
+ const double mu1_factor =
+ computeBarMu1Factor(series_factorT, series_factorRho);
const double dbar_mu1_dbarT =
bar_rho * std::exp(bar_rho * mu1_factor) * dmu1_factor_dbarT;
@@ -155,6 +176,9 @@ double computedBarMu_dbarT(const double barT, double bar_rho)
double computedBarMu_dbarRho(const double barT, double bar_rho)
{
+ const auto& series_factorT = computeSeriesFactorTForMu1(barT, bar_rho);
+ const auto& series_factorRho = computeSeriesFactorRhoForMu1(barT, bar_rho);
+
double dmu1_factor_dbar_rho = 0.0;
for (int i = 0; i < 6; i++)
{
@@ -169,7 +193,8 @@ double computedBarMu_dbarRho(const double barT, double bar_rho)
const double mu0 = 100. * std::sqrt(barT) / computeBarMu0Factor(barT);
- const double mu1_factor = computeBarMu1Factor();
+ const double mu1_factor =
+ computeBarMu1Factor(series_factorT, series_factorRho);
return mu0 * std::exp(bar_rho * mu1_factor) *
(mu1_factor + bar_rho * dmu1_factor_dbar_rho);
}
--
GitLab