correcting some bugs in frictional pressure loss.

This commit is contained in:
Kai Bao 2017-10-03 14:57:41 +02:00
parent f6ad3669f8
commit 92abdc4f23
2 changed files with 11 additions and 16 deletions

View File

@ -71,8 +71,12 @@ namespace mswellhelpers
static double haalandFormular(const double re, const double diameter, const double roughness)
{
const double value = std::exp(10. / 9. * std::log(roughness / (3.7 * diameter)) );
return -3.6 * std::log10( 6.9 / re + value);
const double value = -3.6 * std::log10(6.9 / re + std::pow(roughness / (3.7 * diameter), 10. / 9.) );
// sqrt(1/f) should be non-positive
assert(value >= 0.0);
return 1. / (value * value);
}
@ -110,8 +114,7 @@ namespace mswellhelpers
// TODO: not sure whether mu, density and mass flow rate should be Evaluation
// only use its value for now.
// calculating the friction pressure loss
// l is the segment length
// area is the segment cross area
// diameter is the segment inner diameter
@ -124,7 +127,8 @@ namespace mswellhelpers
const ValueType& density, const ValueType& w, const ValueType& mu)
{
const double f = calculateFrictionFactor(area, diameter, w.value(), roughness, mu.value());
return f * l * w * w / (area * area * diameter * density);
// TODO: a factor of 2 needs to be here based on the dimensional analysis
return 2. * f * l * w * w / (area * area * diameter * density);
}

View File

@ -1360,12 +1360,10 @@ namespace Opm
segment_viscosities_[seg] = 0.;
// calculate the average viscosity
for (int p = 0; p < number_of_phases_; ++p) {
// const EvalWell phase_fraction = mix[p] / b[p] / volrat;
// segment_viscosities_[seg] += visc[p] * phase_fraction;
segment_viscosities_[seg] += visc[p] * mix[p];
const EvalWell phase_fraction = mix[p] / b[p] / volrat;
segment_viscosities_[seg] += visc[p] * phase_fraction;
}
// TODO: not handling solvent for now.
EvalWell density(0.0);
@ -1579,12 +1577,6 @@ namespace Opm
MultisegmentWell<TypeTag>::
assemblePressureEq(const int seg) const
{
// TODO: currently, we only handle the hydrostatic pressure difference.
// We need to add the friction pressure loss and also the acceleration pressure loss
// with the acceleration pressure loss, there will be inlets flow rates (maybe alos the oulet flow)
// not sure whether to handle them implicitly or explicitly
// TODO: we can try to handle them explicitly first, if it does not work, we can handle them
assert(seg != 0); // not top segment
// for top segment, the well control equation will be used.
@ -1595,7 +1587,6 @@ namespace Opm
pressure_equation -= getHydroPressureLoss(seg);
if (frictionalPressureLossConsidered()) {
// TODO: deciding the direction of friction later
pressure_equation -= getFrictionPressureLoss(seg);
}