Skip to content

Commit 74123a4

Browse files
xu-23xu-23
committed
correct scheme files
1 parent aae15ac commit 74123a4

File tree

4 files changed

+408
-496
lines changed

4 files changed

+408
-496
lines changed

‎src/dfGenMatrixIncompressible/GenMatrix_p.C‎

Lines changed: 14 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -15,11 +15,23 @@
1515

1616
namespace Foam{
1717

18-
tmp<fvVectorMatrix>
19-
GenMatrix_p(
2018

19+
template<class Type>
20+
tmp
21+
<
22+
GeometricField
23+
<
24+
typename outerProduct<vector, Type>::type,
25+
fvPatchField,
26+
volMesh
27+
>
28+
>
29+
GenMatrix_p(
30+
const GeometricField<Type, fvPatchField, volMesh>& vsf
2131
){
32+
word name("grad(" + vsf.name() + ')');
2233

34+
const fvMesh& mesh = U.mesh();
2335

2436

2537

‎src/dfGenMatrixIncompressible/gaussGrad.C‎

Lines changed: 260 additions & 56 deletions
Original file line numberDiff line numberDiff line change
@@ -23,100 +23,304 @@ License
2323
2424
\*---------------------------------------------------------------------------*/
2525

26-
#include "GenFvMatrix.H"
27-
#include "fvcSurfaceIntegrate.H"
28-
#include "fvMatrices.H"
29-
#include "gaussConvectionScheme.H"
26+
#include "gaussGrad.H"
27+
#include "extrapolatedCalculatedFvPatchField.H"
3028

3129
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
3230

3331
namespace Foam
3432
{
3533

36-
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
3734

38-
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
3935

4036
template<class Type>
41-
tmp<fvMatrix<Type>>
42-
gaussConvectionSchemeFvmDiv
37+
tmp
38+
<
39+
GeometricField
40+
<
41+
typename outerProduct<vector, Type>::type,
42+
fvPatchField,
43+
volMesh
44+
>
45+
>
46+
gaussGradSchemeGrad
4347
(
44-
const surfaceScalarField& faceFlux,
45-
const GeometricField<Type, fvPatchField, volMesh>& vf,
48+
const GeometricField<Type, fvPatchField, volMesh>& vsf
49+
)
50+
{
51+
return gaussGradSchemeGrad(vsf, "grad(" + vsf.name() + ')');
52+
}
53+
54+
template<class Type>
55+
tmp
56+
<
57+
GeometricField
58+
<
59+
typename outerProduct<vector, Type>::type,
60+
fvPatchField,
61+
volMesh
62+
>
63+
>
64+
gaussGradSchemeGrad
65+
(
66+
const GeometricField<Type, fvPatchField, volMesh>& vsf,
4667
const word& name
4768
)
4869
{
49-
Info << "gaussConvectionSchemeFvmDiv start" << endl;
70+
const fvMesh& mesh = vsf.mesh();
71+
72+
typedef typename outerProduct<vector, Type>::type GradType;
73+
typedef GeometricField<GradType, fvPatchField, volMesh> GradFieldType;
5074

51-
const fvMesh& mesh = vf.mesh();
75+
if (!mesh.changing() && mesh.cache(name))
76+
{
77+
if (!mesh.objectRegistry::template foundObject<GradFieldType>(name))
78+
{
79+
solution::cachePrintMessage("Calculating and caching", name, vsf);
80+
tmp<GradFieldType> tgGrad = gaussGradCalcGrad(vsf, name);
81+
regIOobject::store(tgGrad.ptr());
82+
}
5283

53-
tmp<fv::convectionScheme<Type>> cs = fv::convectionScheme<Type>::New(mesh,faceFlux,mesh.divScheme(name));
54-
fv::gaussConvectionScheme<Type>& gcs = dynamic_cast<fv::gaussConvectionScheme<Type>&>(cs.ref());
84+
solution::cachePrintMessage("Retrieving", name, vsf);
85+
GradFieldType& gGrad =
86+
mesh.objectRegistry::template lookupObjectRef<GradFieldType>
87+
(
88+
name
89+
);
5590

56-
tmp<surfaceScalarField> tweights = gcs.interpScheme().weights(vf);
57-
const surfaceScalarField& weights = tweights();
91+
if (gGrad.upToDate(vsf))
92+
{
93+
return gGrad;
94+
}
95+
else
96+
{
97+
solution::cachePrintMessage("Deleting", name, vsf);
98+
gGrad.release();
99+
delete &gGrad;
58100

59-
tmp<fvMatrix<Type>> tfvm
60-
(
61-
new fvMatrix<Type>
62-
(
63-
vf,
64-
faceFlux.dimensions()*vf.dimensions()
65-
)
66-
);
67-
fvMatrix<Type>& fvm = tfvm.ref();
68-
fvm.lower() = -weights.primitiveField()*faceFlux.primitiveField();
69-
fvm.upper() = fvm.lower() + faceFlux.primitiveField();
70-
fvm.negSumDiag();
71-
forAll(vf.boundaryField(), patchi)
72-
{
73-
const fvPatchField<Type>& psf = vf.boundaryField()[patchi];
74-
const fvsPatchScalarField& patchFlux = faceFlux.boundaryField()[patchi];
75-
const fvsPatchScalarField& pw = weights.boundaryField()[patchi];
101+
solution::cachePrintMessage("Recalculating", name, vsf);
102+
tmp<GradFieldType> tgGrad = gaussGradCalcGrad(vsf, name);
103+
104+
solution::cachePrintMessage("Storing", name, vsf);
105+
regIOobject::store(tgGrad.ptr());
106+
GradFieldType& gGrad =
107+
mesh.objectRegistry::template lookupObjectRef<GradFieldType>
108+
(
109+
name
110+
);
76111

77-
fvm.internalCoeffs()[patchi] = patchFlux*psf.valueInternalCoeffs(pw);
78-
fvm.boundaryCoeffs()[patchi] = -patchFlux*psf.valueBoundaryCoeffs(pw);
112+
return gGrad;
113+
}
79114
}
80-
if (gcs.interpScheme().corrected())
115+
else
81116
{
82-
fvm += fvc::surfaceIntegrate(faceFlux*gcs.interpScheme().correction(vf));
117+
if (mesh.objectRegistry::template foundObject<GradFieldType>(name))
118+
{
119+
GradFieldType& gGrad =
120+
mesh.objectRegistry::template lookupObjectRef<GradFieldType>
121+
(
122+
name
123+
);
124+
125+
if (gGrad.ownedByRegistry())
126+
{
127+
solution::cachePrintMessage("Deleting", name, vsf);
128+
gGrad.release();
129+
delete &gGrad;
130+
}
131+
}
132+
133+
solution::cachePrintMessage("Calculating", name, vsf);
134+
return gaussGradCalcGrad(vsf, name);
83135
}
84-
Info << "gaussConvectionSchemeFvmDiv end" << endl;
85-
return tfvm;
86136
}
87137

88138
template<class Type>
89-
tmp<fvMatrix<Type>>
90-
gaussConvectionSchemeFvmDiv
139+
tmp
140+
<
141+
GeometricField
142+
<
143+
typename outerProduct<vector, Type>::type,
144+
fvPatchField,
145+
volMesh
146+
>
147+
>
148+
gaussGradCalcGrad
91149
(
92-
const surfaceScalarField& faceFlux,
93-
const GeometricField<Type, fvPatchField, volMesh>& vf
150+
const GeometricField<Type, fvPatchField, volMesh>& vsf,
151+
const word& name
94152
)
95153
{
96-
word name("div("+faceFlux.name()+','+vf.name()+')');
97-
return gaussConvectionSchemeFvmDiv(faceFlux,vf,name);
154+
const fvMesh& mesh = vsf.mesh();
155+
156+
tmp<surfaceInterpolationScheme<Type>> tinterpScheme_ =
157+
tmp<surfaceInterpolationScheme<Type>>
158+
(
159+
new linear<Type>(mesh)
160+
);
161+
162+
typedef typename outerProduct<vector, Type>::type GradType;
163+
164+
tmp<GeometricField<Type, fvsPatchField, surfaceMesh>> tinterpolate = tinterpScheme_().interpolate(vsf);
165+
166+
tmp<GeometricField<GradType, fvPatchField, volMesh>> tgGrad
167+
(
168+
gaussGradGradf(tinterpolate.ref(), name)
169+
);
170+
GeometricField<GradType, fvPatchField, volMesh>& gGrad = tgGrad.ref();
171+
172+
gaussGradCorrectBoundaryConditions(vsf, gGrad);
173+
174+
return tgGrad;
175+
}
176+
177+
template<class Type>
178+
void gaussGradCorrectBoundaryConditions
179+
(
180+
const GeometricField<Type, fvPatchField, volMesh>& vsf,
181+
GeometricField
182+
<
183+
typename outerProduct<vector, Type>::type, fvPatchField, volMesh
184+
>& gGrad
185+
)
186+
{
187+
typename GeometricField
188+
<
189+
typename outerProduct<vector, Type>::type, fvPatchField, volMesh
190+
>::Boundary& gGradbf = gGrad.boundaryFieldRef();
191+
192+
forAll(vsf.boundaryField(), patchi)
193+
{
194+
if (!vsf.boundaryField()[patchi].coupled())
195+
{
196+
const vectorField n
197+
(
198+
vsf.mesh().Sf().boundaryField()[patchi]
199+
/ vsf.mesh().magSf().boundaryField()[patchi]
200+
);
201+
202+
gGradbf[patchi] += n *
203+
(
204+
vsf.boundaryField()[patchi].snGrad()
205+
- (n & gGradbf[patchi])
206+
);
207+
}
208+
}
209+
}
210+
211+
template<class Type>
212+
tmp
213+
<
214+
GeometricField
215+
<
216+
typename outerProduct<vector, Type>::type,
217+
fvPatchField,
218+
volMesh
219+
>
220+
>
221+
gaussGradGradf
222+
(
223+
const GeometricField<Type, fvsPatchField, surfaceMesh>& ssf,
224+
const word& name
225+
)
226+
{
227+
typedef typename outerProduct<vector, Type>::type GradType;
228+
229+
const fvMesh& mesh = ssf.mesh();
230+
231+
tmp<GeometricField<GradType, fvPatchField, volMesh>> tgGrad
232+
(
233+
new GeometricField<GradType, fvPatchField, volMesh>
234+
(
235+
IOobject
236+
(
237+
name,
238+
ssf.instance(),
239+
mesh,
240+
IOobject::NO_READ,
241+
IOobject::NO_WRITE
242+
),
243+
mesh,
244+
dimensioned<GradType>
245+
(
246+
"0",
247+
ssf.dimensions()/dimLength,
248+
Zero
249+
),
250+
extrapolatedCalculatedFvPatchField<GradType>::typeName
251+
)
252+
);
253+
GeometricField<GradType, fvPatchField, volMesh>& gGrad = tgGrad.ref();
254+
255+
const labelUList& owner = mesh.owner();
256+
const labelUList& neighbour = mesh.neighbour();
257+
const vectorField& Sf = mesh.Sf();
258+
259+
Field<GradType>& igGrad = gGrad;
260+
const Field<Type>& issf = ssf;
261+
262+
forAll(owner, facei)
263+
{
264+
GradType Sfssf = Sf[facei]*issf[facei];
265+
266+
igGrad[owner[facei]] += Sfssf;
267+
igGrad[neighbour[facei]] -= Sfssf;
268+
}
269+
270+
forAll(mesh.boundary(), patchi)
271+
{
272+
const labelUList& pFaceCells =
273+
mesh.boundary()[patchi].faceCells();
274+
275+
const vectorField& pSf = mesh.Sf().boundaryField()[patchi];
276+
277+
const fvsPatchField<Type>& pssf = ssf.boundaryField()[patchi];
278+
279+
forAll(mesh.boundary()[patchi], facei)
280+
{
281+
igGrad[pFaceCells[facei]] += pSf[facei]*pssf[facei];
282+
}
283+
}
284+
285+
igGrad /= mesh.V();
286+
287+
gGrad.correctBoundaryConditions();
288+
289+
return tgGrad;
98290
}
99291

100-
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
101292

102293
template
103-
tmp<fvMatrix<scalar>>
104-
gaussConvectionSchemeFvmDiv
294+
tmp
295+
<
296+
GeometricField
297+
<
298+
typename outerProduct<vector, scalar>::type,
299+
fvPatchField,
300+
volMesh
301+
>
302+
>
303+
gaussGradSchemeGrad
105304
(
106-
const surfaceScalarField& faceFlux,
107-
const GeometricField<scalar, fvPatchField, volMesh>& vf
305+
const GeometricField<scalar, fvPatchField, volMesh>& vsf
108306
);
109307

308+
110309
template
111-
tmp<fvMatrix<vector>>
112-
gaussConvectionSchemeFvmDiv
310+
tmp
311+
<
312+
GeometricField
313+
<
314+
typename outerProduct<vector, vector>::type,
315+
fvPatchField,
316+
volMesh
317+
>
318+
>
319+
gaussGradSchemeGrad
113320
(
114-
const surfaceScalarField& faceFlux,
115-
const GeometricField<vector, fvPatchField, volMesh>& vf
321+
const GeometricField<vector, fvPatchField, volMesh>& vsf
116322
);
117323

118-
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
119-
120324
} // End namespace Foam
121325

122326
// ************************************************************************* //

0 commit comments

Comments
 (0)