1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 1991-2008 OpenCFD Ltd.
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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 the
13 Free Software Foundation; either version 2 of the License, or (at your
14 option) any later version.
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
21 You should have received a copy of the GNU General Public License
22 along with OpenFOAM; if not, write to the Free Software Foundation,
23 Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
27 \*---------------------------------------------------------------------------*/
29 #include "maxwellSlipUFvPatchVectorField.H"
30 #include "addToRunTimeSelectionTable.H"
31 #include "mathematicalConstants.H"
32 #include "fvPatchFieldMapper.H"
33 #include "volFields.H"
34 #include "surfaceFields.H"
37 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
42 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
44 maxwellSlipUFvPatchVectorField::maxwellSlipUFvPatchVectorField
47 const DimensionedField<vector, volMesh>& iF
50 mixedFixedValueSlipFvPatchVectorField(p, iF),
51 accommodationCoeff_(1.0),
52 Uwall_(p.size(), vector(0.0, 0.0, 0.0)),
58 maxwellSlipUFvPatchVectorField::maxwellSlipUFvPatchVectorField
60 const maxwellSlipUFvPatchVectorField& tdpvf,
62 const DimensionedField<vector, volMesh>& iF,
63 const fvPatchFieldMapper& mapper
66 mixedFixedValueSlipFvPatchVectorField(tdpvf, p, iF, mapper),
67 accommodationCoeff_(tdpvf.accommodationCoeff_),
69 thermalCreep_(tdpvf.thermalCreep_),
70 curvature_(tdpvf.curvature_)
74 maxwellSlipUFvPatchVectorField::maxwellSlipUFvPatchVectorField
77 const DimensionedField<vector, volMesh>& iF,
78 const dictionary& dict
81 mixedFixedValueSlipFvPatchVectorField(p, iF),
82 accommodationCoeff_(readScalar(dict.lookup("accommodationCoeff"))),
83 Uwall_("Uwall", dict, p.size()),
84 thermalCreep_(dict.lookupOrDefault("thermalCreep", true)),
85 curvature_(dict.lookupOrDefault("curvature", true))
89 mag(accommodationCoeff_) < SMALL
91 mag(accommodationCoeff_) > 2.0
96 "maxwellSlipUFvPatchScalarField::"
97 "maxwellSlipUFvPatchScalarField"
98 "(const fvPatch&, const scalarField&, const dictionary&)",
100 ) << "unphysical accommodationCoeff_ specified"
101 << "(0 < accommodationCoeff_ <= 1)" << endl
105 if (dict.found("value"))
107 fvPatchField<vector>::operator=
109 vectorField("value", dict, p.size())
114 mixedFixedValueSlipFvPatchVectorField::evaluate();
119 maxwellSlipUFvPatchVectorField::maxwellSlipUFvPatchVectorField
121 const maxwellSlipUFvPatchVectorField& tdpvf,
122 const DimensionedField<vector, volMesh>& iF
125 mixedFixedValueSlipFvPatchVectorField(tdpvf, iF),
126 accommodationCoeff_(tdpvf.accommodationCoeff_),
127 Uwall_(tdpvf.Uwall_),
128 thermalCreep_(tdpvf.thermalCreep_),
129 curvature_(tdpvf.curvature_)
133 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
135 // Update the coefficients associated with the patch field
136 void maxwellSlipUFvPatchVectorField::updateCoeffs()
143 const fvPatchScalarField& pmu =
144 patch().lookupPatchField<volScalarField, scalar>("mu");
145 const fvPatchScalarField& prho =
146 patch().lookupPatchField<volScalarField, scalar>("rho");
147 const fvPatchField<scalar>& ppsi =
148 patch().lookupPatchField<volScalarField, scalar>("psi");
150 Field<scalar> C1 = sqrt(ppsi*mathematicalConstant::pi/2.0)
151 *(2.0 - accommodationCoeff_)/accommodationCoeff_;
153 Field<scalar> pnu = pmu/prho;
154 valueFraction() = (1.0/(1.0 + patch().deltaCoeffs()*C1*pnu));
160 const volScalarField& vsfT =
161 this->db().objectRegistry::lookupObject<volScalarField>("T");
162 label patchi = this->patch().index();
163 const fvPatchScalarField& pT = vsfT.boundaryField()[patchi];
164 Field<vector> gradpT = fvc::grad(vsfT)().boundaryField()[patchi];
165 vectorField n = patch().nf();
167 refValue() -= 3.0*pnu/(4.0*pT)*transform(I - n*n, gradpT);
172 const fvPatchTensorField& ptauMC =
173 patch().lookupPatchField<volTensorField, tensor>("tauMC");
174 vectorField n = patch().nf();
176 refValue() -= C1/prho*transform(I - n*n, (n & ptauMC));
179 mixedFixedValueSlipFvPatchVectorField::updateCoeffs();
184 void maxwellSlipUFvPatchVectorField::write(Ostream& os) const
186 fvPatchVectorField::write(os);
187 os.writeKeyword("accommodationCoeff")
188 << accommodationCoeff_ << token::END_STATEMENT << nl;
189 Uwall_.writeEntry("Uwall", os);
190 os.writeKeyword("thermalCreep")
191 << thermalCreep_ << token::END_STATEMENT << nl;
192 os.writeKeyword("curvature") << curvature_ << token::END_STATEMENT << nl;
194 os.writeKeyword("refValue")
195 << refValue() << token::END_STATEMENT << nl;
196 os.writeKeyword("valueFraction")
197 << valueFraction() << token::END_STATEMENT << nl;
199 writeEntry("value", os);
203 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
205 makePatchTypeField(fvPatchVectorField, maxwellSlipUFvPatchVectorField);
207 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
209 } // End namespace Foam
211 // ************************************************************************* //