[bvieira]

I’m trying to implement a new solver, essentially adding mesh motion capability into rhoPimpleFoam. I’m using pimpleDyMFoam as a basis to figure out what needs to be changed. My ultimate goal is to be able to run pitching airfoils on a compressible flow solver.

I was able to compile the new solver, but when I try to run, I getting a dimensions error:

  Different dimensions for -=
       dimensions : [1 0 -1 0 0 0 0] = [0 3 -1 0 0 0 0]

This seems to be a difference of density (kg/m3) in one of the sides. Tracking down where the error is coming from, I figured that it happens in the following line of my code:

// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);

Therefore, I suspect there may be something wrong in the correctPhi.H file.
I’d appreciate any thoughts on this issue. I’m attaching the correctPhi.H file, as well as the rhoPimpleDyMFoam.C file.

I'm using the latest OpenFOAM 2.1.0.

Thank you.

correctPhi.H

{
    if (mesh.changing())
    {
        forAll(U.boundaryField(), patchI)
        {
            if (U.boundaryField()[patchI].fixesValue())
            {
                U.boundaryField()[patchI].initEvaluate();
            }
        }

        forAll(U.boundaryField(), patchI)
        {
            if (U.boundaryField()[patchI].fixesValue())
            {
                U.boundaryField()[patchI].evaluate();

                phi.boundaryField()[patchI] =
                    U.boundaryField()[patchI]
                  & mesh.Sf().boundaryField()[patchI];
            }
        }
    }

    wordList pcorrTypes
    (
        p.boundaryField().size(),
        zeroGradientFvPatchScalarField::typeName
    );

    forAll(p.boundaryField(), patchI)
    {
        if (p.boundaryField()[patchI].fixesValue())
        {
            pcorrTypes[patchI] = fixedValueFvPatchScalarField::typeName;
        }
    }

    volScalarField pcorr
    (
        IOobject
        (
            "pcorr",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
        mesh,
        dimensionedScalar("pcorr", p.dimensions(), 0.0),
        pcorrTypes
    );

    while (pimple.correctNonOrthogonal())
    {
        //Changed the correction here to match pEqn.H (rhoPimpleFoam)
        fvScalarMatrix pcorrEqn
        (
            fvm::ddt(psi, pcorr)
          + fvc::div(phi)
          - fvm::laplacian(rho*rAU, pcorr)
        );

        // Removed the pressure referencing
        pcorrEqn.solve();

        if (pimple.finalNonOrthogonalIter())
        {
            phi -= pcorrEqn.flux();
        }
    }
}
//Changed to compressible version (not sure if rhoEqn.H should also be included)
#include "rhoEqn.H"
#include "compressibleContinuityErrs.H"

rhoPimpleDyMFoam.C

#include "fvCFD.H"
#include "basicPsiThermo.H"
#include "turbulenceModel.H"
#include "dynamicFvMesh.H"
#include "bound.H"
#include "pimpleControl.H"
//#include "IObasicSourceList.H"

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

int main(int argc, char *argv[])
{
    #include "setRootCase.H"
    #include "createTime.H"

    #include "createDynamicFvMesh.H"

    pimpleControl pimple(mesh);

    #include "initContinuityErrs.H"
    #include "createFields.H"
    #include "readTimeControls.H"

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

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

    while (runTime.run())
    {
        #include "readControls.H"
        #include "compressibleCourantNo.H"

        // Make the fluxes absolute
        fvc::makeAbsolute(phi, U);

        #include "setDeltaT.H"

        runTime++;

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

        mesh.update();
        
        if (mesh.changing() && correctPhi)
        {
            #include "correctPhi.H" 
        }

        // Make the fluxes relative to the mesh motion
        fvc::makeRelative(phi, U);

        if (mesh.changing() && checkMeshCourantNo)
        {
            #include "meshCourantNo.H"
        }
        
        #include "rhoEqn.H"
 
        // --- Pressure-velocity PIMPLE corrector loop
        while (pimple.loop())
        {
            #include "UEqn.H"
            #include "hEqn.H"

            // --- Pressure corrector loop
            while (pimple.correct())
            {
                #include "pEqn.H"
            }

            if (pimple.turbCorr())
            {
                turbulence->correct();
            }
        }

        runTime.write();

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

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

    return 0;
}


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

[kalle]

Interesting work. I recently wrote a low mach number combustion solver for OF2.1 capable of dynamic meshes (for adaptive refining) Here are is my correctPhi.H and another file which is called before the time loop starts.

correctPhi.H

{
    if (mesh.changing())
    {
        forAll(U.boundaryField(), patchi)
        {
            if (U.boundaryField()[patchi].fixesValue())
            {
                U.boundaryField()[patchi].initEvaluate();
            }
        }

        forAll(U.boundaryField(), patchi)
        {
            if (U.boundaryField()[patchi].fixesValue())
            {
                U.boundaryField()[patchi].evaluate();

                phi.boundaryField()[patchi] =
                U.boundaryField()[patchi] 
                & mesh.Sf().boundaryField()[patchi];
            }
        }
    }

#include "compressibleContinuityErrs.H"

    volScalarField pcorr
    (
        IOobject
        (
            "pcorr",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
        mesh,
        dimensionedScalar("pcorr", p.dimensions(), 0.0),
        pcorrTypes
    );

    dimensionedScalar rAUf("(1|A(U))", dimTime, 1.0);

    while (pimple.correctNonOrthogonal())
    {
        fvScalarMatrix pcorrEqn
        (
            fvm::laplacian(rAUf, pcorr) == fvc::div(phi) - divU
        );

        pcorrEqn.solve();

        if (pimple.finalNonOrthogonalIter())
        {
            phi -= pcorrEqn.flux();
        }
    }

#include "compressibleContinuityErrs.H"
}

And a set up file createPcorrTypes.H

   
 wordList pcorrTypes
    (
        p.boundaryField().size(),
        zeroGradientFvPatchScalarField::typeName
    );

    for (label i=0; i<p.boundaryField().size(); i++)
    {
        if (p.boundaryField()[i].fixesValue())
        {
            pcorrTypes[i] = fixedValueFvPatchScalarField::typeName;
        }
    }

Remember that this flux reconstruction is only necessary for complex mesh changes, that may flip faces. Just stretching meshes wont need this:

http://www.cfd-online.com/Forums/ope...tml#post222397

Regards,
Kalle

[bvieira]

Thank you for the reply!

I'd like to mention a few things, and I'd appreciate any comments.

1) I'm trying to model an airfoil that will be performing a pitching motion about the 1/4 chord point. My understanding is that my mesh doesn't need to deform at all, only rotate (I'm assuming the mesh motion routine will be able to translate the effects of both angle-of-attack changes, but also due to pitch rate (effective-camber effects that depend on pivot-point location)).

2) Let's say I don't need the correctPhi.H routine, should I comment it out, or it will be called only if the there are any flipping faces ?

3) About your implementation of correctPhi.H, the equation you’re using for the nonOrthogonal Mesh correction is the same as the one used in compressibleInterDyMFoam. I used the equation from rhoPimpleFoam solver under pEqn.H. Do you have comments on which one should be used ?

4) By inspecting where my error really came from, I realized that it isn't related to the correctPhi.H file. The "fvc::makeRelative(phi, U)" line gives a dimension error only after updating the mesh "mesh.update()" . If placed before it doesn't. Replacing the makeRelative command by,

