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
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037 #ifndef cubic_H
00038 #define cubic_H
00039
00040 #include "linear.H"
00041 #include "gaussGrad.H"
00042
00043
00044
00045 namespace Foam
00046 {
00047
00048
00049
00050
00051
00052 template<class Type>
00053 class cubic
00054 :
00055 public linear<Type>
00056 {
00057
00058
00059
00060 cubic(const cubic&);
00061
00062
00063 void operator=(const cubic&);
00064
00065
00066 public:
00067
00068
00069 TypeName("cubic");
00070
00071
00072
00073
00074
00075 cubic(const fvMesh& mesh)
00076 :
00077 linear<Type>(mesh)
00078 {}
00079
00080
00081 cubic
00082 (
00083 const fvMesh& mesh,
00084 Istream&
00085 )
00086 :
00087 linear<Type>(mesh)
00088 {}
00089
00090
00091 cubic
00092 (
00093 const fvMesh& mesh,
00094 const surfaceScalarField&,
00095 Istream&
00096 )
00097 :
00098 linear<Type>(mesh)
00099 {}
00100
00101
00102
00103
00104
00105 virtual bool corrected() const
00106 {
00107 return true;
00108 }
00109
00110
00111 virtual tmp<GeometricField<Type, fvPatchField, surfaceMesh> >
00112 correction
00113 (
00114 const GeometricField<Type, fvPatchField, volMesh>& vf
00115 ) const
00116 {
00117 const fvMesh& mesh = this->mesh();
00118
00119
00120 const surfaceScalarField& lambda = mesh.weights();
00121
00122 surfaceScalarField kSc =
00123 lambda*(scalar(1) - lambda*(scalar(3) - scalar(2)*lambda));
00124
00125 surfaceScalarField kVecP = sqr(scalar(1) - lambda)*lambda;
00126 surfaceScalarField kVecN = sqr(lambda)*(lambda - scalar(1));
00127
00128 tmp<GeometricField<Type, fvPatchField, surfaceMesh> > tsfCorr
00129 (
00130 new GeometricField<Type, fvPatchField, surfaceMesh>
00131 (
00132 IOobject
00133 (
00134 vf.name(),
00135 mesh.time().timeName(),
00136 mesh
00137 ),
00138 surfaceInterpolationScheme<Type>::interpolate(vf, kSc, -kSc)
00139 )
00140 );
00141
00142 GeometricField<Type, fvPatchField, surfaceMesh>& sfCorr = tsfCorr();
00143
00144 for (direction cmpt=0; cmpt<pTraits<Type>::nComponents; cmpt++)
00145 {
00146 sfCorr.replace
00147 (
00148 cmpt,
00149 sfCorr.component(cmpt)
00150 + (
00151 surfaceInterpolationScheme
00152 <
00153 typename outerProduct
00154 <
00155 vector,
00156 typename pTraits<Type>::cmptType
00157 >::type
00158 >::interpolate
00159 (
00160 fv::gaussGrad
00161 <typename pTraits<Type>::cmptType>(mesh)
00162 .grad(vf.component(cmpt)),
00163 kVecP,
00164 kVecN
00165 ) & mesh.Sf()
00166 )/mesh.magSf()/mesh.surfaceInterpolation::deltaCoeffs()
00167 );
00168 }
00169
00170 forAll (sfCorr.boundaryField(), pi)
00171 {
00172 if (!sfCorr.boundaryField()[pi].coupled())
00173 {
00174 sfCorr.boundaryField()[pi] = pTraits<Type>::zero;
00175 }
00176 }
00177
00178 return tsfCorr;
00179 }
00180 };
00181
00182
00183
00184
00185 }
00186
00187
00188
00189 #endif
00190
00191