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

pointPatchField.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 "pointPatchField.H"
00027 #include <OpenFOAM/pointMesh.H>
00028 #include <OpenFOAM/dictionary.H>
00029 
00030 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00031 
00032 namespace Foam
00033 {
00034 
00035 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
00036 
00037 template<class Type>
00038 pointPatchField<Type>::pointPatchField
00039 (
00040     const pointPatch& p,
00041     const DimensionedField<Type, pointMesh>& iF
00042 )
00043 :
00044     patch_(p),
00045     internalField_(iF),
00046     updated_(false),
00047     patchType_(word::null)
00048 {}
00049 
00050 
00051 template<class Type>
00052 pointPatchField<Type>::pointPatchField
00053 (
00054     const pointPatch& p,
00055     const DimensionedField<Type, pointMesh>& iF,
00056     const dictionary& dict
00057 )
00058 :
00059     patch_(p),
00060     internalField_(iF),
00061     updated_(false),
00062     patchType_(dict.lookupOrDefault<word>("patchType", word::null))
00063 {}
00064 
00065 
00066 template<class Type>
00067 pointPatchField<Type>::pointPatchField
00068 (
00069     const pointPatchField<Type>& ptf
00070 )
00071 :
00072     patch_(ptf.patch_),
00073     internalField_(ptf.internalField_),
00074     updated_(false),
00075     patchType_(ptf.patchType_)
00076 {}
00077 
00078 
00079 template<class Type>
00080 pointPatchField<Type>::pointPatchField
00081 (
00082     const pointPatchField<Type>& ptf,
00083     const DimensionedField<Type, pointMesh>& iF
00084 )
00085 :
00086     patch_(ptf.patch_),
00087     internalField_(iF),
00088     updated_(false),
00089     patchType_(ptf.patchType_)
00090 {}
00091 
00092 
00093 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00094 
00095 template<class Type>
00096 const objectRegistry& pointPatchField<Type>::db() const
00097 {
00098     return patch_.boundaryMesh().mesh()();
00099 }
00100 
00101 
00102 template<class Type>
00103 void pointPatchField<Type>::write(Ostream& os) const
00104 {
00105     os.writeKeyword("type") << type() << token::END_STATEMENT << nl;
00106 
00107     if (patchType_.size())
00108     {
00109         os.writeKeyword("patchType") << patchType_
00110             << token::END_STATEMENT << nl;
00111     }
00112 }
00113 
00114 
00115 template<class Type>
00116 tmp<Field<Type> > pointPatchField<Type>::patchInternalField() const
00117 {
00118     return patchInternalField(internalField());
00119 }
00120 
00121 
00122 template<class Type>
00123 template<class Type1>
00124 tmp<Field<Type1> > pointPatchField<Type>::patchInternalField
00125 (
00126     const Field<Type1>& iF
00127 ) const
00128 {
00129     // Check size
00130     if (iF.size() != internalField().size())
00131     {
00132         FatalErrorIn
00133         (
00134             "tmp<Field<Type1> > pointPatchField<"
00135             "Type>::"
00136             "patchInternalField(const Field<Type1>& iF) const"
00137         )   << "given internal field does not correspond to the mesh. "
00138             << "Field size: " << iF.size()
00139             << " mesh size: " << internalField().size()
00140             << abort(FatalError);
00141     }
00142 
00143     // get addressing
00144     const labelList& meshPoints = patch().meshPoints();
00145 
00146     tmp<Field<Type1> > tvalues(new Field<Type1>(meshPoints.size()));
00147     Field<Type1>& values = tvalues();
00148 
00149     forAll (meshPoints, pointI)
00150     {
00151         values[pointI] = iF[meshPoints[pointI]];
00152     }
00153 
00154     return tvalues;
00155 }
00156 
00157 
00158 template<class Type>
00159 template<class Type1>
00160 void pointPatchField<Type>::addToInternalField
00161 (
00162     Field<Type1>& iF,
00163     const Field<Type1>& pF
00164 ) const
00165 {
00166     // Check size
00167     if (iF.size() != internalField().size())
00168     {
00169         FatalErrorIn
00170         (
00171             "void pointPatchField<Type>::"
00172             "addToInternalField("
00173             "Field<Type1>& iF, const Field<Type1>& iF) const"
00174         )   << "given internal field does not correspond to the mesh. "
00175             << "Field size: " << iF.size()
00176             << " mesh size: " << internalField().size()
00177             << abort(FatalError);
00178     }
00179 
00180     if (pF.size() != size())
00181     {
00182         FatalErrorIn
00183         (
00184             "void pointPatchField<Type>::"
00185             "addToInternalField("
00186             "Field<Type1>& iF, const Field<Type1>& iF) const"
00187         )   << "given patch field does not correspond to the mesh. "
00188             << "Field size: " << pF.size()
00189             << " mesh size: " << size()
00190             << abort(FatalError);
00191     }
00192 
00193     // Get the addressing
00194     const labelList& mp = patch().meshPoints();
00195 
00196     forAll (mp, pointI)
00197     {
00198         iF[mp[pointI]] += pF[pointI];
00199     }
00200 }
00201 
00202 
00203 template<class Type>
00204 template<class Type1>
00205 void pointPatchField<Type>::setInInternalField
00206 (
00207     Field<Type1>& iF,
00208     const Field<Type1>& pF,
00209     const labelList& meshPoints
00210 ) const
00211 {
00212     // Check size
00213     if (iF.size() != internalField().size())
00214     {
00215         FatalErrorIn
00216         (
00217             "void pointPatchField<Type>::"
00218             "setInInternalField("
00219             "Field<Type1>& iF, const Field<Type1>& iF) const"
00220         )   << "given internal field does not correspond to the mesh. "
00221             << "Field size: " << iF.size()
00222             << " mesh size: " << internalField().size()
00223             << abort(FatalError);
00224     }
00225 
00226     if (pF.size() != meshPoints.size())
00227     {
00228         FatalErrorIn
00229         (
00230             "void pointPatchField<Type>::"
00231             "setInInternalField("
00232             "Field<Type1>& iF, const Field<Type1>& iF) const"
00233         )   << "given patch field does not correspond to the meshPoints. "
00234             << "Field size: " << pF.size()
00235             << " meshPoints size: " << size()
00236             << abort(FatalError);
00237     }
00238 
00239     forAll (meshPoints, pointI)
00240     {
00241         iF[meshPoints[pointI]] = pF[pointI];
00242     }
00243 }
00244 
00245 
00246 template<class Type>
00247 template<class Type1>
00248 void pointPatchField<Type>::setInInternalField
00249 (
00250     Field<Type1>& iF,
00251     const Field<Type1>& pF
00252 ) const
00253 {
00254     setInInternalField(iF, pF, patch().meshPoints());
00255 }
00256 
00257 
00258 template<class Type>
00259 void pointPatchField<Type>::evaluate(const Pstream::commsTypes)
00260 {
00261     if (!updated_)
00262     {
00263         updateCoeffs();
00264     }
00265 
00266     updated_ = false;
00267 }
00268 
00269 
00270 // * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
00271 
00272 template<class Type>
00273 Ostream& operator<<
00274 (
00275     Ostream& os,
00276     const pointPatchField<Type>& ptf
00277 )
00278 {
00279     ptf.write(os);
00280 
00281     os.check("Ostream& operator<<(Ostream&, const pointPatchField<Type>&)");
00282 
00283     return os;
00284 }
00285 
00286 
00287 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00288 
00289 } // End namespace Foam
00290 
00291 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00292 
00293 #include "newPointPatchField.C"
00294 
00295 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines