[rieuk]

Hey guys, real simple question. My moving reference frame is accelerating -0.5*pi*A*T*sin(2*pi*t/T) in x-dirn and -sqrt(3)/2*pi*A*T*sin(2*pi*t/T) in y-dirn. A, T, pi are constants, t is time. I'm trying to add these to the momentum equation but it doesn't compile. Here's what I've got:

Ueqn.H

volScalarField pi = 3.1415926535897932384626433832795028841971693993751;
volVectorField acc = -pi*0.1091911*sin(2*pi*runtime.time().value()/1.091911)*(0.5*vector(1,0,0)+pow(3,0.5)*vector(0,1,0));
    
fvVectorMatrix UEqn
    (
        fvm::ddt(U)
      + fvm::div(phi, U)
      - fvm::laplacian(nu, U)
      - acc
    );

    if (momentumPredictor)
    {
        solve(UEqn == -fvc::grad(p));
    }

What's wrong with it??

[rieuk]

Ok....put it this way,

1. How can I create an empty volScalarField with dimensions acceleration?
2. How do I introduce the runtime class into Ueqn.H? (What header file do I need to introduce)

I couldn't find any of this basic stuff in the programmer's guide...

[ozgur]

Hi Pavan,

I am not expert but if I make some comments on your code, I can say that:

1. I wouldn't define pi as volScalarField, because it is a constant scalar, not a field variable which has different values in time and space. So I think,

scalar pi 3.1415926535897932;

would be sufficient.

There might (/or not) also be an already defined macro for the number pi in OpenFOAM but I haven't used it before.

2. I wouldn't define acc as volVectorField neither (if I won't need to use it over and over, or do other calculations with it). Because, if you want to define a volVectorField, then (as far as I understood up to now) then also you will need to construct an IOobjet and make mesh association, defining the dimensions etc., such as:

    volVectorField acc
    (
        IOobject
        (
            "acc",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::AUTO_WRITE
        ),
        mesh
    );

So the code gets more complicated, extra memory is used etc, etc...

What can be done may be writing the acc formula directly into the Ueqn.

However, I am not sure about how will the expressions like 0.5*vector(1,0,0) compiled. They might need to be verified separately.

3. I haven't used icoDyMFoam so I don't know anything about it, however, I have been using interDyMFoam for a while and I guess the "DyM" stands for dynamic mesh (and it maybe the same in icoDyMFoam). In interDyMFoam the mesh movement is imposed by an additional utility called gen6DoF, which is in turn equivalent to an extra acceleration term to UEqn. Since your acceleration term is also only a function of time, and is independent of space, use of gen6DoF explicitly might also be considered.

Regards,

Ozgur

[Simon Lapointe]

Hi Pavan,

Since your source term is constant in space you don't need to define a volVectorField, a simple vector would do it.

You could try something like:

dimensionedVector accX
(
"accX",dimensionSet(0, 1, -2, 0, 0),vector(1, 0, 0)
);

dimensionedVector accY
(
"accY",dimensionSet(0, 1, -2, 0, 0),vector(0, 1, 0)
);

scalar Pi = mathematicalConstant::pi;

dimensionedVector acc = -Pi*0.1091911*sin(2*Pi*runTime.value()/1.091911)*(0.5*accX+pow(3,0.5)*accY);

[rieuk]

Thanks, Ozgur!

1. Agreed about the scalar pi thing.

2. I don't know I was referring to this thread: http://www.cfd-online.com/Forums/ope...nce-frame.html (3rd post by Henry) where Fcent is defined as a volVectorField without being a field variable.

3. Yeah icoDyMFoam is just the dynamic mesh version of icoFoam, the incompressible fluid solver with no turbulence.

[rieuk]

Also, thanks Simon! Will try this asap

[rieuk]

Well, with the extra help of an experienced OF user I got the following to compile:

scalar Pi = mathematicalConstant::pi;
dimensionedScalar accCt("accCt",dimensionSet(0,1,-2,0,0,0,0),10000*Foam::cos(2*Pi*runTime.value()/0.5));
dimensionedVector acc("acc",dimensionSet(0,1,-2,0,0,0,0), vector::zero);

acc.component(vector::X) = 0.5 * accCt;
acc.component(vector::Y) = 0.5 * Foam::sqrt(3.0) * accCt;

        fvVectorMatrix UEqn
        (
            fvm::ddt(U)
          + fvm::div(phi, U)
          - fvm::laplacian(nu, U)
        );


        solve(UEqn == acc-fvc::grad(p));

However the acceleration term added does not seem to have any effect whatsoever on the solution. I just start with a stationary 2d cylinder in a rectangular domain with no flow. Solution stays zero the whole time even with the sinusoidally varying acceleration term -- why?? Is there some part of the PISO loop I have to edit?

[rieuk]

I notice in dnsFoam the force term is implemented like so:

        fvVectorMatrix UEqn
        (
            fvm::ddt(U)
          + fvm::div(phi, U)
          - fvm::laplacian(nu, U)
    ==
      acc
        );

        solve(UEqn == -fvc::grad(p));

As opposed to:

        fvVectorMatrix UEqn
        (
            fvm::ddt(U)
          + fvm::div(phi, U)
          - fvm::laplacian(nu, U)
        );


        solve(UEqn == (acc-fvc::grad(p)));

What is the difference? They give vastly differing results...

[su_junwei]

They are different.
In the first, acc is introduced to momentum equation, and in the second acc was introduced into pressure equation and pressure equation should also be altered accordingly (if you don't include acc relative term in pressure, the results should differ much). If acc has effect on the velocity, doing in the second way is more implicit.
please see derivation of pressure equation in Henrik Rusche's phd thesis.

regards, Junwei

[rieuk]

Thanks for the reply Su Junwei! So I take it the first approach is valid but not as implicit as the second approach (with modification to the pressure equation)?