thermoSingleLayer::q: Use a local "alpha" without hydrophilic/phobic adjustment
Resolves bug-report https://bugs.openfoam.org/view.php?id=2605
This commit is contained in:
@ -204,7 +204,7 @@ tmp<volScalarField> kinematicSingleLayer::pp()
|
||||
|
||||
void kinematicSingleLayer::correctAlpha()
|
||||
{
|
||||
alpha_ == pos0(delta_ - deltaSmall_);
|
||||
alpha_ == pos(delta_ - deltaSmall_);
|
||||
}
|
||||
|
||||
|
||||
|
||||
@ -268,15 +268,17 @@ void thermoSingleLayer::updateSubmodels()
|
||||
|
||||
tmp<fvScalarMatrix> thermoSingleLayer::q(volScalarField& hs) const
|
||||
{
|
||||
const volScalarField alpha(pos(delta_ - deltaSmall_));
|
||||
|
||||
return
|
||||
(
|
||||
// Heat-transfer to the primary region
|
||||
- fvm::Sp(htcs_->h()/Cp_, hs)
|
||||
+ htcs_->h()*(hs/Cp_ + alpha_*(TPrimary_ - T_))
|
||||
+ htcs_->h()*(hs/Cp_ + alpha*(TPrimary_ - T_))
|
||||
|
||||
// Heat-transfer to the wall
|
||||
- fvm::Sp(htcw_->h()/Cp_, hs)
|
||||
+ htcw_->h()*(hs/Cp_ + alpha_*(Tw_- T_))
|
||||
+ htcw_->h()*(hs/Cp_ + alpha*(Tw_- T_))
|
||||
);
|
||||
}
|
||||
|
||||
|
||||
Reference in New Issue
Block a user