
[Sponsors] 
March 9, 2016, 17:08 
Minimal hacking to interFoam to add funny gravity

#1 
New Member
Daniel Duque
Join Date: Jan 2011
Location: ETSIN, Madrid
Posts: 16
Rep Power: 8 
Hi,
I had this idea to model a RayleighTaylor instability for the simplest possible case. I needed a funny gravity field that would pull up on one of the phases, and down on the other. I hacked the interFoam solver minimally, and I think I have got correct results. My first trial was actually succesful. UEqn.H was changed to: Code:
MRF.correctBoundaryVelocity(U); fvVectorMatrix UEqn ( fvm::ddt(rho, U) + fvm::div(rhoPhi, U) + MRF.DDt(rho, U) + turbulence>divDevRhoReff(rho, U) == fvOptions(rho, U) // fake gravity field + (2*alpha11) * rho * g ); UEqn.relax(); fvOptions.constrain(UEqn); if (pimple.momentumPredictor()) { solve ( UEqn == fvc::reconstruct ( ( mixture.surfaceTensionForce()  ghf*fvc::snGrad(rho) //  fvc::snGrad(p_rgh)  fvc::snGrad(p) ) * mesh.magSf() ) ); fvOptions.correct(U); } But then I thought adding the gravity on the righthand side of the equation would be more correct: Code:
MRF.correctBoundaryVelocity(U); fvVectorMatrix UEqn ( fvm::ddt(rho, U) + fvm::div(rhoPhi, U) + MRF.DDt(rho, U) + turbulence>divDevRhoReff(rho, U) == fvOptions(rho, U) ); UEqn.relax(); fvOptions.constrain(UEqn); if (pimple.momentumPredictor()) { solve ( UEqn == fvc::reconstruct ( ( mixture.surfaceTensionForce()  ghf*fvc::snGrad(rho) //  fvc::snGrad(p_rgh)  fvc::snGrad(p) ) * mesh.magSf() ) // fake gravity field + (2*alpha11) * rho * g ); fvOptions.correct(U); } comparison_step_60.jpg Thanks in advance for any help! 

March 9, 2016, 21:54 

#2 
Senior Member
Pablo Higuera
Join Date: Jan 2011
Location: Singapore
Posts: 379
Rep Power: 10 
Hi Daniel,
probably you have not switched on the momentum predictor in fvSolution, it is usually switched off by default. Best, Pablo 

March 10, 2016, 06:09 

#3  
New Member
Daniel Duque
Join Date: Jan 2011
Location: ETSIN, Madrid
Posts: 16
Rep Power: 8 
Quote:
You are right, it would seem that " if (pimple.momentumPredictor()) " would require the momentum predictor step. However, I have tried enabling it and now the velocity seems to evolve ... but not the pressure, which does not change, or the color field. Aside from " momentumPredictor yes;" at the PIMPLE section, I had to add the following to fvSolutions. Code:
UFinal { $U; solver smoothSolver; smoother symGaussSeidel; tolerance 1e06; relTol 0; } [It's actually grad( p  rho g z ) + rho z grad(rho), but it's the same ] 

March 10, 2016, 08:12 

#4 
Senior Member
Pablo Higuera
Join Date: Jan 2011
Location: Singapore
Posts: 379
Rep Power: 10 
Hi again,
not really, you should always include it in the first one (to assemble the system of equations properly), and if you want a momentum predictor, also in the second equation. Momentum predictor just solves an intermediate step and updates U. Then it applies that value to obtain the independent terms (explicit) and solve the PPE. Otherwise, the U used is the one from the beginning of the time step. Best, Pablo 

March 10, 2016, 09:03 

#5 
New Member
Daniel Duque
Join Date: Jan 2011
Location: ETSIN, Madrid
Posts: 16
Rep Power: 8 
Thanks, that seems to work. I find little difference with the previous results, but that may be ok. To check incompressibility, I am running foamCalc div phi, and the resulting field is between about 3e4 and 3e4. I guess that's small. Perhaps all this could go in another thread.


Tags 
interfoam;ueqn, rayleightaylor 
Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
InterFoam: add a source term in alpha eq.  Alucard  OpenFOAM Programming & Development  11  August 3, 2017 02:23 
InterFoam: Add an equation to only one phase  voingiappone  OpenFOAM Programming & Development  36  July 30, 2016 20:39 
Add source term in alphaEqn.H of interFoam  tayo  OpenFOAM  1  October 23, 2013 03:40 
how to add sourceterm in momentum equation for interFoam?  anon_g  OpenFOAM  9  October 18, 2011 12:47 
add temperature equation to interFoam  nygbook  OpenFOAM Running, Solving & CFD  3  September 4, 2011 07:28 