FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

turbulentTemperatureCoupledBaffleFvPatchScalarField.C

Go to the documentation of this file.
00001 /*---------------------------------------------------------------------------*\
00002   =========                 |
00003   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
00004    \\    /   O peration     |
00005     \\  /    A nd           | Copyright (C) 1991-2010 OpenCFD Ltd.
00006      \\/     M anipulation  |
00007 -------------------------------------------------------------------------------
00008 License
00009     This file is part of OpenFOAM.
00010 
00011     OpenFOAM is free software: you can redistribute it and/or modify it
00012     under the terms of the GNU General Public License as published by
00013     the Free Software Foundation, either version 3 of the License, or
00014     (at your option) any later version.
00015 
00016     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
00017     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
00018     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
00019     for more details.
00020 
00021     You should have received a copy of the GNU General Public License
00022     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
00023 
00024 \*---------------------------------------------------------------------------*/
00025 
00026 #include "turbulentTemperatureCoupledBaffleFvPatchScalarField.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 #include <finiteVolume/fvPatchFieldMapper.H>
00029 #include <finiteVolume/volFields.H>
00030 #include <meshTools/directMappedPatchBase.H>
00031 #include "regionProperties.H"
00032 #include <basicThermophysicalModels/basicThermo.H>
00033 #include <compressibleRASModels/RASModel.H>
00034 #include <OpenFOAM/mapDistribute.H>
00035 
00036 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
00037 
00038 bool Foam::turbulentTemperatureCoupledBaffleFvPatchScalarField::interfaceOwner
00039 (
00040     const polyMesh& nbrRegion,
00041     const polyPatch& nbrPatch
00042 ) const
00043 {
00044     const fvMesh& myRegion = patch().boundaryMesh().mesh();
00045 
00046     if (nbrRegion.name() == myRegion.name())
00047     {
00048         return patch().index() < nbrPatch.index();
00049     }
00050     else
00051     {
00052         const regionProperties& props =
00053             myRegion.objectRegistry::parent().lookupObject<regionProperties>
00054             (
00055                 "regionProperties"
00056             );
00057 
00058         label myIndex = findIndex(props.fluidRegionNames(), myRegion.name());
00059         if (myIndex == -1)
00060         {
00061             label i = findIndex(props.solidRegionNames(), myRegion.name());
00062 
00063             if (i == -1)
00064             {
00065                 FatalErrorIn
00066                 (
00067                     "turbulentTemperatureCoupledBaffleFvPatchScalarField"
00068                     "::interfaceOwner(const polyMesh&"
00069                     ", const polyPatch&)const"
00070                 )   << "Cannot find region " << myRegion.name()
00071                     << " neither in fluids " << props.fluidRegionNames()
00072                     << " nor in solids " << props.solidRegionNames()
00073                     << exit(FatalError);
00074             }
00075             myIndex = props.fluidRegionNames().size() + i;
00076         }
00077         label nbrIndex = findIndex
00078         (
00079             props.fluidRegionNames(),
00080             nbrRegion.name()
00081         );
00082         if (nbrIndex == -1)
00083         {
00084             label i = findIndex(props.solidRegionNames(), nbrRegion.name());
00085 
00086             if (i == -1)
00087             {
00088                 FatalErrorIn
00089                 (
00090                     "coupleManager::interfaceOwner"
00091                     "(const polyMesh&, const polyPatch&) const"
00092                 )   << "Cannot find region " << nbrRegion.name()
00093                     << " neither in fluids " << props.fluidRegionNames()
00094                     << " nor in solids " << props.solidRegionNames()
00095                     << exit(FatalError);
00096             }
00097             nbrIndex = props.fluidRegionNames().size() + i;
00098         }
00099 
00100         return myIndex < nbrIndex;
00101     }
00102 }
00103 
00104 
00105 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00106 
00107 Foam::turbulentTemperatureCoupledBaffleFvPatchScalarField::
00108 turbulentTemperatureCoupledBaffleFvPatchScalarField
00109 (
00110     const fvPatch& p,
00111     const DimensionedField<scalar, volMesh>& iF
00112 )
00113 :
00114     fixedValueFvPatchScalarField(p, iF),
00115     neighbourFieldName_("undefined-neighbourFieldName"),
00116     KName_("undefined-K")
00117 {}
00118 
00119 
00120 Foam::turbulentTemperatureCoupledBaffleFvPatchScalarField::
00121 turbulentTemperatureCoupledBaffleFvPatchScalarField
00122 (
00123     const turbulentTemperatureCoupledBaffleFvPatchScalarField& ptf,
00124     const fvPatch& p,
00125     const DimensionedField<scalar, volMesh>& iF,
00126     const fvPatchFieldMapper& mapper
00127 )
00128 :
00129     fixedValueFvPatchScalarField(ptf, p, iF, mapper),
00130     neighbourFieldName_(ptf.neighbourFieldName_),
00131     KName_(ptf.KName_)
00132 {}
00133 
00134 
00135 Foam::turbulentTemperatureCoupledBaffleFvPatchScalarField::
00136 turbulentTemperatureCoupledBaffleFvPatchScalarField
00137 (
00138     const fvPatch& p,
00139     const DimensionedField<scalar, volMesh>& iF,
00140     const dictionary& dict
00141 )
00142 :
00143     fixedValueFvPatchScalarField(p, iF, dict),
00144     neighbourFieldName_(dict.lookup("neighbourFieldName")),
00145     KName_(dict.lookup("Kcond"))
00146 {
00147     if (!isA<directMappedPatchBase>(this->patch().patch()))
00148     {
00149         FatalErrorIn
00150         (
00151             "turbulentTemperatureCoupledBaffleFvPatchScalarField::"
00152             "turbulentTemperatureCoupledBaffleFvPatchScalarField\n"
00153             "(\n"
00154             "    const fvPatch& p,\n"
00155             "    const DimensionedField<scalar, volMesh>& iF,\n"
00156             "    const dictionary& dict\n"
00157             ")\n"
00158         )   << "\n    patch type '" << p.type()
00159             << "' not type '" << directMappedPatchBase::typeName << "'"
00160             << "\n    for patch " << p.name()
00161             << " of field " << dimensionedInternalField().name()
00162             << " in file " << dimensionedInternalField().objectPath()
00163             << exit(FatalError);
00164     }
00165 }
00166 
00167 
00168 Foam::turbulentTemperatureCoupledBaffleFvPatchScalarField::
00169 turbulentTemperatureCoupledBaffleFvPatchScalarField
00170 (
00171     const turbulentTemperatureCoupledBaffleFvPatchScalarField& wtcsf,
00172     const DimensionedField<scalar, volMesh>& iF
00173 )
00174 :
00175     fixedValueFvPatchScalarField(wtcsf, iF),
00176     neighbourFieldName_(wtcsf.neighbourFieldName_),
00177     KName_(wtcsf.KName_)
00178 {}
00179 
00180 
00181 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00182 
00183 Foam::tmp<Foam::scalarField>
00184 Foam::turbulentTemperatureCoupledBaffleFvPatchScalarField::K() const
00185 {
00186     const fvMesh& mesh = patch().boundaryMesh().mesh();
00187 
00188     if (KName_ == "none")
00189     {
00190         const compressible::RASModel& model =
00191             db().lookupObject<compressible::RASModel>("RASProperties");
00192 
00193         tmp<volScalarField> talpha = model.alphaEff();
00194 
00195         const basicThermo& thermo =
00196             db().lookupObject<basicThermo>("thermophysicalProperties");
00197 
00198         return
00199             talpha().boundaryField()[patch().index()]
00200            *thermo.Cp()().boundaryField()[patch().index()];
00201     }
00202     else if (mesh.objectRegistry::foundObject<volScalarField>(KName_))
00203     {
00204         return patch().lookupPatchField<volScalarField, scalar>(KName_);
00205     }
00206     else if (mesh.objectRegistry::foundObject<volSymmTensorField>(KName_))
00207     {
00208         const symmTensorField& KWall =
00209             patch().lookupPatchField<volSymmTensorField, scalar>(KName_);
00210 
00211         vectorField n = patch().nf();
00212 
00213         return n & KWall & n;
00214     }
00215     else
00216     {
00217         FatalErrorIn
00218         (
00219             "turbulentTemperatureCoupledBaffleFvPatchScalarField::K() const"
00220         )   << "Did not find field " << KName_
00221             << " on mesh " << mesh.name() << " patch " << patch().name()
00222             << endl
00223             << "Please set 'K' to 'none', a valid volScalarField"
00224             << " or a valid volSymmTensorField." << exit(FatalError);
00225 
00226         return scalarField(0);
00227     }
00228 }
00229 
00230 
00231 void Foam::turbulentTemperatureCoupledBaffleFvPatchScalarField::updateCoeffs()
00232 {
00233     if (updated())
00234     {
00235         return;
00236     }
00237 
00238     // Get the coupling information from the directMappedPatchBase
00239     const directMappedPatchBase& mpp = refCast<const directMappedPatchBase>
00240     (
00241         patch().patch()
00242     );
00243     const polyMesh& nbrMesh = mpp.sampleMesh();
00244     const fvPatch& nbrPatch = refCast<const fvMesh>
00245     (
00246         nbrMesh
00247     ).boundary()[mpp.samplePolyPatch().index()];
00248 
00249     // Force recalculation of mapping and schedule
00250     const mapDistribute& distMap = mpp.map();
00251     (void)distMap.schedule();
00252 
00253     tmp<scalarField> intFld = patchInternalField();
00254 
00255     if (interfaceOwner(nbrMesh, nbrPatch.patch()))
00256     {
00257         // Note: other side information could be cached - it only needs
00258         // to be updated the first time round the iteration (i.e. when
00259         // switching regions) but unfortunately we don't have this information.
00260 
00261 
00262         // Calculate the temperature by harmonic averaging
00263         // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
00264 
00265         const turbulentTemperatureCoupledBaffleFvPatchScalarField& nbrField =
00266         refCast<const turbulentTemperatureCoupledBaffleFvPatchScalarField>
00267         (
00268             nbrPatch.lookupPatchField<volScalarField, scalar>
00269             (
00270                 neighbourFieldName_
00271             )
00272         );
00273 
00274         // Swap to obtain full local values of neighbour internal field
00275         scalarField nbrIntFld = nbrField.patchInternalField();
00276         mapDistribute::distribute
00277         (
00278             Pstream::defaultCommsType,
00279             distMap.schedule(),
00280             distMap.constructSize(),
00281             distMap.subMap(),           // what to send
00282             distMap.constructMap(),     // what to receive
00283             nbrIntFld
00284         );
00285 
00286         // Swap to obtain full local values of neighbour K*delta
00287         scalarField nbrKDelta = nbrField.K()*nbrPatch.deltaCoeffs();
00288         mapDistribute::distribute
00289         (
00290             Pstream::defaultCommsType,
00291             distMap.schedule(),
00292             distMap.constructSize(),
00293             distMap.subMap(),           // what to send
00294             distMap.constructMap(),     // what to receive
00295             nbrKDelta
00296         );
00297 
00298         tmp<scalarField> myKDelta = K()*patch().deltaCoeffs();
00299 
00300         // Calculate common wall temperature. Reuse *this to store common value.
00301         scalarField Twall
00302         (
00303             (myKDelta()*intFld() + nbrKDelta*nbrIntFld)
00304           / (myKDelta() + nbrKDelta)
00305         );
00306         // Assign to me
00307         fvPatchScalarField::operator=(Twall);
00308         // Distribute back and assign to neighbour
00309         mapDistribute::distribute
00310         (
00311             Pstream::defaultCommsType,
00312             distMap.schedule(),
00313             nbrField.size(),
00314             distMap.constructMap(),     // reverse : what to send
00315             distMap.subMap(),
00316             Twall
00317         );
00318         const_cast<turbulentTemperatureCoupledBaffleFvPatchScalarField&>
00319         (
00320             nbrField
00321         ).fvPatchScalarField::operator=(Twall);
00322     }
00323 
00324     if (debug)
00325     {
00326         //tmp<scalarField> normalGradient =
00327         //    (*this-intFld())
00328         //  * patch().deltaCoeffs();
00329 
00330         scalar Q = gSum(K()*patch().magSf()*snGrad());
00331 
00332         Info<< patch().boundaryMesh().mesh().name() << ':'
00333             << patch().name() << ':'
00334             << this->dimensionedInternalField().name() << " <- "
00335             << nbrMesh.name() << ':'
00336             << nbrPatch.name() << ':'
00337             << this->dimensionedInternalField().name() << " :"
00338             << " heat[W]:" << Q
00339             << " walltemperature "
00340             << " min:" << gMin(*this)
00341             << " max:" << gMax(*this)
00342             << " avg:" << gAverage(*this)
00343             << endl;
00344     }
00345 
00346     fixedValueFvPatchScalarField::updateCoeffs();
00347 }
00348 
00349 
00350 void Foam::turbulentTemperatureCoupledBaffleFvPatchScalarField::write
00351 (
00352     Ostream& os
00353 ) const
00354 {
00355     fixedValueFvPatchScalarField::write(os);
00356     os.writeKeyword("neighbourFieldName")<< neighbourFieldName_
00357         << token::END_STATEMENT << nl;
00358     os.writeKeyword("Kcond") << KName_ << token::END_STATEMENT << nl;
00359 }
00360 
00361 
00362 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00363 
00364 namespace Foam
00365 {
00366 
00367 makePatchTypeField
00368 (
00369     fvPatchScalarField,
00370     turbulentTemperatureCoupledBaffleFvPatchScalarField
00371 );
00372 
00373 } // End namespace Foam
00374 
00375 // ************************************************************************* //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines