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 <errorEstimation/resErrorLaplacian.H>
00027
00028
00029
00030 namespace Foam
00031 {
00032
00033
00034
00035 namespace resError
00036 {
00037
00038 template<class Type>
00039 tmp<errorEstimate<Type> >
00040 laplacian
00041 (
00042 const GeometricField<Type, fvPatchField, volMesh>& vf
00043 )
00044 {
00045 surfaceScalarField Gamma
00046 (
00047 IOobject
00048 (
00049 "gamma",
00050 vf.time().constant(),
00051 vf.db(),
00052 IOobject::NO_READ
00053 ),
00054 vf.mesh(),
00055 dimensionedScalar("1", dimless, 1.0)
00056 );
00057
00058 return resError::laplacian(Gamma, vf);
00059 }
00060
00061
00062 template<class Type>
00063 tmp<errorEstimate<Type> >
00064 laplacian
00065 (
00066 const dimensionedScalar& gamma,
00067 const GeometricField<Type, fvPatchField, volMesh>& vf
00068 )
00069 {
00070 surfaceScalarField Gamma
00071 (
00072 IOobject
00073 (
00074 gamma.name(),
00075 vf.time().timeName(),
00076 vf.db(),
00077 IOobject::NO_READ
00078 ),
00079 vf.mesh(),
00080 gamma
00081 );
00082
00083 return resError::laplacian(Gamma, vf);
00084 }
00085
00086
00087 template<class Type>
00088 tmp<errorEstimate<Type> >
00089 laplacian
00090 (
00091 const volScalarField& gamma,
00092 const GeometricField<Type, fvPatchField, volMesh>& vf
00093 )
00094 {
00095 return resError::laplacian(fvc::interpolate(gamma), vf);
00096 }
00097
00098
00099 template<class Type>
00100 tmp<errorEstimate<Type> >
00101 laplacian
00102 (
00103 const tmp<volScalarField>& tgamma,
00104 const GeometricField<Type, fvPatchField, volMesh>& vf
00105 )
00106 {
00107 tmp<errorEstimate<Type> > Laplacian(resError::laplacian(tgamma(), vf));
00108 tgamma.clear();
00109 return Laplacian;
00110 }
00111
00112
00113 template<class Type>
00114 tmp<errorEstimate<Type> >
00115 laplacian
00116 (
00117 const surfaceScalarField& gamma,
00118 const GeometricField<Type, fvPatchField, volMesh>& vf
00119 )
00120 {
00121 const fvMesh& mesh = vf.mesh();
00122
00123 const scalarField& vols = mesh.V();
00124 const surfaceVectorField& Sf = mesh.Sf();
00125 const surfaceScalarField magSf = mesh.magSf();
00126 const fvPatchList& patches = mesh.boundary();
00127 const unallocLabelList& owner = mesh.owner();
00128 const unallocLabelList& neighbour = mesh.neighbour();
00129
00130 const surfaceScalarField& delta =
00131 mesh.surfaceInterpolation::deltaCoeffs();
00132
00133 Field<Type> res(vols.size(), pTraits<Type>::zero);
00134 scalarField aNorm(vols.size(), 0.0);
00135
00136
00137 GeometricField
00138 <
00139 typename outerProduct<vector, Type>::type, fvPatchField, volMesh
00140 > gradVf = fvc::grad(vf);
00141
00142
00143 forAll (owner, faceI)
00144 {
00145
00146
00147
00148 res[owner[faceI]] -=
00149 gamma[faceI]*(Sf[faceI] & gradVf[owner[faceI]]);
00150
00151 aNorm[owner[faceI]] += delta[faceI]*gamma[faceI]*magSf[faceI];
00152
00153
00154
00155
00156 res[neighbour[faceI]] +=
00157 gamma[faceI]*(Sf[faceI] & gradVf[neighbour[faceI]]);
00158
00159 aNorm[neighbour[faceI]] += delta[faceI]*gamma[faceI]*magSf[faceI];
00160
00161 }
00162
00163 forAll (patches, patchI)
00164 {
00165 const vectorField& patchSf = Sf.boundaryField()[patchI];
00166 const scalarField& patchMagSf = magSf.boundaryField()[patchI];
00167 const scalarField& patchGamma = gamma.boundaryField()[patchI];
00168 const scalarField& patchDelta = delta.boundaryField()[patchI];
00169
00170 const labelList& fCells = patches[patchI].faceCells();
00171
00172 forAll (fCells, faceI)
00173 {
00174
00175 res[fCells[faceI]] -=
00176 patchGamma[faceI]*
00177 (
00178 patchSf[faceI] & gradVf[fCells[faceI]]
00179 );
00180
00181 aNorm[fCells[faceI]] +=
00182 patchDelta[faceI]*patchGamma[faceI]*patchMagSf[faceI];
00183 }
00184 }
00185
00186 res /= vols;
00187 aNorm /= vols;
00188
00189 return tmp<errorEstimate<Type> >
00190 (
00191 new errorEstimate<Type>
00192 (
00193 vf,
00194 delta.dimensions()*gamma.dimensions()*magSf.dimensions()
00195 *vf.dimensions(),
00196 res,
00197 aNorm
00198 )
00199 );
00200 }
00201
00202 template<class Type>
00203 tmp<errorEstimate<Type> >
00204 laplacian
00205 (
00206 const tmp<surfaceScalarField>& tgamma,
00207 const GeometricField<Type, fvPatchField, volMesh>& vf
00208 )
00209 {
00210 tmp<errorEstimate<Type> > tresError(resError::laplacian(tgamma(), vf));
00211 tgamma.clear();
00212 return tresError;
00213 }
00214
00215
00216 template<class Type>
00217 tmp<errorEstimate<Type> >
00218 laplacian
00219 (
00220 const volTensorField& gamma,
00221 const GeometricField<Type, fvPatchField, volMesh>& vf
00222 )
00223 {
00224 const fvMesh& mesh = vf.mesh();
00225
00226 return resError::laplacian
00227 (
00228 (mesh.Sf() & fvc::interpolate(gamma) & mesh.Sf())
00229 /sqr(mesh.magSf()),
00230 vf
00231 );
00232 }
00233
00234 template<class Type>
00235 tmp<errorEstimate<Type> >
00236 laplacian
00237 (
00238 const tmp<volTensorField>& tgamma,
00239 const GeometricField<Type, fvPatchField, volMesh>& vf
00240 )
00241 {
00242 tmp<errorEstimate<Type> > Laplacian = resError::laplacian(tgamma(), vf);
00243 tgamma.clear();
00244 return Laplacian;
00245 }
00246
00247
00248 template<class Type>
00249 tmp<errorEstimate<Type> >
00250 laplacian
00251 (
00252 const surfaceTensorField& gamma,
00253 const GeometricField<Type, fvPatchField, volMesh>& vf
00254 )
00255 {
00256 const fvMesh& mesh = vf.mesh();
00257
00258 return resError::laplacian
00259 (
00260 (mesh.Sf() & gamma & mesh.Sf())/sqr(mesh.magSf()),
00261 vf
00262 );
00263 }
00264
00265 template<class Type>
00266 tmp<errorEstimate<Type> >
00267 laplacian
00268 (
00269 const tmp<surfaceTensorField>& tgamma,
00270 const GeometricField<Type, fvPatchField, volMesh>& vf
00271 )
00272 {
00273 tmp<errorEstimate<Type> > Laplacian = resError::laplacian(tgamma(), vf);
00274 tgamma.clear();
00275 return Laplacian;
00276 }
00277
00278
00279
00280
00281 }
00282
00283
00284
00285 }
00286
00287