00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 #include "smoluchowskiJumpTFvPatchScalarField.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 #include <finiteVolume/fvPatchFieldMapper.H>
00029 #include <finiteVolume/volFields.H>
00030 #include <OpenFOAM/mathematicalConstants.H>
00031
00032
00033
00034 namespace Foam
00035 {
00036
00037
00038
00039 smoluchowskiJumpTFvPatchScalarField::smoluchowskiJumpTFvPatchScalarField
00040 (
00041 const fvPatch& p,
00042 const DimensionedField<scalar, volMesh>& iF
00043 )
00044 :
00045 mixedFvPatchScalarField(p, iF),
00046 accommodationCoeff_(1.0),
00047 Twall_(p.size(), 0.0),
00048 gamma_(1.4)
00049 {
00050 refValue() = 0.0;
00051 refGrad() = 0.0;
00052 valueFraction() = 0.0;
00053 }
00054
00055 smoluchowskiJumpTFvPatchScalarField::smoluchowskiJumpTFvPatchScalarField
00056 (
00057 const smoluchowskiJumpTFvPatchScalarField& ptf,
00058 const fvPatch& p,
00059 const DimensionedField<scalar, volMesh>& iF,
00060 const fvPatchFieldMapper& mapper
00061 )
00062 :
00063 mixedFvPatchScalarField(ptf, p, iF, mapper),
00064 accommodationCoeff_(ptf.accommodationCoeff_),
00065 Twall_(ptf.Twall_),
00066 gamma_(ptf.gamma_)
00067 {}
00068
00069
00070 smoluchowskiJumpTFvPatchScalarField::smoluchowskiJumpTFvPatchScalarField
00071 (
00072 const fvPatch& p,
00073 const DimensionedField<scalar, volMesh>& iF,
00074 const dictionary& dict
00075 )
00076 :
00077 mixedFvPatchScalarField(p, iF),
00078 accommodationCoeff_(readScalar(dict.lookup("accommodationCoeff"))),
00079 Twall_("Twall", dict, p.size()),
00080 gamma_(dict.lookupOrDefault<scalar>("gamma", 1.4))
00081 {
00082 if
00083 (
00084 mag(accommodationCoeff_) < SMALL
00085 || mag(accommodationCoeff_) > 2.0
00086 )
00087 {
00088 FatalIOErrorIn
00089 (
00090 "smoluchowskiJumpTFvPatchScalarField::"
00091 "smoluchowskiJumpTFvPatchScalarField"
00092 "("
00093 " const fvPatch&,"
00094 " const DimensionedField<scalar, volMesh>&,"
00095 " const dictionary&"
00096 ")",
00097 dict
00098 ) << "unphysical accommodationCoeff specified"
00099 << "(0 < accommodationCoeff <= 1)" << endl
00100 << exit(FatalError);
00101 }
00102
00103 if (dict.found("value"))
00104 {
00105 fvPatchField<scalar>::operator=
00106 (
00107 scalarField("value", dict, p.size())
00108 );
00109 }
00110 else
00111 {
00112 fvPatchField<scalar>::operator=(patchInternalField());
00113 }
00114
00115 refValue() = *this;
00116 refGrad() = 0.0;
00117 valueFraction() = 0.0;
00118 }
00119
00120
00121 smoluchowskiJumpTFvPatchScalarField::smoluchowskiJumpTFvPatchScalarField
00122 (
00123 const smoluchowskiJumpTFvPatchScalarField& ptpsf,
00124 const DimensionedField<scalar, volMesh>& iF
00125 )
00126 :
00127 mixedFvPatchScalarField(ptpsf, iF),
00128 accommodationCoeff_(ptpsf.accommodationCoeff_),
00129 Twall_(ptpsf.Twall_),
00130 gamma_(ptpsf.gamma_)
00131 {}
00132
00133
00134
00135
00136
00137 void smoluchowskiJumpTFvPatchScalarField::autoMap
00138 (
00139 const fvPatchFieldMapper& m
00140 )
00141 {
00142 mixedFvPatchScalarField::autoMap(m);
00143 }
00144
00145
00146
00147 void smoluchowskiJumpTFvPatchScalarField::rmap
00148 (
00149 const fvPatchField<scalar>& ptf,
00150 const labelList& addr
00151 )
00152 {
00153 mixedFvPatchField<scalar>::rmap(ptf, addr);
00154 }
00155
00156
00157
00158 void smoluchowskiJumpTFvPatchScalarField::updateCoeffs()
00159 {
00160 if (updated())
00161 {
00162 return;
00163 }
00164
00165 const fvPatchScalarField& pmu =
00166 patch().lookupPatchField<volScalarField, scalar>("mu");
00167 const fvPatchScalarField& prho =
00168 patch().lookupPatchField<volScalarField, scalar>("rho");
00169 const fvPatchField<scalar>& ppsi =
00170 patch().lookupPatchField<volScalarField, scalar>("psi");
00171 const fvPatchVectorField& pU =
00172 patch().lookupPatchField<volVectorField, vector>("U");
00173
00174
00175 const dictionary& thermophysicalProperties =
00176 db().lookupObject<IOdictionary>("thermophysicalProperties");
00177 dimensionedScalar Pr = dimensionedScalar("Pr", dimless, 1.0);
00178 if (thermophysicalProperties.found("Pr"))
00179 {
00180 Pr = thermophysicalProperties.lookup("Pr");
00181 }
00182
00183 Field<scalar> C2 = pmu/prho
00184 *sqrt(ppsi*mathematicalConstant::pi/2.0)
00185 *2.0*gamma_/Pr.value()/(gamma_ + 1.0)
00186 *(2.0 - accommodationCoeff_)/accommodationCoeff_;
00187
00188 Field<scalar> aCoeff = prho.snGrad() - prho/C2;
00189 Field<scalar> KEbyRho = 0.5*magSqr(pU);
00190
00191 valueFraction() = (1.0/(1.0 + patch().deltaCoeffs()*C2));
00192 refValue() = Twall_;
00193 refGrad() = 0.0;
00194
00195 mixedFvPatchScalarField::updateCoeffs();
00196 }
00197
00198
00199
00200 void smoluchowskiJumpTFvPatchScalarField::write(Ostream& os) const
00201 {
00202 fvPatchScalarField::write(os);
00203 os.writeKeyword("accommodationCoeff")
00204 << accommodationCoeff_ << token::END_STATEMENT << nl;
00205 Twall_.writeEntry("Twall", os);
00206 os.writeKeyword("gamma")
00207 << gamma_ << token::END_STATEMENT << nl;
00208 writeEntry("value", os);
00209 }
00210
00211
00212
00213
00214 makePatchTypeField(fvPatchScalarField, smoluchowskiJumpTFvPatchScalarField);
00215
00216
00217
00218 }
00219
00220