CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (http://www.cfd-online.com/Forums/openfoam-solving/)
-   -   implementing this equation (http://www.cfd-online.com/Forums/openfoam-solving/66272-implementing-equation.html)

sambatra July 10, 2009 05:27

implementing this equation
 
Hi there, It's been a whlile I haven't logged in, Please,Can anybody show me how to implement this equation of species in a chemical model? d/dx(rho*Ux*xi)+d/dy(rho*Uy*xi) = Sxi d/dx and d/dy stand for the derivative with respect to x and y Ux,Uy are the x-velocity and y-velocity components xi is the mass fraction of a species Sxi is the source term that generates the species at this stage the pressure,the temperature,U and rho have just been solved Here is what I tried but it diverged: fvm::div(phiv,xi) = Sxi; phiv is surfaceScalarField ("fvc::interpolate(rho*U) & mesh.Sf()") when I run it,it crashes.(error... prinstack...sigpFe) So my question is,how to implement this in the code? Cause I'm quite sure that problem lies in the way I implement it but I just don't know how.I'm not very experienced in OpenFoam Thanks!!

henrik July 12, 2009 14:10

Dear sambatra,

1) Make sure the flux is conservative
2) Make the source implicit
3) Underrelax

Read a book on Finite Volumes if this does not ring a bell.

Henrik

sambatra July 14, 2009 05:26

Dear Henrik,
Thank you so much for the advice,it's kind of you.
Sorry for the delay in replying you.

isabel August 14, 2009 14:28

Sorry for bother you.

I am studying boiling (growing of a vapor bubble in a liquid medium). I have to add a mass source to the continuity equation:

div(u) = cte*grad(T)*grad(gamma)

gamma is the liquid volume fraction, and cte=1.065e-3
So I have programmed the pEqn.H as this:

volScalarField rUA = 1.0/UEqn.A();
surfaceScalarField rUAf = fvc::interpolate(rUA);

U = rUA*UEqn.H();

surfaceScalarField phiU
(
"phiU",
(fvc::interpolate(U) & mesh.Sf()) + fvc::ddtPhiCorr(rUA, rho, U, phi)
);

phi = phiU +
(
fvc::interpolate(interface.sigmaK())*fvc::snGrad(g amma)
- ghf*fvc::snGrad(rho)
)*rUAf*mesh.magSf();

adjustPhi(phi, U, pd);

volVectorField gradT = fvc::grad(T);
volVectorField gradgamma = fvc::grad(gamma);
sourcemasa = -1.065e-3*(gradT & gradgamma);

for(int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
{
fvScalarMatrix pdEqn
(
fvm::laplacian(rUAf, pd) - fvc::div(phi) - sourcemasa
);

pdEqn.setReference(pdRefCell, pdRefValue);

if (corr == nCorr-1 && nonOrth == nNonOrthCorr)
{
pdEqn.solve(mesh.solver(pd.name() + "Final"));
}

else
{
pdEqn.solve(mesh.solver(pd.name()));
}


if (nonOrth == nNonOrthCorr)
{
phi -= pdEqn.flux();
}
}

U += rUA*fvc::reconstruct((phi - phiU)/rUAf);
U.correctBoundaryConditions();





The solver runs Ok, but I have convergence problems. Does anybody know why??

sandy August 16, 2009 08:36

Maybe ....
 
Hi isabel, maybe you can try to add an artifical compression term into this equation as similar as gamma*(1-gamma)*U of the gammaEqn.H? I don't know. But, why you not try it?

sandy August 20, 2009 06:48

Hi isabel, I just remind it. If you add a source term to pdEqn, you think, whether or not this equation can still keep symmetrical?? If not, you should choose PBiCG method to solve it (in system / fvSolution) but not PCG; if yes, you can still use PCG. I can not make sure, because I don't know how to analyze the coefficient matrix structure to the momentum, continuity and gamma fraction volume equations in OpenFOAM, f. e. their eigenvalues, diagnal predominance et. at. Do you have some direct methods? Thanks.

sandy August 20, 2009 10:05

Hi isabel, in fact, I ever tested in my own code (not OpenFOAM) that, whether the pressure equation is symmetrical or not, it can been solved very well by PBiCG method. So, please don't worry about anything and try it. :)

But my preconditioning method is just ILU(0) but not Diagonal-ILU or Multi-grid. What is the matter with them? You think.

isabel August 26, 2009 06:42

Thank you very much, sandy, but when I try to change PCG to PBiCG, I obtain this error:


Unknown symmetric matrix solver PBiCG
Valid symmetric matrix solvers are :
4
(
ICCG
smoothSolver
PCG
GAMG
)


It means that my matrix is still symmetric despite I have added a source.

sandy August 26, 2009 10:50

Hi isabel, I wake up again. Source term = Sp * psi +Sc, so Sp actually is acutally added into diagonal elements in the ldumatrix. So, the pressure equation is still symmetric. If you add a ddt term as like as the sonicTurbFoam, the pressure equation will loss symmetric, I think. Sorry...

sandy August 27, 2009 08:55

Hi isabel, I wake up again. Source term = Sp * psi +Sc, so Sp actually is acutally added into diagonal elements in the ldumatrix. So, the pressure equation is still symmetric. If you add a ddt term as like as the sonicTurbFoam, the pressure equation will loss symmetric, I think. Sorry... :o


All times are GMT -4. The time now is 19:48.