Hi OpenFOAMers Wish you all
 

Hi OpenFOAMers

Wish you all a nice evening

I have attempted to modify PISO based interfoam solver to the Transient interFoam solver which uses MRF as well. I am posting here the code and hope to get some expert feedback. Please correct me if I am wrong somewhere.

--\*---------------------------------------------------------------------------* /

#include "fvCFD.H"
#include "MULES.H"
#include "subCycle.H"
#include "interfaceProperties.H"
#include "twoPhaseMixture.H"
#include "MRFZones.H"
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

int main(int argc, char *argv[])
{

# include "setRootCase.H"
# include "createTime.H"
# include "createMesh.H"
# include "readEnvironmentalProperties.H"
# include "readPISOControls.H"
# include "initContinuityErrs.H"
# include "createFields.H"
# include "readTimeControls.H"
# include "correctPhi.H"
# include "setInitialDeltaT.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

Info<< "\nStarting time loop\n" << endl;

for (runTime++; !runTime.end(); runTime++)
{
# include "readPISOControls.H"
# include "readTimeControls.H"
# include "CourantNo.H"

Info<< "Time = " << runTime.timeName() << nl << endl;

twoPhaseProperties.correct();

# include "gammaEqnSubCycle.H"

surfaceScalarField muf = twoPhaseProperties.muf();

// pressure velocity Transient SIMPLE corrector loop

for (int corr=0; corr<nCorr; corr++)
{
tmp<fvvectormatrix> UEqn
(
fvm::ddt(rho, U)
+ fvm::div(rhoPhi, U)
- fvm::laplacian(muf, U)
- (fvc::grad(U) & fvc::grad(muf))
);


mrfZone_1.addCoriolis(UEqn());

UEqn().relax();

solve
(
UEqn()
==
fvc::reconstruct
(
(
fvc::interpolate(interface.sigmaK())*fvc::snGrad(g amma)
- ghf*fvc::snGrad(rho)
- fvc::snGrad(pd)
) * mesh.magSf()
)
);

pd.boundaryField().updateCoeffs();

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

U = rUA*UEqn().H();
UEqn.clear();

surfaceScalarField phiU
(
"phiU",
(fvc::interpolate(U) & mesh.Sf())
);
phi = phiU +
(
fvc::interpolate(interface.sigmaK())*fvc::snGrad(g amma)
- ghf*fvc::snGrad(rho)
)*rUAf*mesh.magSf();


mrfZone_1.relativeFlux(phi);
adjustPhi(phi, U, pd);

pd.storePrevIter();

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

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();
}
}

# include "movingMeshRhoUContinuityErrs.H"

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

}


runTime.write();

Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
<< " ClockTime = " << runTime.elapsedClockTime() << " s"
<< nl << endl;
}

Info<< "End\n" << endl;

return(0);
}

--------------------------------
Thanks in advance
Best Regards
Jaswinder

MRF interFoam solver
 

I realize this thread is a bit stale, but I wanted to give an update in this area. Based on the code posted by Jaswinder, I made a MRFinterFoam solver and have been doing some testing. The first test case was to take the setup from the MRFSimpleFoam tutorial on the wiki (http://openfoamwiki.net/index.php/MRFSimpleFoam) which consists of a rotating cube in a cylinder. For this test I started with the liquid level just below the top of the rotating box. Originally I tried an MRF omega value that I though was reasonable (~10-100) and got problems so I changed it to 1.0 (I assume this is radians per second, but the result make me suspicious). Here is the gamma and U distributions on a couple cross-sections at t=0.5s :
http://sites.google.com/site/kwardle..._0p5s_1rps.tif

http://sites.google.com/site/kwardle..._0p5s_1rps.tif

1 rad/s is pretty slow rotation, the maximum velocity should be very small, but as you can see from the the U image things are up to ~0.5 m/s.

Is there some problem with using MRFZones with interfoam that is causing unphysical results?

MRF interFoam solver #2
 

Second example, same problem. For another test, I constructed a bottom-fed centrifuge type geometry where there are two feeds coming into a T junction which flows upward into the centrifuge bowl (the rotating domain). The centrifuge has two small baffles and a weir above which the fluid flows out radially. The radius of the centrifuge cylinder is only 1cm.
Here you can see the geometry and the initially gamma field:

http://sites.google.com/site/kwardle...mega0p1_0s.tif

The inlet velocity on the two inlets was 0.1m/s.

For this MRFZones was set to:
1
(
rotor
{
patches (rotwalls outlet);
origin origin [0 1 0 0 0 0 0] (0 0 0);
axis axis [0 0 0 0 0 0 0] (0 0 1);
omega omega [0 0 -1 0 0 0 0] 0.1;
}
)
Incidentally, I also tried leaving out the two patches and didn't see any difference in the solution. Again, I also tried initially a larger omega, but got extremely small time steps (indicating large velocity).

Here is the gamma field and U vectors for t=1s
http://sites.google.com/site/kwardle...ga0p1_1s-2.tif

You can see here that the max velocities are up to 1.19 m/s. For this radius (r=1cm) and omega = 0.1 rad/s, the outer wall velocity should only be U_tan = r*omega = 0.1*0.01 = 0.001 m/s = 1mm/s = snail speed. Something is wrong. Also I didn't show the axes, but the rotation here is clockwise--perhaps I am missing something, but I would have thought that omega would have the standard definition and the rotation of the flow should be counter-clockwise.

Can anyone offer help? I would be happy to send and/or post the solver I am using and the files for this test case.

Hi Kent

I got your mail. Send me the solver and test case . My guess is that there
is some very basic mistake in the solver. Back then I just left it like that as there was no feedback. We can try to bring this solver back to life :-)

