fixedShearStressFvPatchVectorField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2020 OpenFOAM Foundation
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
28 #include "fvPatchFieldMapper.H"
29 #include "volFields.H"
30 #include "surfaceFields.H"
31 #include "momentumTransportModel.H"
32 
33 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
34 
36 (
37  const fvPatch& p,
39 )
40 :
41  fixedValueFvPatchVectorField(p, iF),
42  tau0_(Zero)
43 {}
44 
45 
47 (
48  const fvPatch& p,
50  const dictionary& dict
51 )
52 :
53  fixedValueFvPatchVectorField(p, iF, dict, false),
54  tau0_(dict.lookupOrDefault<vector>("tau", Zero))
55 {
56  fvPatchField<vector>::operator=(patchInternalField());
57 }
58 
59 
61 (
63  const fvPatch& p,
65  const fvPatchFieldMapper& mapper
66 )
67 :
68  fixedValueFvPatchVectorField(ptf, p, iF, mapper),
69  tau0_(ptf.tau0_)
70 {}
71 
72 
74 (
77 )
78 :
79  fixedValueFvPatchVectorField(ptf, iF),
80  tau0_(ptf.tau0_)
81 {}
82 
83 
84 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
85 
87 {
88  if (updated())
89  {
90  return;
91  }
92 
93  const momentumTransportModel& turbModel =
94  db().lookupObject<momentumTransportModel>
95  (
97  (
98  momentumTransportModel::typeName,
99  internalField().group()
100  )
101  );
102 
103  scalarField nuEff(turbModel.nuEff(patch().index()));
104 
105  const vectorField Uc(patchInternalField());
106 
107  vector tauHat = tau0_/(mag(tau0_) + rootVSmall);
108 
109  const scalarField& ry = patch().deltaCoeffs();
110 
111  operator==(tauHat*(tauHat & (tau0_*(1.0/(ry*nuEff)) + Uc)));
112 
113  fixedValueFvPatchVectorField::updateCoeffs();
114 }
115 
116 
118 {
120  writeEntry(os, "tau", tau0_);
121  writeEntry(os, "value", *this);
122 }
123 
124 
125 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
126 
127 namespace Foam
128 {
130  (
133  );
134 }
135 
136 // ************************************************************************* //
Foam::surfaceFields.
const char *const group
Group name for atomic constants.
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
Set a constant shear stress as tau0 = -nuEff dU/dn.
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:62
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:260
Macros for easy insertion into run-time selection tables.
tmp< fvMatrix< Type > > operator==(const fvMatrix< Type > &, const fvMatrix< Type > &)
static word groupName(Name name, const word &group)
Foam::fvPatchFieldMapper.
static const zero Zero
Definition: zero.H:97
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
Abstract base class for turbulence models (RAS, LES and laminar).
T lookupOrDefault(const word &, const T &, bool recursive=false, bool patternMatch=true) const
Find and return a T,.
fixedShearStressFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
dimensioned< scalar > mag(const dimensioned< Type > &)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
virtual tmp< volScalarField > nuEff() const =0
Return the effective viscosity.
Namespace for OpenFOAM.