FreeFOAM The Cross-Platform CFD Toolkit
waveTransmissiveFvPatchField.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd.
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 
29 #include <finiteVolume/volFields.H>
33 
34 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 
39 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
40 
41 template<class Type>
43 (
44  const fvPatch& p,
46 )
47 :
49  psiName_("psi"),
50  gamma_(0.0)
51 {}
52 
53 
54 template<class Type>
56 (
58  const fvPatch& p,
60  const fvPatchFieldMapper& mapper
61 )
62 :
63  advectiveFvPatchField<Type>(ptf, p, iF, mapper),
64  psiName_(ptf.psiName_),
65  gamma_(ptf.gamma_)
66 {}
67 
68 
69 template<class Type>
71 (
72  const fvPatch& p,
74  const dictionary& dict
75 )
76 :
77  advectiveFvPatchField<Type>(p, iF, dict),
78  psiName_(dict.lookupOrDefault<word>("psi", "psi")),
79  gamma_(readScalar(dict.lookup("gamma")))
80 {}
81 
82 
83 template<class Type>
85 (
86  const waveTransmissiveFvPatchField& ptpsf
87 )
88 :
90  psiName_(ptpsf.psiName_),
91  gamma_(ptpsf.gamma_)
92 {}
93 
94 
95 template<class Type>
97 (
98  const waveTransmissiveFvPatchField& ptpsf,
100 )
101 :
102  advectiveFvPatchField<Type>(ptpsf, iF),
103  psiName_(ptpsf.psiName_),
104  gamma_(ptpsf.gamma_)
105 {}
106 
107 
108 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
109 
110 template<class Type>
112 {
113  // Lookup the velocity and compressibility of the patch
114  const fvPatchField<scalar>& psip = this->patch().lookupPatchField
115  (
116  psiName_,
117  reinterpret_cast<const volScalarField*>(0),
118  reinterpret_cast<const scalar*>(0)
119  );
120 
121  const surfaceScalarField& phi =
122  this->db().objectRegistry::lookupObject<surfaceScalarField>
123  (this->phiName_);
124 
125  fvsPatchField<scalar> phip = this->patch().lookupPatchField
126  (
127  this->phiName_,
128  reinterpret_cast<const surfaceScalarField*>(0),
129  reinterpret_cast<const scalar*>(0)
130  );
131 
133  {
134  const fvPatchScalarField& rhop = this->patch().lookupPatchField
135  (
136  this->rhoName_,
137  reinterpret_cast<const volScalarField*>(0),
138  reinterpret_cast<const scalar*>(0)
139  );
140 
141  phip /= rhop;
142  }
143 
144  // Calculate the speed of the field wave w
145  // by summing the component of the velocity normal to the boundary
146  // and the speed of sound (sqrt(gamma_/psi)).
147  return phip/this->patch().magSf() + sqrt(gamma_/psip);
148 }
149 
150 
151 template<class Type>
153 {
155 
156  if (this->phiName_ != "phi")
157  {
158  os.writeKeyword("phi") << this->phiName_ << token::END_STATEMENT << nl;
159  }
160  if (this->rhoName_ != "rho")
161  {
162  os.writeKeyword("rho") << this->rhoName_ << token::END_STATEMENT << nl;
163  }
164  if (psiName_ != "psi")
165  {
166  os.writeKeyword("psi") << psiName_ << token::END_STATEMENT << nl;
167  }
168  os.writeKeyword("gamma") << gamma_ << token::END_STATEMENT << nl;
169 
170  if (this->lInf_ > SMALL)
171  {
172  os.writeKeyword("fieldInf") << this->fieldInf_
173  << token::END_STATEMENT << nl;
174  os.writeKeyword("lInf") << this->lInf_
175  << token::END_STATEMENT << nl;
176  }
177 
178  this->writeEntry("value", os);
179 }
180 
181 
182 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
183 
184 } // End namespace Foam
185 
186 // ************************ vim: set sw=4 sts=4 et: ************************ //