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 #ifndef localBlended_H
00037 #define localBlended_H
00038
00039 #include "surfaceInterpolationScheme.H"
00040
00041
00042
00043 namespace Foam
00044 {
00045
00046
00047
00048
00049
00050 template<class Type>
00051 class localBlended
00052 :
00053 public surfaceInterpolationScheme<Type>
00054 {
00055
00056
00057
00058 tmp<surfaceInterpolationScheme<Type> > tScheme1_;
00059
00060
00061 tmp<surfaceInterpolationScheme<Type> > tScheme2_;
00062
00063
00064
00065 localBlended(const localBlended&);
00066
00067
00068 void operator=(const localBlended&);
00069
00070
00071 public:
00072
00073
00074 TypeName("localBlended");
00075
00076
00077
00078
00079
00080
00081
00082 localBlended
00083 (
00084 const fvMesh& mesh,
00085 Istream& is
00086 )
00087 :
00088 surfaceInterpolationScheme<Type>(mesh),
00089 tScheme1_
00090 (
00091 surfaceInterpolationScheme<Type>::New(mesh, is)
00092 ),
00093 tScheme2_
00094 (
00095 surfaceInterpolationScheme<Type>::New(mesh, is)
00096 )
00097 {}
00098
00099
00100 localBlended
00101 (
00102 const fvMesh& mesh,
00103 const surfaceScalarField& faceFlux,
00104 Istream& is
00105 )
00106 :
00107 surfaceInterpolationScheme<Type>(mesh),
00108 tScheme1_
00109 (
00110 surfaceInterpolationScheme<Type>::New(mesh, faceFlux, is)
00111 ),
00112 tScheme2_
00113 (
00114 surfaceInterpolationScheme<Type>::New(mesh, faceFlux, is)
00115 )
00116 {}
00117
00118
00119
00120
00121
00122 tmp<surfaceScalarField> weights
00123 (
00124 const GeometricField<Type, fvPatchField, volMesh>& vf
00125 ) const
00126 {
00127 const surfaceScalarField& blendingFactor =
00128 this->mesh().objectRegistry::
lookupObject<const surfaceScalarField>
00129 (
00130 word(vf.name() + "BlendingFactor")
00131 );
00132
00133 return
00134 blendingFactor*tScheme1_().weights(vf)
00135 + (scalar(1) - blendingFactor)*tScheme2_().weights(vf);
00136 }
00137
00138
00139
00140 tmp<GeometricField<Type, fvPatchField, surfaceMesh> >
00141 interpolate(const GeometricField<Type, fvPatchField, volMesh>& vf) const
00142 {
00143 const surfaceScalarField& blendingFactor =
00144 (
00145 this->mesh().objectRegistry::
lookupObject<const surfaceScalarField>
00146 (
00147 word(vf.name() + "BlendingFactor")
00148 )
00149 );
00150
00151 return
00152 blendingFactor*tScheme1_().interpolate(vf)
00153 + (scalar(1) - blendingFactor)*tScheme2_().interpolate(vf);
00154 }
00155 };
00156
00157
00158
00159
00160 }
00161
00162
00163
00164 #endif
00165
00166
00167