FreeFOAM The Cross-Platform CFD Toolkit
pEqn.H
Go to the documentation of this file.
1 {
2  surfaceScalarField alphaf = fvc::interpolate(alpha);
3  surfaceScalarField betaf = scalar(1) - alphaf;
4 
5  volScalarField rUaA = 1.0/UaEqn.A();
6  volScalarField rUbA = 1.0/UbEqn.A();
7 
10 
11  Ua = rUaA*UaEqn.H();
12  Ub = rUbA*UbEqn.H();
13 
14  surfaceScalarField phiDraga =
15  fvc::interpolate(beta/rhoa*dragCoef*rUaA)*phib + rUaAf*(g & mesh.Sf());
16  surfaceScalarField phiDragb =
17  fvc::interpolate(alpha/rhob*dragCoef*rUbA)*phia + rUbAf*(g & mesh.Sf());
18 
19  forAll(p.boundaryField(), patchi)
20  {
21  if (isA<zeroGradientFvPatchScalarField>(p.boundaryField()[patchi]))
22  {
23  phiDraga.boundaryField()[patchi] = 0.0;
24  phiDragb.boundaryField()[patchi] = 0.0;
25  }
26  }
27 
28  phia = (fvc::interpolate(Ua) & mesh.Sf()) + fvc::ddtPhiCorr(rUaA, Ua, phia)
29  + phiDraga;
30  phib = (fvc::interpolate(Ub) & mesh.Sf()) + fvc::ddtPhiCorr(rUbA, Ub, phib)
31  + phiDragb;
32 
33  phi = alphaf*phia + betaf*phib;
34 
35  surfaceScalarField Dp("(rho*(1|A(U)))", alphaf*rUaAf/rhoa + betaf*rUbAf/rhob);
36 
37  for(int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
38  {
39  fvScalarMatrix pEqn
40  (
41  fvm::laplacian(Dp, p) == fvc::div(phi)
42  );
43 
44  pEqn.setReference(pRefCell, pRefValue);
45  pEqn.solve();
46 
47  if (nonOrth == nNonOrthCorr)
48  {
49  surfaceScalarField SfGradp = pEqn.flux()/Dp;
50 
51  phia -= rUaAf*SfGradp/rhoa;
52  phib -= rUbAf*SfGradp/rhob;
53  phi = alphaf*phia + betaf*phib;
54 
55  p.relax();
56  SfGradp = pEqn.flux()/Dp;
57 
58  Ua += (fvc::reconstruct(phiDraga - rUaAf*SfGradp/rhoa));
59  //Ua += rUaA*(fvc::reconstruct(phiDraga/rUaAf - SfGradp/rhoa));
60  Ua.correctBoundaryConditions();
61 
62  Ub += (fvc::reconstruct(phiDragb - rUbAf*SfGradp/rhob));
63  //Ub += rUbA*(fvc::reconstruct(phiDragb/rUbAf - SfGradp/rhob));
64  Ub.correctBoundaryConditions();
65 
66  U = alpha*Ua + beta*Ub;
67  }
68  }
69 }
70