diff --git a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUSpaldingWallFunction/nutUSpaldingWallFunctionFvPatchScalarField.C b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUSpaldingWallFunction/nutUSpaldingWallFunctionFvPatchScalarField.C index 6cb1d4b832..fb5d2cba2b 100644 --- a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUSpaldingWallFunction/nutUSpaldingWallFunctionFvPatchScalarField.C +++ b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUSpaldingWallFunction/nutUSpaldingWallFunctionFvPatchScalarField.C @@ -53,11 +53,39 @@ tmp nutUSpaldingWallFunctionFvPatchScalarField::calcNut() const const tmp tnuw = turbModel.nu(patchi); const scalarField& nuw = tnuw(); - return max + + // Calculate new viscosity + tmp tnutw ( - scalar(0), - sqr(calcUTau(magGradU))/(magGradU + ROOTVSMALL) - nuw + max + ( + scalar(0), + sqr(calcUTau(magGradU))/(magGradU + ROOTVSMALL) - nuw + ) ); + + if (tolerance_ != 0.01) + { + // User-specified tolerance. Find out if current nut already satisfies + // eqns. + + // Run ut for one iteration + scalarField err; + tmp UTau(calcUTau(magGradU, 1, err)); + + // Preserve nutw if the initial error (err) already lower than the + // tolerance. + + scalarField& nutw = tnutw.ref(); + forAll(err, facei) + { + if (err[facei] < tolerance_) + { + nutw[facei] = this->operator[](facei); + } + } + } + return tnutw; } @@ -65,6 +93,18 @@ tmp nutUSpaldingWallFunctionFvPatchScalarField::calcUTau ( const scalarField& magGradU ) const +{ + scalarField err; + return calcUTau(magGradU, maxIter_, err); +} + + +tmp nutUSpaldingWallFunctionFvPatchScalarField::calcUTau +( + const scalarField& magGradU, + const label maxIter, + scalarField& err +) const { const label patchi = patch().index(); @@ -89,6 +129,9 @@ tmp nutUSpaldingWallFunctionFvPatchScalarField::calcUTau tmp tuTau(new scalarField(patch().size(), Zero)); scalarField& uTau = tuTau.ref(); + err.setSize(uTau.size()); + err = 0.0; + forAll(uTau, facei) { scalar ut = sqrt((nutw[facei] + nuw[facei])*magGradU[facei]); @@ -98,7 +141,6 @@ tmp nutUSpaldingWallFunctionFvPatchScalarField::calcUTau if (ut > ROOTVSMALL) { int iter = 0; - scalar err = GREAT; do { @@ -116,12 +158,17 @@ tmp nutUSpaldingWallFunctionFvPatchScalarField::calcUTau + 1/E_*kUu*fkUu/ut; scalar uTauNew = ut + f/df; - err = mag((ut - uTauNew)/ut); + err[facei] = mag((ut - uTauNew)/ut); ut = uTauNew; //iterations_++; - } while (ut > ROOTVSMALL && err > tolerance_ && ++iter < maxIter_); + } while + ( + ut > ROOTVSMALL + && err[facei] > tolerance_ + && ++iter < maxIter + ); uTau[facei] = max(0.0, ut); diff --git a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUSpaldingWallFunction/nutUSpaldingWallFunctionFvPatchScalarField.H b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUSpaldingWallFunction/nutUSpaldingWallFunctionFvPatchScalarField.H index d330b1d14d..d28a799693 100644 --- a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUSpaldingWallFunction/nutUSpaldingWallFunctionFvPatchScalarField.H +++ b/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutUSpaldingWallFunction/nutUSpaldingWallFunctionFvPatchScalarField.H @@ -71,9 +71,18 @@ Note Suffers from non-exact restart since correctNut() (called through turbulence->validate) returns a slightly different value every time it is called. This is since the seed for the Newton-Raphson iteration - uses the current value of *this (= nut). Can be avoided by seeding the - NR with e.g. the laminar viscosity or tightening the convergence tolerance - to e.g. 1e-7 and the max number of iterations to 100. + uses the current value of *this (= nut). + + This can be avoided by overriding the tolerance. This also switches on + a pre-detection whether the current nut already satisfies the turbulence + conditions and if so does not change it at all. This means that the nut + only changes if it 'has really changed'. This probably should be used with + a tight tolerance, e.g. + + maxIter 100; + tolerance 1e-7; + + to make sure to kick every iteration. SourceFiles nutUSpaldingWallFunctionFvPatchScalarField.C @@ -123,6 +132,15 @@ protected: //- Calculate the friction velocity virtual tmp calcUTau(const scalarField& magGradU) const; + //- Calculate the friction velocity and number of iterations for + // convergence + virtual tmp calcUTau + ( + const scalarField& magGradU, + const label maxIter, + scalarField& err + ) const; + //- Write local wall function variables virtual void writeLocalEntries(Ostream&) const;