movingWallVelocityFvPatchVectorField.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-2022 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 
27 #include "volFields.H"
28 #include "surfaceFields.H"
29 #include "fvcMeshPhi.H"
31 
32 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
33 
36 (
37  const fvPatch& p,
39 )
40 :
41  fixedValueFvPatchVectorField(p, iF)
42 {}
43 
44 
47 (
48  const fvPatch& p,
50  const dictionary& dict
51 )
52 :
53  fixedValueFvPatchVectorField(p, iF, dict)
54 {}
55 
56 
59 (
61  const fvPatch& p,
63  const fvPatchFieldMapper& mapper
64 )
65 :
66  fixedValueFvPatchVectorField(ptf, p, iF, mapper)
67 {}
68 
69 
72 (
75 )
76 :
77  fixedValueFvPatchVectorField(mwvpvf, iF)
78 {}
79 
80 
81 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
82 
84 {
85  if (updated())
86  {
87  return;
88  }
89 
90  const fvMesh& mesh = internalField().mesh();
91 
92  if (mesh.moving())
93  {
94  const fvPatch& p = patch();
95  const polyPatch& pp = p.patch();
96  const pointField& oldPoints = mesh.oldPoints();
97 
98  vectorField oldFc(pp.size());
99 
100  forAll(oldFc, i)
101  {
102  oldFc[i] = pp[i].centre(oldPoints);
103  }
104 
105  const scalar deltaT = mesh.time().deltaTValue();
106 
107  const vectorField Up((pp.faceCentres() - oldFc)/deltaT);
108 
109  const volVectorField& U =
110  static_cast<const volVectorField&>(internalField());
111 
112  const scalarField phip(fvc::meshPhi(U, p.index()));
113 
114  const vectorField n(p.nf());
115  const scalarField& magSf = p.magSf();
116  tmp<scalarField> Un = phip/(magSf + vSmall);
117 
118  vectorField::operator=(Up + n*(Un - (n & Up)));
119  }
120 
121  fixedValueFvPatchVectorField::updateCoeffs();
122 }
123 
124 
126 {
128  writeEntry(os, "value", *this);
129 }
130 
131 
132 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
133 
134 namespace Foam
135 {
137  (
140  );
141 }
142 
143 // ************************************************************************* //
Foam::surfaceFields.
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
This boundary condition provides a velocity condition for cases with moving walls.
bool moving() const
Is mesh moving.
Definition: polyMesh.H:525
U
Definition: pEqn.H:72
const scalarField & magSf() const
Return face area magnitudes.
Definition: fvPatch.C:142
label index() const
Return the index of this patch in the fvBoundaryMesh.
Definition: fvPatch.H:175
A list of keyword definitions, which are a keyword followed by any number of values (e...
Definition: dictionary.H:156
tmp< vectorField > nf() const
Return face normals.
Definition: fvPatch.C:130
A finiteVolume patch using a polyPatch and a fvBoundaryMesh.
Definition: fvPatch.H:63
Abstract base class with a fat-interface to all derived classes covering all possible ways in which t...
Definition: fvPatchField.H:66
const Time & time() const
Return the top-level database.
Definition: fvMesh.H:372
virtual void write(Ostream &) const
Write.
Definition: fvPatchField.C:243
fvMesh & mesh
Macros for easy insertion into run-time selection tables.
const polyPatch & patch() const
Return the polyPatch.
Definition: fvPatch.H:139
Foam::fvPatchFieldMapper.
scalar deltaTValue() const
Return time step value.
Definition: TimeStateI.H:34
virtual const pointField & oldPoints() const
Return old points for mesh motion.
Definition: polyMesh.C:1249
An Ostream is an abstract base class for all output systems (streams, files, token lists...
Definition: Ostream.H:54
void operator=(const Field< vector > &)
Definition: Field.C:526
void writeEntry(Ostream &os, const HashTable< T, Key, Hash > &ht)
Definition: HashTableIO.C:96
Calculate the mesh motion flux and convert fluxes from absolute to relative and back.
Mesh data needed to do the Finite Volume discretisation.
Definition: fvMesh.H:95
Field with dimensions and associated with geometry type GeoMesh which is used to size the field and a...
tmp< surfaceScalarField > meshPhi(const volVectorField &U)
Definition: fvcMeshPhi.C:34
label n
virtual void updateCoeffs()
Update the coefficients associated with the patch field.
movingWallVelocityFvPatchVectorField(const fvPatch &, const DimensionedField< vector, volMesh > &)
Construct from patch and internal field.
A class for managing temporary objects.
Definition: PtrList.H:53
A patch is a list of labels that address the faces in the global face list.
Definition: polyPatch.H:66
makePatchTypeField(fvPatchScalarField, atmBoundaryLayerInletEpsilonFvPatchScalarField)
Namespace for OpenFOAM.
const vectorField::subField faceCentres() const
Return face centres.
Definition: polyPatch.C:303