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 <finiteVolume/volFields.H>
00027 #include <finiteVolume/surfaceFields.H>
00028 #include <finiteVolume/fvMatrix.H>
00029
00030
00031
00032 template<class Type>
00033 Foam::tmp<Foam::fvMatrix<Type> >
00034 Foam::fvm::Su
00035 (
00036 const DimensionedField<Type, volMesh>& su,
00037 GeometricField<Type, fvPatchField, volMesh>& vf
00038 )
00039 {
00040 const fvMesh& mesh = vf.mesh();
00041
00042 tmp<fvMatrix<Type> > tfvm
00043 (
00044 new fvMatrix<Type>
00045 (
00046 vf,
00047 dimVol*su.dimensions()
00048 )
00049 );
00050 fvMatrix<Type>& fvm = tfvm();
00051
00052 fvm.source() -= mesh.V()*su.field();
00053
00054 return tfvm;
00055 }
00056
00057 template<class Type>
00058 Foam::tmp<Foam::fvMatrix<Type> >
00059 Foam::fvm::Su
00060 (
00061 const tmp<DimensionedField<Type, volMesh> >& tsu,
00062 GeometricField<Type, fvPatchField, volMesh>& vf
00063 )
00064 {
00065 tmp<fvMatrix<Type> > tfvm = fvm::Su(tsu(), vf);
00066 tsu.clear();
00067 return tfvm;
00068 }
00069
00070 template<class Type>
00071 Foam::tmp<Foam::fvMatrix<Type> >
00072 Foam::fvm::Su
00073 (
00074 const tmp<GeometricField<Type, fvPatchField, volMesh> >& tsu,
00075 GeometricField<Type, fvPatchField, volMesh>& vf
00076 )
00077 {
00078 tmp<fvMatrix<Type> > tfvm = fvm::Su(tsu(), vf);
00079 tsu.clear();
00080 return tfvm;
00081 }
00082
00083 template<class Type>
00084 Foam::zeroField
00085 Foam::fvm::Su
00086 (
00087 const zeroField&,
00088 GeometricField<Type, fvPatchField, volMesh>& vf
00089 )
00090 {
00091 return zeroField();
00092 }
00093
00094
00095 template<class Type>
00096 Foam::tmp<Foam::fvMatrix<Type> >
00097 Foam::fvm::Sp
00098 (
00099 const DimensionedField<scalar, volMesh>& sp,
00100 GeometricField<Type, fvPatchField, volMesh>& vf
00101 )
00102 {
00103 const fvMesh& mesh = vf.mesh();
00104
00105 tmp<fvMatrix<Type> > tfvm
00106 (
00107 new fvMatrix<Type>
00108 (
00109 vf,
00110 dimVol*sp.dimensions()*vf.dimensions()
00111 )
00112 );
00113 fvMatrix<Type>& fvm = tfvm();
00114
00115 fvm.diag() += mesh.V()*sp.field();
00116
00117 return tfvm;
00118 }
00119
00120 template<class Type>
00121 Foam::tmp<Foam::fvMatrix<Type> >
00122 Foam::fvm::Sp
00123 (
00124 const tmp<DimensionedField<scalar, volMesh> >& tsp,
00125 GeometricField<Type, fvPatchField, volMesh>& vf
00126 )
00127 {
00128 tmp<fvMatrix<Type> > tfvm = fvm::Sp(tsp(), vf);
00129 tsp.clear();
00130 return tfvm;
00131 }
00132
00133 template<class Type>
00134 Foam::tmp<Foam::fvMatrix<Type> >
00135 Foam::fvm::Sp
00136 (
00137 const tmp<volScalarField>& tsp,
00138 GeometricField<Type, fvPatchField, volMesh>& vf
00139 )
00140 {
00141 tmp<fvMatrix<Type> > tfvm = fvm::Sp(tsp(), vf);
00142 tsp.clear();
00143 return tfvm;
00144 }
00145
00146
00147 template<class Type>
00148 Foam::tmp<Foam::fvMatrix<Type> >
00149 Foam::fvm::Sp
00150 (
00151 const dimensionedScalar& sp,
00152 GeometricField<Type, fvPatchField, volMesh>& vf
00153 )
00154 {
00155 const fvMesh& mesh = vf.mesh();
00156
00157 tmp<fvMatrix<Type> > tfvm
00158 (
00159 new fvMatrix<Type>
00160 (
00161 vf,
00162 dimVol*sp.dimensions()*vf.dimensions()
00163 )
00164 );
00165 fvMatrix<Type>& fvm = tfvm();
00166
00167 fvm.diag() += mesh.V()*sp.value();
00168
00169 return tfvm;
00170 }
00171
00172 template<class Type>
00173 Foam::zeroField
00174 Foam::fvm::Sp
00175 (
00176 const zeroField&,
00177 GeometricField<Type, fvPatchField, volMesh>&
00178 )
00179 {
00180 return zeroField();
00181 }
00182
00183
00184 template<class Type>
00185 Foam::tmp<Foam::fvMatrix<Type> >
00186 Foam::fvm::SuSp
00187 (
00188 const DimensionedField<scalar, volMesh>& susp,
00189 GeometricField<Type, fvPatchField, volMesh>& vf
00190 )
00191 {
00192 const fvMesh& mesh = vf.mesh();
00193
00194 tmp<fvMatrix<Type> > tfvm
00195 (
00196 new fvMatrix<Type>
00197 (
00198 vf,
00199 dimVol*susp.dimensions()*vf.dimensions()
00200 )
00201 );
00202 fvMatrix<Type>& fvm = tfvm();
00203
00204 fvm.diag() += mesh.V()*max(susp.field(), scalar(0));
00205
00206 fvm.source() -= mesh.V()*min(susp.field(), scalar(0))
00207 *vf.internalField();
00208
00209 return tfvm;
00210 }
00211
00212 template<class Type>
00213 Foam::tmp<Foam::fvMatrix<Type> >
00214 Foam::fvm::SuSp
00215 (
00216 const tmp<DimensionedField<scalar, volMesh> >& tsusp,
00217 GeometricField<Type, fvPatchField, volMesh>& vf
00218 )
00219 {
00220 tmp<fvMatrix<Type> > tfvm = fvm::SuSp(tsusp(), vf);
00221 tsusp.clear();
00222 return tfvm;
00223 }
00224
00225 template<class Type>
00226 Foam::tmp<Foam::fvMatrix<Type> >
00227 Foam::fvm::SuSp
00228 (
00229 const tmp<volScalarField>& tsusp,
00230 GeometricField<Type, fvPatchField, volMesh>& vf
00231 )
00232 {
00233 tmp<fvMatrix<Type> > tfvm = fvm::SuSp(tsusp(), vf);
00234 tsusp.clear();
00235 return tfvm;
00236 }
00237
00238 template<class Type>
00239 Foam::zeroField
00240 Foam::fvm::SuSp
00241 (
00242 const zeroField&,
00243 GeometricField<Type, fvPatchField, volMesh>& vf
00244 )
00245 {
00246 return zeroField();
00247 }
00248
00249
00250