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

surfaceInterpolationScheme.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 Description
00025     Abstract base class for surface interpolation schemes.
00026 
00027 \*---------------------------------------------------------------------------*/
00028 
00029 #include "surfaceInterpolationScheme.H"
00030 #include <finiteVolume/volFields.H>
00031 #include <finiteVolume/surfaceFields.H>
00032 #include <finiteVolume/coupledFvPatchField.H>
00033 
00034 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00035 
00036 namespace Foam
00037 {
00038 
00039 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
00040 
00041 // Return weighting factors for scheme given by name in dictionary
00042 template<class Type>
00043 tmp<surfaceInterpolationScheme<Type> > surfaceInterpolationScheme<Type>::New
00044 (
00045     const fvMesh& mesh,
00046     Istream& schemeData
00047 )
00048 {
00049     if (schemeData.eof())
00050     {
00051         FatalIOErrorIn
00052         (
00053             "surfaceInterpolationScheme<Type>::New(const fvMesh&, Istream&)",
00054             schemeData
00055         )   << "Discretisation scheme not specified"
00056             << endl << endl
00057             << "Valid schemes are :" << endl
00058             << MeshConstructorTablePtr_->sortedToc()
00059             << exit(FatalIOError);
00060     }
00061 
00062     word schemeName(schemeData);
00063 
00064     if (surfaceInterpolation::debug || surfaceInterpolationScheme<Type>::debug)
00065     {
00066         Info<< "surfaceInterpolationScheme<Type>::New"
00067                "(const fvMesh&, Istream&)"
00068                " : discretisation scheme = "
00069             << schemeName
00070             << endl;
00071     }
00072 
00073     typename MeshConstructorTable::iterator constructorIter =
00074         MeshConstructorTablePtr_->find(schemeName);
00075 
00076     if (constructorIter == MeshConstructorTablePtr_->end())
00077     {
00078         FatalIOErrorIn
00079         (
00080             "surfaceInterpolationScheme<Type>::New(const fvMesh&, Istream&)",
00081             schemeData
00082         )   << "Unknown discretisation scheme " << schemeName
00083             << endl << endl
00084             << "Valid schemes are :" << endl
00085             << MeshConstructorTablePtr_->sortedToc()
00086             << exit(FatalIOError);
00087     }
00088 
00089     return constructorIter()(mesh, schemeData);
00090 }
00091 
00092 
00093 // Return weighting factors for scheme given by name in dictionary
00094 template<class Type>
00095 tmp<surfaceInterpolationScheme<Type> > surfaceInterpolationScheme<Type>::New
00096 (
00097     const fvMesh& mesh,
00098     const surfaceScalarField& faceFlux,
00099     Istream& schemeData
00100 )
00101 {
00102     if (schemeData.eof())
00103     {
00104         FatalIOErrorIn
00105         (
00106             "surfaceInterpolationScheme<Type>::New"
00107             "(const fvMesh&, const surfaceScalarField&, Istream&)",
00108             schemeData
00109         )   << "Discretisation scheme not specified"
00110             << endl << endl
00111             << "Valid schemes are :" << endl
00112             << MeshConstructorTablePtr_->sortedToc()
00113             << exit(FatalIOError);
00114     }
00115 
00116     word schemeName(schemeData);
00117 
00118     if (surfaceInterpolation::debug || surfaceInterpolationScheme<Type>::debug)
00119     {
00120         Info<< "surfaceInterpolationScheme<Type>::New"
00121                "(const fvMesh&, const surfaceScalarField&, Istream&)"
00122                " : discretisation scheme = "
00123             << schemeName
00124             << endl;
00125     }
00126 
00127     typename MeshFluxConstructorTable::iterator constructorIter =
00128         MeshFluxConstructorTablePtr_->find(schemeName);
00129 
00130     if (constructorIter == MeshFluxConstructorTablePtr_->end())
00131     {
00132         FatalIOErrorIn
00133         (
00134             "surfaceInterpolationScheme<Type>::New"
00135             "(const fvMesh&, const surfaceScalarField&, Istream&)",
00136             schemeData
00137         )   << "Unknown discretisation scheme " << schemeName
00138             << endl << endl
00139             << "Valid schemes are :" << endl
00140             << MeshFluxConstructorTablePtr_->sortedToc()
00141             << exit(FatalIOError);
00142     }
00143 
00144     return constructorIter()(mesh, faceFlux, schemeData);
00145 }
00146 
00147 
00148 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
00149 
00150 template<class Type>
00151 surfaceInterpolationScheme<Type>::~surfaceInterpolationScheme()
00152 {}
00153 
00154 
00155 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
00156 
00157 //- Return the face-interpolate of the given cell field
00158 //  with the given owner and neighbour weighting factors
00159 template<class Type>
00160 tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
00161 surfaceInterpolationScheme<Type>::interpolate
00162 (
00163     const GeometricField<Type, fvPatchField, volMesh>& vf,
00164     const tmp<surfaceScalarField>& tlambdas,
00165     const tmp<surfaceScalarField>& tys
00166 )
00167 {
00168     if (surfaceInterpolation::debug)
00169     {
00170         Info<< "surfaceInterpolationScheme<Type>::uncorrectedInterpolate"
00171                "(const GeometricField<Type, fvPatchField, volMesh>&, "
00172                "const tmp<surfaceScalarField>&, "
00173                "const tmp<surfaceScalarField>&) : "
00174                "interpolating volTypeField from cells to faces "
00175                "without explicit correction"
00176             << endl;
00177     }
00178 
00179     const surfaceScalarField& lambdas = tlambdas();
00180     const surfaceScalarField& ys = tys();
00181 
00182     const Field<Type>& vfi = vf.internalField();
00183     const scalarField& lambda = lambdas.internalField();
00184     const scalarField& y = ys.internalField();
00185 
00186     const fvMesh& mesh = vf.mesh();
00187     const unallocLabelList& P = mesh.owner();
00188     const unallocLabelList& N = mesh.neighbour();
00189 
00190     tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsf
00191     (
00192         new GeometricField<Type, fvsPatchField, surfaceMesh>
00193         (
00194             IOobject
00195             (
00196                 "interpolate("+vf.name()+')',
00197                 vf.instance(),
00198                 vf.db()
00199             ),
00200             mesh,
00201             vf.dimensions()
00202         )
00203     );
00204     GeometricField<Type, fvsPatchField, surfaceMesh>& sf = tsf();
00205 
00206     Field<Type>& sfi = sf.internalField();
00207 
00208     for (register label fi=0; fi<P.size(); fi++)
00209     {
00210         sfi[fi] = lambda[fi]*vfi[P[fi]] + y[fi]*vfi[N[fi]];
00211     }
00212 
00213 
00214     // Interpolate across coupled patches using given lambdas and ys
00215 
00216     forAll (lambdas.boundaryField(), pi)
00217     {
00218         const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
00219         const fvsPatchScalarField& pY = ys.boundaryField()[pi];
00220 
00221         if (vf.boundaryField()[pi].coupled())
00222         {
00223             sf.boundaryField()[pi] =
00224                 pLambda*vf.boundaryField()[pi].patchInternalField()
00225               + pY*vf.boundaryField()[pi].patchNeighbourField();
00226         }
00227         else
00228         {
00229             sf.boundaryField()[pi] = vf.boundaryField()[pi];
00230         }
00231     }
00232 
00233     tlambdas.clear();
00234     tys.clear();
00235 
00236     return tsf;
00237 }
00238 
00239 
00240 //- Return the face-interpolate of the given cell field
00241 //  with the given weigting factors
00242 template<class Type>
00243 tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
00244 surfaceInterpolationScheme<Type>::interpolate
00245 (
00246     const GeometricField<Type, fvPatchField, volMesh>& vf,
00247     const tmp<surfaceScalarField>& tlambdas
00248 )
00249 {
00250     if (surfaceInterpolation::debug)
00251     {
00252         Info<< "surfaceInterpolationScheme<Type>::interpolate"
00253                "(const GeometricField<Type, fvPatchField, volMesh>&, "
00254                "const tmp<surfaceScalarField>&) : "
00255                "interpolating volTypeField from cells to faces "
00256                "without explicit correction"
00257             << endl;
00258     }
00259 
00260     const surfaceScalarField& lambdas = tlambdas();
00261 
00262     const Field<Type>& vfi = vf.internalField();
00263     const scalarField& lambda = lambdas.internalField();
00264 
00265     const fvMesh& mesh = vf.mesh();
00266     const unallocLabelList& P = mesh.owner();
00267     const unallocLabelList& N = mesh.neighbour();
00268 
00269     tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsf
00270     (
00271         new GeometricField<Type, fvsPatchField, surfaceMesh>
00272         (
00273             IOobject
00274             (
00275                 "interpolate("+vf.name()+')',
00276                 vf.instance(),
00277                 vf.db()
00278             ),
00279             mesh,
00280             vf.dimensions()
00281         )
00282     );
00283     GeometricField<Type, fvsPatchField, surfaceMesh>& sf = tsf();
00284 
00285     Field<Type>& sfi = sf.internalField();
00286 
00287     for (register label fi=0; fi<P.size(); fi++)
00288     {
00289         sfi[fi] = lambda[fi]*(vfi[P[fi]] - vfi[N[fi]]) + vfi[N[fi]];
00290     }
00291 
00292     // Interpolate across coupled patches using given lambdas
00293 
00294     forAll (lambdas.boundaryField(), pi)
00295     {
00296         const fvsPatchScalarField& pLambda = lambdas.boundaryField()[pi];
00297 
00298         if (vf.boundaryField()[pi].coupled())
00299         {
00300             tsf().boundaryField()[pi] =
00301                 pLambda*vf.boundaryField()[pi].patchInternalField()
00302              + (1.0 - pLambda)*vf.boundaryField()[pi].patchNeighbourField();
00303         }
00304         else
00305         {
00306             sf.boundaryField()[pi] = vf.boundaryField()[pi];
00307         }
00308     }
00309 
00310     tlambdas.clear();
00311 
00312     return tsf;
00313 }
00314 
00315 
00316 //- Return the face-interpolate of the given cell field
00317 //  with explicit correction
00318 template<class Type>
00319 tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
00320 surfaceInterpolationScheme<Type>::interpolate
00321 (
00322     const GeometricField<Type, fvPatchField, volMesh>& vf
00323 ) const
00324 {
00325     if (surfaceInterpolation::debug)
00326     {
00327         Info<< "surfaceInterpolationScheme<Type>::interpolate"
00328                "(const GeometricField<Type, fvPatchField, volMesh>&) : "
00329             << "interpolating volTypeField from cells to faces"
00330             << endl;
00331     }
00332 
00333     tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tsf
00334         = interpolate(vf, weights(vf));
00335 
00336     if (corrected())
00337     {
00338         tsf() += correction(vf);
00339     }
00340 
00341     return tsf;
00342 }
00343 
00344 
00345 //- Return the face-interpolate of the given cell field
00346 //  with explicit correction
00347 template<class Type>
00348 tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
00349 surfaceInterpolationScheme<Type>::interpolate
00350 (
00351     const tmp<GeometricField<Type, fvPatchField, volMesh> >& tvf
00352 ) const
00353 {
00354     tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > tinterpVf
00355         = interpolate(tvf());
00356     tvf.clear();
00357     return tinterpVf;
00358 }
00359 
00360 
00361 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
00362 
00363 } // End namespace Foam
00364 
00365 // ************************ vim: set sw=4 sts=4 et: ************************ //
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines