[adambarfi]

Hi everybody,

I'm trying to modify the PisoFoam in a way that all dU/dz=0.
actually it's a method to reach a fully developed flow by exerting dU/dz=0 and setting the pressure gradient as a negative constant.

so in momentum equation I should get rid of any dU/dz, can I do something like:

            fvVectorMatrix UEqn
            (
                fvm::ddt(U)
//              + fvm::div(phi, U)
	      - fvm::div(phi, myU)
              + turbulence->divDevReff(U)
            );

            UEqn.relax();

            if (momentumPredictor)
            {
                solve(UEqn== -theGradP);
            }

which myU is the velocity vector that its z-component is zero???

the solver complied successfully by when I run it the following error appears

--> FOAM FATAL ERROR: 
incompatible fields for operation 
    [U] - [myU]

    From function checkMethod(const fvMatrix<Type>&, const fvMatrix<Type>&)
    in file /home/mostafa/OpenFOAM/OpenFOAM-2.1.0/src/finiteVolume/lnInclude/fvMatrix.C at line 1300.

FOAM aborting

I think that's because I have both U and myU in one equation, but I don't know how to fixed this problem!!!!

any hints or tips would be appreciated.

Regards,
Mostafa

[marupio]

You could try by creating an fvVectorMatrix with fvm::div(phi, U) first, then setting the z-component to zero, then bringing it into the full equation.

[adambarfi]

thanks for the reply.

I tried it before but I couldn't find anything that I can define a new divergence scheme. I found the below address but I got nothing from that!
OpenFOAM/OpenFOAM-2.2.2/src/finiteVolume/finiteVolume/divSchemes

can you please explain more about it?

Regards,
Mostafa

p.s:
should I use something like:

84 tmp<fvVectorMatrix> Smagorinsky2::divDevReff
   85 (
   86     volVectorField& U
   87 ) const
   88 {
   89     volTensorField gradU(fvc::grad(U));
   90 
   91     volSymmTensorField aniNuEff
   92     (
   93         "aniNuEff",
   94         I*nuEff() + cD2_*delta()*symm(gradU)
   95     );
   96 
   97     return
   98     (
   99       - fvm::laplacian(aniNuEff, U) - fvc::div(nuEff()*dev(T(fvc::grad(U))))
  100     );
  101 }

but where and how?

[marupio]

I don't think it's necessary to create a custom divergence scheme. What I meant was to do this all at the top-level, in your solver. Kind of like this (pseudo code):

// Create a regular matrix
fvVectorMatrix divUWithoutZ(fvm::div(phi, U));

// Set the z-component to zero
// Not sure exactly how off the top of my head, but it can be done

// Put it in your full equation
fvVectorMatrix UEqn
(
fvm::ddt(U)
+ divUWithoutZ
+ turbulence->divDevReff(U)
);

And so on. Would take some tweaking - just an idea.

[adambarfi]

Dear David,

That's a great idea. meanwhile, I need to exert this condition for all the dU/dz in the momentum eqn including the LES turbulent model, so I'm trying to combine the fvm::div(phi, U) and turbulence->divDevReff(U) into only turbulence->mydivDevReff() which will be defined in turbulenceModel and the LESModel. I hope using this method give me the correct solution, I doesn't I'll be here again

again thanks for your consideration,

Regards,
Mostafa

[adambarfi]

hmmm..., I did it but again the same error:

#0  Foam::error::printStack(Foam::Ostream&) in "/home/mostafa/OpenFOAM/OpenFOAM-2.1.0/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#1  Foam::sigSegv::sigHandler(int) in "/home/mostafa/OpenFOAM/OpenFOAM-2.1.0/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#2   in "/lib/x86_64-linux-gnu/libc.so.6"
#3  Foam::OSstream::write(Foam::word const&) in "/home/mostafa/OpenFOAM/OpenFOAM-2.1.0/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#4  Foam::operator<<(Foam::Ostream&, Foam::word const&) in "/home/mostafa/OpenFOAM/OpenFOAM-2.1.0/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#5  void Foam::checkMethod<Foam::Vector<double> >(Foam::fvMatrix<Foam::Vector<double> > const&, Foam::fvMatrix<Foam::Vector<double> > const&, char const*) in "/home/mostafa/OpenFOAM/OpenFOAM-2.1.0/platforms/linux64GccDPOpt/bin/pisoFoamCS2"
#6  
 in "/home/mostafa/OpenFOAM/OpenFOAM-2.1.0/platforms/linux64GccDPOpt/bin/pisoFoamCS2"