fvc::makeRelative(phi, rho, U)

the error went away, but I'm not sure if that makes sense or not... Do you have any ideas?
I know that because this solver is compressible my "phi" definition is different than in an incompressible solver by a rho factor (it's using compressibleCreatePhi.H instead of CreatePhi.H)

I'd really appreciate any comments you may have on this...

[kalle]

Hi!

1,2) For such a case you could write a code that simple oscillates the complete mesh. I would guess that you can take a lot of ideas on how the AMI propeller case is done, even though I never looked at it. Just make sure you call makeAbsolute and makeRelative correctly. That approach would not flip any faces, and you can do without this correctPhi code.

3) I was using the code from interDyMFoam. Thinking more about it, is might be more relevant to look at the code in compressibleInterDyMFoam. I all cases, one should validate the approach taken by inspecting phi fields after the correctPhi.H. I did not get that far yet.

4) Have a look at the implementation of the makeRelative/Absolute. They have different methods depending on compressiblity.

Regards,
Kalle

[bvieira]

Thanks for the comments!

I was able to make the case run without any errors.
I've included rho in all MakeAbsolute and MakeRelative to be consistent and some other minor corrections and the solver started to work.

By inspecting the mesh during the motion, it seems that my mesh is rotating just as I desired, without any flipping or distorting faces.

