FreeFOAM The Cross-Platform CFD Toolkit
Hosted by SourceForge:
Get FreeFOAM at SourceForge.net.
            Fast, secure and Free Open Source software downloads

pEqn.H

Go to the documentation of this file.
00001 {
00002     surfaceScalarField alphaf = fvc::interpolate(alpha);
00003     surfaceScalarField betaf = scalar(1) - alphaf;
00004 
00005     volScalarField rUaA = 1.0/UaEqn.A();
00006     volScalarField rUbA = 1.0/UbEqn.A();
00007 
00008     surfaceScalarField rUaAf = fvc::interpolate(rUaA);
00009     surfaceScalarField rUbAf = fvc::interpolate(rUbA);
00010 
00011     Ua = rUaA*UaEqn.H();
00012     Ub = rUbA*UbEqn.H();
00013 
00014     surfaceScalarField phiDraga =
00015         fvc::interpolate(beta/rhoa*dragCoef*rUaA)*phib + rUaAf*(g & mesh.Sf());
00016     surfaceScalarField phiDragb =
00017         fvc::interpolate(alpha/rhob*dragCoef*rUbA)*phia + rUbAf*(g & mesh.Sf());
00018 
00019     forAll(p.boundaryField(), patchi)
00020     {
00021         if (isA<zeroGradientFvPatchScalarField>(p.boundaryField()[patchi]))
00022         {
00023             phiDraga.boundaryField()[patchi] = 0.0;
00024             phiDragb.boundaryField()[patchi] = 0.0;
00025         }
00026     }
00027 
00028     phia = (fvc::interpolate(Ua) & mesh.Sf()) + fvc::ddtPhiCorr(rUaA, Ua, phia)
00029         + phiDraga;
00030     phib = (fvc::interpolate(Ub) & mesh.Sf()) + fvc::ddtPhiCorr(rUbA, Ub, phib)
00031         + phiDragb;
00032 
00033     phi = alphaf*phia + betaf*phib;
00034 
00035     surfaceScalarField Dp("(rho*(1|A(U)))", alphaf*rUaAf/rhoa + betaf*rUbAf/rhob);
00036 
00037     for(int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
00038     {
00039         fvScalarMatrix pEqn
00040         (
00041             fvm::laplacian(Dp, p) == fvc::div(phi)
00042         );
00043 
00044         pEqn.setReference(pRefCell, pRefValue);
00045         pEqn.solve();
00046 
00047         if (nonOrth == nNonOrthCorr)
00048         {
00049             surfaceScalarField SfGradp = pEqn.flux()/Dp;
00050 
00051             phia -= rUaAf*SfGradp/rhoa;
00052             phib -= rUbAf*SfGradp/rhob;
00053             phi = alphaf*phia + betaf*phib;
00054 
00055             p.relax();
00056             SfGradp = pEqn.flux()/Dp;
00057 
00058             Ua += (fvc::reconstruct(phiDraga - rUaAf*SfGradp/rhoa));
00059             //Ua += rUaA*(fvc::reconstruct(phiDraga/rUaAf - SfGradp/rhoa));
00060             Ua.correctBoundaryConditions();
00061 
00062             Ub += (fvc::reconstruct(phiDragb - rUbAf*SfGradp/rhob));
00063             //Ub += rUbA*(fvc::reconstruct(phiDragb/rUbAf - SfGradp/rhob));
00064             Ub.correctBoundaryConditions();
00065 
00066             U = alpha*Ua + beta*Ub;
00067         }
00068     }
00069 }
00070 
00071 #include <finiteVolume/continuityErrs.H>
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines