pressureInletOutletVelocityFvPatchVectorField.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-2023 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 "volFields.H"
29 #include "surfaceFields.H"
30 
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
36 (
37  const fvPatch& p,
39  const dictionary& dict
40 )
41 :
42  directionMixedFvPatchVectorField(p, iF),
43  phiName_(dict.lookupOrDefault<word>("phi", "phi"))
44 {
46 
47  if (dict.found("tangentialVelocity"))
48  {
49  tangentialVelocity_ =
50  Function1<vector>::New("tangentialVelocity", dict);
51  }
52 
53  refValue() = Zero;
54  refGrad() = Zero;
55  valueFraction() = Zero;
56 }
57 
58 
61 (
63  const fvPatch& p,
65  const fvPatchFieldMapper& mapper
66 )
67 :
68  directionMixedFvPatchVectorField(ptf, p, iF, mapper),
69  phiName_(ptf.phiName_),
70  tangentialVelocity_(ptf.tangentialVelocity_, false)
71 {}
72 
73 
76 (
79 )
80 :
81  directionMixedFvPatchVectorField(pivpvf, iF),
82  phiName_(pivpvf.phiName_),
83  tangentialVelocity_(pivpvf.tangentialVelocity_, false)
84 {}
85 
86 
87 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
88 
90 {
91  if (updated())
92  {
93  return;
94  }
95 
96  if (tangentialVelocity_.valid())
97  {
98  const scalar t = this->db().time().userTimeValue();
99  const vector tangentialVelocity = tangentialVelocity_->value(t);
100  const vectorField n(patch().nf());
101  refValue() = tangentialVelocity - n*(n & tangentialVelocity);
102  }
103 
104  const fvsPatchField<scalar>& phip =
105  patch().lookupPatchField<surfaceScalarField, scalar>(phiName_);
106 
107  valueFraction() = neg(phip)*(I - sqr(patch().nf()));
108 
109  directionMixedFvPatchVectorField::updateCoeffs();
111 }
112 
113 
115 (
116  Ostream& os
117 )
118 const
119 {
121  writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
122  if (tangentialVelocity_.valid())
123  {
124  writeEntry(os, tangentialVelocity_());
125  }
126  writeEntry(os, "value", *this);
127 }
128 
129 
130 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
131 
132 void Foam::pressureInletOutletVelocityFvPatchVectorField::operator=
133 (
134  const fvPatchField<vector>& pvf
135 )
136 {
137  tmp<vectorField> normalValue = transform(valueFraction(), refValue());
138  tmp<vectorField> transformGradValue = transform(I - valueFraction(), pvf);
139  fvPatchField<vector>::operator=(normalValue + transformGradValue);
140 }
141 
142 
143 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
144 
145 namespace Foam
146 {
148  (
151  );
152 }
153 
154 // ************************************************************************* //
label n
Macros for easy insertion into run-time selection tables.
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
static autoPtr< Function1< Type > > New(const word &name, const dictionary &dict)
Selector.
Definition: Function1New.C:32
An Ostream is an abstract base class for all output systems (streams, files, token lists,...
Definition: Ostream.H:57
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
Foam::fvPatchFieldMapper.
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:87
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:231
virtual void operator=(const UList< Type > &)
Definition: fvPatchField.C:251
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:64
An abstract base class with a fat-interface to all derived classes covering all possible ways in whic...
Definition: fvsPatchField.H:82
Velocity inlet/outlet boundary condition for patches where the pressure is specified in some manner,...
pressureInletOutletVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &, const dictionary &)
Construct from patch, internal field and dictionary.
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
A class for managing temporary objects.
Definition: tmp.H:55
A class for handling words, derived from string.
Definition: word.H:62
Namespace for OpenFOAM.
static const zero Zero
Definition: zero.H:97
dimensionedSymmTensor sqr(const dimensionedVector &dv)
SurfaceField< scalar > surfaceScalarField
static const Identity< scalar > I
Definition: Identity.H:93
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
dimensionSet transform(const dimensionSet &)
Definition: dimensionSet.C:483
Field< vector > vectorField
Specialisation of Field<T> for vector.
dimensionedScalar neg(const dimensionedScalar &ds)
void evaluate(GeometricField< Type, PatchField, GeoMesh > &result, const Function1< Type > &func, const GeometricField< Type, PatchField, GeoMesh > &x)
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
dictionary dict
volScalarField & p
Foam::surfaceFields.