But now, I'm facing some instability problems with the run. Residuals don't diverge, but I get instabilities in the field variables (p, rho) near the leading edge which propagate downstream and messes with the solution.
Thinking that this could be related to the new solver that I put together, I decided to run a case in rhoPimpleFoam with no moving mesh, just a fixed airfoil but also transient. I'm also getting some weird instabilities in the flow after some time. I'm still experimenting with relaxation factors and subiterations, but i may post a question on this soon. Again any ideas would be useful.

Thanks.

[fredo490]

Hello,
First I have to say that your work is really nice ! I found your topic because I also need to implement a mesh motion in rhoPimpleFoam.

I am very curious, did you succeed to solve your instability problem ? And may I ask you to share the final version of your solver ?

Thx, Fred

[prasant]

Hello All,

Did you ever success with the rhoPimpleDymFoam Solver. I had compiled the code successfully. But while running the solver, I am getting error. Could you please help me regarding this?

Regards
Prasant.

[fredo490]

How is your case with a simple rhoPimpleFoam ? Is your case stable with a static mesh ?

You can run a checkMesh at each step to check if the divergence come from the mesh Morphing.

[bieshuxuhe]

hi !
could your rhoPimpleDyMFoam be used to simulate the water ?
or it only could be used for gas computing?

thanks

xuhe

[bieshuxuhe]

hello everyone!
which do you think is easier ：
1 developing a rhoPimpleDyMFoam from rhoPimpleFoam
2 developing a compressiblePimpleDyMFoam from PimpleDyMFoam

[prasant]

Hello,

The solver is already availalble in the current version OpenFOAM
You can check the compressible solvers in OpenFOAM-2.3.0
Its well working.

Regards
Prasant.

Originally Posted by bieshuxuhe

hello everyone!
which do you think is easier ：
1 developing a rhoPimpleDyMFoam from rhoPimpleFoam
2 developing a compressiblePimpleDyMFoam from PimpleDyMFoam

[bieshuxuhe]

thanks for your reply !
I have check the compressible solvers in OpenFOAM-2.3.0 , but I didn't find the solver .
http://www.openfoam.org/docs/user/st...p#x13-890003.5
could you help me?

[prasant]

Hello,

It is available. Download OpenFOAM-2.3.0. and See the compressible solvers.
It is located in the rhoPimpleFOAM. there is a solver called rhoPimpleDymFoam.
It was implemented from OpenFOAM-2.2.1 onwards and working fine.

The page which you are seeing in the OpenFOAM website is just for information only. We need to download and see each and every thing at every latest Release.

Since this is a open source code, don't expect more information from the site.
We need to explore.....

Regards
Prasant

[bieshuxuhe]

thanks ！
I will explore !

[ChristianR1988]

Hello everybody!
I'm trying to simulate a rotating compressor-blade with rhoPimpleDyMFoam.
I use cyclicAMI and a dynamic mesh.
The solver I use is GAMG in general.
Due to the high Ma (Ma=0,95) I need to run the simulation transonic, but it doesnt work. It seems, that the pressure explodes after a little while.
First I thought its based on bad startup-conditions, so I decided to make a starting solution with transonic = no.
If I switch off the transonic-option everything works fine and I get a Solution. When I turn it back on the same problem appears again.
What can I do? Is it simply not possible or am I doing something wrong?
Thanks for your support.

Best regards, Christian.

[jvd.mechanic]

hi every body
I'm trying to solve my case with rhoPimpleDymFoam
When I write rhoPimpleDymFoam in terminal , my case solves correctly but when I want to solve by parallel situation,it give me this error in rhoPimpleDymFoam.log :

Create time

Create mesh for time = 0

Selecting dynamicFvMesh dynamicMotionSolverFvMesh
Selecting motion solver: displacementLaplacian
Selecting motion diffusion: inverseDistance

PIMPLE: Operating solver in PISO mode

Reading thermophysical properties

Selecting thermodynamics package
{
type hePsiThermo;
mixture pureMixture;
transport sutherland;
thermo hConst;
equationOfState perfectGas;
specie specie;
energy sensibleEnthalpy;
}

Reading field U

Reading/calculating face flux field phi

Creating turbulence model