Kind Regards
Jaswi

MRFinterFoam solver
 

Atached is the MRFinterFoam solver that I constructed by modifying the interFoam solver based on the code posted by Jaswinder above along with one based directly on his code (MRFinterFoamJPS) with only a few modifications needed to make it compile and run without errors. Also, here are the two test cases described in my previous posts. I get essentially the same behavior with either solver.

http://sites.google.com/site/kwardle..._rotbox.tar.gz

http://sites.google.com/site/kwardle...tBaffle.tar.gz

Thank you for taking a look at this!
-Kent

Hi Kent

To start with I have first done some first principle checks using the MRFSimpleFOAM on the rotating box case setup. The velocity values and direction for this setup is correct.

Here is a picture of velocity field. As expected the maximum velocity is at the diagonals of the cube. The cube is 1.5 cm (LxBxH) and at 10 rad/sec the velocity calculates to 0.021 m/sec

The inner cylinder shows the size of the rotating zone - rotor.


I will analyze the MRFSimpleFOAM code once again and report back tomorrow about my findings. On that basis we will able to figure out what are we doing wrong in MRFInterFAOM. My wildest guess is that water's density is having some effect here as we have mass flux term rhoPhi in this case and a quick nudge is creating those velocities. I have to check what happens If I increase the velocity gradually .

Best Regards
Jaswi

Hello Guys,

besides phi, there is another flux in this solver: phiU. I guess you should insert a

mrfZones.relativeFlux(phiU);

Best regards, Hannes

Hi Hannes

Thanks for the tip.

I guess the correct formulation will turn up once we have the VOF equations in rotating frame of reference.

Any idea or reference will be really great. :)

Regards
Jaswi

MRFInterFoam solver
 

Jaswi,
Thought I would post an update on this since you have given of your time to look into this.

Here is a working MRFInterFoam solver written by Henry Weller (I finally got a support contract).
MRFInterFoam.tgz

Apparently, it requires you to download and recompile the most recent 1.5.x distribution. I have done this and it works great. Note that this solver also solves directly for p (rather than pd--a big plus for me if you have read some of my other posts about issues with static pressure and pressure outlets) a change that will be part of the interFoam solvers for 1.6. Not sure if this solver will be included in the 1.6 release (not sure when the release will be either, but I was told 'soon').

BTW,
I think Hannes is right--it should have been "mrfZones.relativeFlux(phiU);" in pEqn.H (not phi). I am not sure what other changes were needed to the main libraries that requires recompilation and use of 1.5.x. Maybe that was just for the p solve?

Hi Kent,

thanks for the link.
I think 1.5.x is needed because changes in the "MRFZone" class. There is a variant of addCoriolis() in 1.5.x which takes the density.

Regards, Hannes

MRFInterFoam
 

Hello,

kwardle, could you, please, post the link to MRFInterFoam.tgz again? The link you gave in your post on June 23, 2009, 09:25 does not work any more.

I have just started diving into openFoam to simulate a stirred tank and it would be very helpful to look at this solver. As far as I understood, it works in OF 1.5.x only. Am I right?

Thank you in advance.

included in 1.6
 

Since it was written by the guys at OpenCFD, the MRFinterFoam solver is included as part of 1.6--look in the tutorials under multiphase/interFoam/MRFinterFoam. If you want the 1.5.x version, the link from the previous post should work now.

Thank you, kwardle, very much. It works great!

unit

Hallo,


Does anyone know how to set time varying angular speed (omega)?
I'd like to simulate rotation with given acceleration. I'm using MRFInterFoam.


I couldn't find on the forum any hints but please give me a link if similar problem has been already somewhere discussed.


Thanks a lot!


Przemek

I have not done this before so I am just taking a stab, but you might take a look at this post which talks about implementing a time varying g.

http://www.cfd-online.com/Forums/ope...vity-time.html

Perhaps you would do something similar for omega--that is, define a new field for omega and write an equation to change omega within your time loop. You can look in src/finiteVolume/cfdTools/general/MRF/MRFZone.* to get an idea of how omega is set. It looks like "omega" is read from the dict as a dimensionedScalar but is immediately changed to a dimensionedVector "Omega" by multiplying with the "axis" vector. I assume "Omega" is what is used later on.

Hope this helps in some way...

MRFInterFoam
 

Very good morning kent

i am trying to solve a airpurifier by using MRFInterFoam whether it suits good or not.,
if not which solver is best to solve for it

thanking you in advance

MRFInterFoam
 

dear foamers

i am trying to solve a airpurifier by using mrfinterfoam so please, if there if any tutorial on mrfinterfoam can you please send it to me.
and also please suggest me if the solver can solve the airpurifier or not, and if not which is the best solver to solve the problem

MRFInterFoam
 

dear foamers

i am trying to solve a airpurifier by using mrfinterfoam so please, if there if any tutorial on mrfinterfoam can you please send it to me.
and also please suggest me if the solver can solve the airpurifier or not, and if not which is the best solver to solve the problem