#7  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
#8  
 in "/home/mostafa/OpenFOAM/OpenFOAM-2.1.0/platforms/linux64GccDPOpt/bin/pisoFoamCS2"
Segmentation fault (core dumped)

the part of error which is bolded says me that again the problem is "incompatible fields for operation". for eliminating the dU/dz in the eqns I defined a myU vector which its z-component is zero throughout the grids using below code:

    volVectorField myU("myU", U);
    forAll(myU,cellI) 
    {
         scalar Ux=U[cellI].component(vector::X);
         scalar Uy=U[cellI].component(vector::Y);
         myU[cellI]=vector(Ux,Uy,0); 
    }

so in all my equation I have two separate fields: U and myU, so this error appears!

now I'm wondering how can I get rid of this component without using a new field? is there any way to turn this problem around?

Regards,
Mostafa


P.S:
the bigger problem is that I have grad(U) in my LES model, so I should exert this condition on it, too. therefore I probably should call something except U for my turbulence model and it's the biggest crisis!!!! because my model is unsteady and I have ddt(U) in my momentum eqn and again It would lead to having two different fields in one equation.

[marupio]

You posted the stack, not the error.

Don't define myU. Use the solution I suggested.

[adambarfi]

dear David,

It was the only thing I got!! No error message
What do you mean "Use the solution"? Can you illuminate it more?

[marupio]

Sorry I'm not clear. It's just an idea.

The first bit of code you posted, the top bit with "fvVectorMatrix UEqn"... replace that constructor with the bit of code I posted above (after tweaking). Don't define a myU.

As for the LES, you may have to define a custom model because you don't have access to it at the solver level. But again, do the same thing there.

[adambarfi]

thanks David,

I know what you mean and I understand where you're going, my problem is how to set the z-component to zero?
I need more time to think about it and brainstorming!

Regards,
Mostafa

[wyldckat]

Greetings to all!

At Mostafa's request, yesterday I finally took a look into this and after about 1h I didn't manage to figure it out. And unfortunately I currently do not have more time to look into this.

Nonetheless, Mostafa still asked me this after that:

is there any way to modify the div, grad, and laplacian schemes in openFoam?
I searched the src but I got nothing!

This part of OpenFOAM's source code feels almost like having to go straight to a bee's nest with no protection. It's a part of the source code that is very heavily dependant on C++ templates, which can make it quite hard to comprehend how things work. Honestly, this makes me want to ask: Mostafa, are you reaaaaally certain that you want to force "dU/dz = 0"? Why not simply run the case in pseudo-2D?

Anyway... using the User Guide as a reference, the names provided on this chapter should help us find the right keywords in the source code: http://www.openfoam.org/docs/user/fvSchemes.php
The thing is that, for example, the "div" schemes are essentially just "Gauss" using an interpolation scheme, therefore one should not look for a "div" scheme to hack/modify, but rather it's the interpolation method that you should look for.
Therefore, for example, you can look for a reference file that defined "linearUpwind", like this:

find $FOAM_SRC -name "*.H" -type f | xargs grep linearUpwind

And what should you do when you find it? Simple: you have to study all of the interpolation schemes that are present in the parent folder(s) and search for ideas on how a new interpolation scheme can be implemented. Essentially, study the complete folder of "$FOAM_SRC/finiteVolume/interpolation/surfaceInterpolation/schemes".

Beyond this, there is already a lot of studying material here: http://openfoamwiki.net/index.php/OpenFOAM_guide

Good luck! Best regards,
Bruno

[adambarfi]

Hi Bruno,

Wooooow! that's a lot of work to do, I thought that's easy to do! everything is tied together
You know, I want to solve a fully-developed turbulence flow in a channel using by a LES model. So I thought that if I can exert the dU/dx in all the equations and setting a fixed pressure gradient I can just solve a section instead of whole channel. In this way, I shouldn't wait 2 days to solve just a geometry. But it seems that solving the whole geometry would be more quicker than exerting this condition.
anyway, due to the lack of time, I'll try to modify the schemes later and I certainly going to do that. Because I need the same condition for my next project.

Thanks Bruno for the time you spent to find a solution for my problem and I appreciate it. It gives me a better prospective which how can I do that. Hope I can modify the schemes finally.

Best Regards,
Mostafa


P.S:
How can I run a pseudo-2D case?

[wyldckat]

Originally Posted by adambarfi

How can I run a pseudo-2D case?

An example is the very first tutorial case on the user guide: http://www.openfoam.org/docs/user/cavity.php#x5-40002.1