Selecting turbulence model type RASModel
Selecting RAS turbulence model kEpsilon
kEpsilonCoeffs
{
Cmu 0.09;
C1 1.44;
C2 1.92;
C3 -0.33;
sigmak 1;
sigmaEps 1.3;
Prt 1;
}

Creating field dpdt

Creating field kinetic energy K

No finite volume options present

Courant Number mean: 4.35221 max: 91.4995

Starting time loop
Courant Number mean: 0.0460965 max: 1.04908
deltaT = 5.20885e-06
Time = 5.20885e-06

DICPCG: Solving for cellDisplacementx, Initial residual = 0, Final residual = 0, No Iterations 0
DICPCG: Solving for cellDisplacementy, Initial residual = 1, Final residual = 8.12622e-09, No Iterations 102
GAMG: Solving for pcorr, Initial residual = 1, Final residual = 0.00856315, No Iterations 13
diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
rhoEqn max/min : 1.50713 0.482044
smoothSolver: Solving for Ux, Initial residual = 0.0041289, Final residual = 8.5165e-08, No Iterations 2
smoothSolver: Solving for Uy, Initial residual = 0.00335382, Final residual = 1.45029e-07, No Iterations 2
smoothSolver: Solving for h, Initial residual = 0.0345707, Final residual = 7.37076e-07, No Iterations 2
GAMG: Solving for p, Initial residual = 0.00353276, Final residual = 1.39692e-08, No Iterations 1
diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 0.733222, global = -0.520471, cumulative = -0.520471
rho max/min : 2 0.5
GAMG: Solving for p, Initial residual = 0.000194151, Final residual = 9.01511e-10, No Iterations 1
diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 0.73311, global = -0.520391, cumulative = -1.04086
rho max/min : 2 0.5
GAMG: Solving for p, Initial residual = 2.18972e-06, Final residual = 2.5789e-10, No Iterations 1
diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
time step continuity errors : sum local = 0.73311, global = -0.520391, cumulative = -1.56125
rho max/min : 2 0.5
smoothSolver: Solving for epsilon, Initial residual = 0.0295855, Final residual = 9.76084e-08, No Iterations 3
smoothSolver: Solving for k, Initial residual = 0.0882087, Final residual = 1.24819e-07, No Iterations 3
ExecutionTime = 0.25 s ClockTime = 0 s
.
.
.
Courant Number mean: 0.0453958 max: 0.998186
deltaT = 5.34283e-06
Time = 8.4506e-05

DICPCG: Solving for cellDisplacementx, Initial residual = 0, Final residual = 0, No Iterations 0
DICPCG: Solving for cellDisplacementy, Initial residual = 0.101147, Final residual = 8.99411e-09, No Iterations 90
GAMG: Solving for pcorr, Initial residual = 1, Final residual = 0.00702641, No Iterations 9
diagonal: Solving for rho, Initial residual = 0, Final residual = 0, No Iterations 0
rhoEqn max/min : 2.03222 0.464785
smoothSolver: Solving for Ux, Initial residual = 0.0042863, Final residual = 9.48459e-08, No Iterations 2
smoothSolver: Solving for Uy, Initial residual = 0.00329678, Final residual = 1.42144e-07, No Iterations 2
smoothSolver: Solving for h, Initial residual = 0.0230587, Final residual = 3.37243e-07, No Iterations 2
[2] #0 Foam::error: :rintStack(Foam::Ostream&) at ??:?
[2] #1 Foam::sigFpe::sigHandler(int) at ??:?
[2] #2 in "/lib/x86_64-linux-gnu/libc.so.6"
[2] #3 Foam::hePsiThermo<Foam: :siThermo, Foam: : ureMixture<Foam::sutherlandTransport<Foam::species ::thermo<Foam::hConstThermo<Foam: : erfectGas<Foam::specie> >, Foam::sensibleEnthalpy> > > >::calculate() at ??:?
[2] #4 Foam::hePsiThermo<Foam: :siThermo, Foam: : ureMixture<Foam::sutherlandTransport<Foam::species ::thermo<Foam::hConstThermo<Foam: : erfectGas<Foam::specie> >, Foam::sensibleEnthalpy> > > >::correct() at ??:?
[2] #5
[2] at ??:?
[2] #6 __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
[2] #7
[2] at ??:?
*** Process received signal ***
Signal: Floating point exception (8)
Signal code: (-6)
Failing at address: 0x3e800000dd7
[ 0] /lib/x86_64-linux-gnu/libc.so.6(+0x370b0) [0x7f71007a70b0]
[ 1] /lib/x86_64-linux-gnu/libc.so.6(gsignal+0x37) [0x7f71007a7037]
[ 2] /lib/x86_64-linux-gnu/libc.so.6(+0x370b0) [0x7f71007a70b0]
[ 3] /opt/openfoam221/platforms/linux64GccDPOpt/lib/libfluidThermophysicalModels.so(_ZN4Foam11hePsiThe rmoINS_9psiThermoENS_11pureMixtureINS_19sutherland TransportINS_7species6thermoINS_12hConstThermoINS_ 10perfectGasINS_6specieEEEEENS_16sensibleEnthalpyE EEEEEEE9calculateEv+0x20f) [0x7f7105b513bf]
[ 4] /opt/openfoam221/platforms/linux64GccDPOpt/lib/libfluidThermophysicalModels.so(_ZN4Foam11hePsiThe rmoINS_9psiThermoENS_11pureMixtureINS_19sutherland TransportINS_7species6thermoINS_12hConstThermoINS_ 10perfectGasINS_6specieEEEEENS_16sensibleEnthalpyE EEEEEEE7correctEv+0x32) [0x7f7105b5d9f2]
[ 5] rhoPimpleDyMFoam() [0x42a05b]
[ 6] /lib/x86_64-linux-gnu/libc.so.6(__libc_start_main+0xf5) [0x7f7100791ea5]
[ 7] rhoPimpleDyMFoam() [0x42feb5]
*** End of error message ***
--------------------------------------------------------------------------
mpirun noticed that process rank 2 with PID 3543 on node jvd-K53SV exited on signal 8 (Floating point exception).

I don't understand any thing of this error.Can every body help me to solve this error?
Thanks before

[prasant]

Hello,

As per your log file, you are not using AMI approach, correct???
Then there may be a data transfer from one processor to another processor is not doing well.
check your decomposition settings. Be sure that if you are using hierarchical decomposition means, use the partitions equally in all directions. for example, 8 processors means use (2 2 2). So that data transfer between the processors will be uniform and it will not diverge.

let me know if you are still facing the issue.

Regards
Prasanth.

[jvd.mechanic]

hi prasant
Thank u for your replay.
Can u please say me what is the AMI approach? I don't know about that
my method for decomposition is simple and it's my decomposeParDict :

numberOfSubdomains 6;
method simple;
simpleCoeffs
{
n ( 3 2 1 );
delta 0.001;
}
hierarchicalCoeffs
{
n ( 3 2 1 );
delta 0.001;
order xyz;
}


thanks a lot.

[crixman]

Dear all,
I am having problems starting a rhoPimpleDyMFoam simulation.
The Error is the following

#0 Foam::error:rintStack(Foam::Ostream&) at ??:?
#1 Foam::sigFpe::sigHandler(int) at ??:?
#2 in "/lib/x86_64-linux-gnu/libc.so.6"
#3 Foam::hePsiThermo<Foam:siThermo, Foam:ureMixture<Foam::sutherlandTransport<Foam:: species::thermo<Foam::hConstThermo<Foam:erfectGa s<Foam::specie> >, Foam::sensibleEnthalpy> > > >::calculate() at ??:?




It appears at the beginning. Does anyone had a similar mistake? should I correct something in the source code?
I can provide the case if you want.
Thanks in advance!

[crixman]

((continuing the error message )

#4 Foam:siThermo::addfvMeshConstructorToTable<Foam: :hePsiThermo<Foam:siThermo, Foam:ureMixture<Foam::sutherlandTransport<Foam:: species::thermo<Foam::hConstThermo<Foam:erfectGa s<Foam::specie> >, Foam::sensibleEnthalpy> > > > >::New(Foam::fvMesh const&, Foam::word const&) at ??:?
#5 Foam::autoPtr<Foam:siThermo> Foam::basicThermo::New<Foam:siThermo>(Foam::fvMe sh const&, Foam::word const&) at ??:?
#6 Foam:siThermo::New(Foam::fvMesh const&, Foam::word const&) at ??:?
#7
at ??:?
#8 __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
#9
at ??:?
Floating point exception (core dumped)