CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM

Creating a solver for the nondimensionalized Navier-Stokes equation

Register Blogs Community New Posts Updated Threads Search

Like Tree12Likes

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   November 15, 2012, 10:52
Question Creating a solver for the nondimensionalized Navier-Stokes equation
  #1
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
Hi everyone,

I have recently started using OpenFoam.
I want to resolve the nondimensionalized Navier-Stakes equation (in which appears the Reynolds number).

My equation is ( you can also find it attached ) : Re . (U . div) U + grad(meanP) - laplacian(U) = - grad(fluctP)


meanP is a unit vector which is constant and P=meanP+fluctP
Re is the reynolds number



I want to resolve my equation for U and fluctP


For this, I have modified the IcoFoam solver as follows, but it seems that there is an error ( may be with the '' meanP '' )



fvVectorMatrix UEqn
(
fvm::ddt(U)
+ Re*(fvm::div(phi, U))
- fvm::laplacian(U)
+ magSqr(meanP)
);

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


I will be thankfull If someone can help me to create this solver or share a web site which treated the implementation of the nondimensionalized NS equation.


Best regards

Mehrez
Attached Files
File Type: gz equations.bmp.tar.gz (2.9 KB, 33 views)

Last edited by Mehrez; November 16, 2012 at 05:04.
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 04:01
Default
  #2
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
Hi

Can someone help me please.

Thank you
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 04:49
Default
  #3
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
please find attached the OpenFoam files of the solver (I have proceeded by modifying the IcoFoam solver to myIcoFoamB solver).
Attached Files
File Type: gz myIcoFoamB.tar.gz (2.1 KB, 162 views)
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 05:39
Default
  #4
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Hi,

there are many errors:

Code:
fvVectorMatrix UEqn
          (
              fvm::ddt(U)
            + Re*(fvm::div(phi, U))
            - fvm::laplacian(U)
            + magSqr(meanP)
          );
First: laplacian(U) is not a class. If you use laplacian you have to a coefficent like:

Code:
laplacian(nu, U)
You set Re as reynoldsnumber. But the solver do not know what "Re" is. Its a none declared variable. The same like magSqr(meanP)

OF do not know what meanP is.
So you have to define

Re
meanP

and the other variables you implemented!

The following example is working:
Code:

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

        solve(UEqn == -fvc::grad(p));
Tobi is offline   Reply With Quote

Old   November 16, 2012, 05:50
Default
  #5
Senior Member
 
Hisham's Avatar
 
Hisham Elsafti
Join Date: Apr 2011
Location: Braunschweig, Germany
Posts: 257
Blog Entries: 10
Rep Power: 17
Hisham is on a distinguished road
Hi Mehrez,

First it would be better if you modify the Make/files as in the user guide chapter 3 to be:
Code:
myIcoFoamB.C

EXE = $(FOAM_USER_APPBIN)/myIcoFoamB
When I try to compile your solver (initially) two errors occur from different places. This points out that you have made a punch of modifications then decided to compile ... Not a good idea .. It is always a good practice to compile regularly after small code changes to keep track of code bugs early on. It is a little bit tedious but very efficient. Back to the errors:

Code:
createFields.H: In function ‘int main(int, char**)’:
createFields.H:23:5: error: no matching function for call to ‘Foam::GeometricField<Foam::Vector<double>, Foam::fvPatchField\
, Foam::volMesh>::GeometricField(Foam::ITstream&)’
createFields.H:23:5: note: candidates are:
Then the compiler gives you some ideas about what could be the problem. Anyway the code that caused this error is:
Code:
volVectorField gradP
    (
        transportProperties.lookup("gradP")
    );
Essentially this is not how you declare a volVectorField. The compiler is complaining that there is no lookup function for the volVectorField and if there were one (which isn't the case) you are doing it wrong and hence the error. volVectorField are defined like the p if you'll read as user input but if you need to calculate:
Code:
volVectorField gradP
    (
        IOobject
        (
            "gradP",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ,
            IOobject::NO_WRITE
        ),
    fvc::grad(p)        
    );
The second error:
Code:
myIcoFoamB.C:61:18: error: no match for ‘operator+’ in ‘Foam::operator-(const Foam::tmp<Foam::fvMatrix<Type> >&, const Foam\
::tmp<Foam::fvMatrix<Type> >&) [with Type = Foam::Vector<double>]((*(const Foam::tmp<Foam::fvMatrix<Foam::Vector<double> > \
>*)(& Foam::fvm::laplacian(const Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>&) [with Type = Foam::Vector<\
double>]()))) + Foam::magSqr(const Foam::GeometricField<Type, PatchField, GeoMesh>&) [with Type = Foam::Vector<double>, Pat\
chField = Foam::fvPatchField, GeoMesh = Foam::volMesh]()’
myIcoFoamB.C:61:18: note: candidates are:
Comes from the eq:
Code:
fvVectorMatrix UEqn
        (
            fvm::ddt(U)
          + Re*(fvm::div(phi, U))
          - fvm::laplacian(U)
	  + magSqr(gradP)
        );
Maybe it will be solved if you fix the first one or not ...

Best regards,
Hisham
Tobi likes this.
Hisham is offline   Reply With Quote

Old   November 16, 2012, 06:34
Default
  #6
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Quote:
Originally Posted by Hisham View Post
Hi Mehrez,
Code:
fvVectorMatrix UEqn
        (
            fvm::ddt(U)
          + Re*(fvm::div(phi, U))
          - fvm::laplacian(U)
      + magSqr(gradP)
        );
Maybe it will be solved if you fix the first one or not ...

Best regards,
Hisham
Hi Hisham, the laplacian class is wrong like I posted befor.
I just had a very short look into the code and dont see the declaration of the other variables

Well nice replay!
Code:
        (
            fvm::ddt(U)
          + Re*(fvm::div(phi, U))
          - fvm::laplacian(nu, U)
      + magSqr(gradP)
        );
Tobi is offline   Reply With Quote

Old   November 16, 2012, 06:43
Default
  #7
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
Hi guys

Thank you for your help.

I think that it is better if you take a look to my equation (which I have attached in my first message). You can see that I don't have : nu . laplacian (U) but I just have : laplacian (U)

Concerning the declaration of the variables : I have declared " Re " and " gradP " in the file createFields.H

dimensionedScalar Re
(
transportProperties.lookup("Re")
);

volVectorField gradP
(
transportProperties.lookup("gradP")
);
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 06:50
Default
  #8
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Oh sorry,

the laplacian(U) class existis! So its working
Mehrez likes this.
Tobi is offline   Reply With Quote

Old   November 16, 2012, 06:58
Default
  #9
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
Hi Hisham

Thank you for your answer.

I think that I don't have to declare " gradP " like this :

volVectorField gradP ( IOobject ( "gradP", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE ), fvc::grad(p) );

As I said before, the '' gradP '' is a unit vector which is constant and which I put it like an input (so this is why OF don't have to compute it).
I have just to declare it like a vector which I will specify the value in constantProperties file.
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 07:01
Default
  #10
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Quote:
Originally Posted by Mehrez View Post
Hi Hisham

Thank you for your answer.

I think that I don't have to declare " gradP " like this :

volVectorField gradP ( IOobject ( "gradP", runTime.timeName(), mesh, IOobject::NO_READ, IOobject::NO_WRITE ), fvc::grad(p) );

As I said before, the '' gradP '' is a unit vector which constant and which I put it like an input (so this is why OF don't have to compute it).
I have just to declare it like a vector which I will specify the value in constantProperties file.
Hmmm ... you calculate your gradP field with grad(p) or?
But how whould you get the information of the "vector" couse the pressure is just a scalar without a direction.
Tobi is offline   Reply With Quote

Old   November 16, 2012, 07:06
Default
  #11
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
Hi Tobi

As you can see in my equation gradP is the gradient of the nondimensionalized mean pressure. I have it like an input (I will study my flow in function of given values of Re and gradP)

I think that the problem is how to define this vector (gradP) in the equation
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 07:15
Default
  #12
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
I think that it is more clear like this

Equation : Re . (U . div) U + gradP - laplacian(U) = - grad(fluctP)

OF syntax :

fvVectorMatrix UEqn
(
fvm::ddt(U)
+ Re*(fvm::div(phi, U))
- fvm::laplacian(U)
+ magSqr(gradP)
);

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


Input :
Re : scalar
gradP : vector

output (after computation)
U : vector field
fluctP : vector field
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 10:22
Default
  #13
Senior Member
 
Bernhard
Join Date: Sep 2009
Location: Delft
Posts: 790
Rep Power: 21
Bernhard is on a distinguished road
If gradP is a constant vector, you should not store it in a volVectorField, but in a dimensionedVector from the transportproperties.

Two more things:
- Are you sure Reynolds should appear like you wrote it. As far as I know, in the dimensionless NS equations, you have 1/Re in the diffusive terms.
- You're now adding magSqr(gradP), which is a scalar, to a vector equation. This can not be correct.
Mehrez likes this.
Bernhard is offline   Reply With Quote

Old   November 16, 2012, 10:52
Default
  #14
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
Dear Bernhard,

I'm sure that I have a right equation (attached here)

Thank you very much for your help. Actually I can compile my solver with your corrections.

fvVectorMatrix UEqn
(
fvm::ddt(U)
Re*(fvm::div(phi, U))
- fvm::laplacian(U)
+ gradP
);

and :

dimensionedVector gradP
(
transportProperties.lookup("gradP")
);


To test this solver, I've taken the attached example "cavity" after putting U, p, Re, and gradP in a dimensionless form. I've done the same with the domain (I've commented the line " conversion to meters " ).

Executing "cavity" with "myIcoFoamB" solver, OpenFoam returns an error message on the dimensions:

-> FOAM FATAL ERROR:
incompatible dimensions for operation
[U [0 -1 0 0 0 0 0]] - [U [0 -2 0 0 0 0 0]]

I think it comes from the operators (grad gives : L ^ -1 and Laplacian gives : L ^ -2). Is there a way to put the components (x, y, z) of the domain in a dimensionless form in order to get dimensionless values ​​by applying one of the operators?


Thank you very much for your help.

Best regards.

Mehrez
Attached Files
File Type: zip equations.bmp.zip (3.0 KB, 11 views)
File Type: zip solveur+cavity.zip (31.6 KB, 15 views)
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 11:24
Default
  #15
Senior Member
 
Bernhard
Join Date: Sep 2009
Location: Delft
Posts: 790
Rep Power: 21
Bernhard is on a distinguished road
Complicated one, you can choose to ignore the dimension-checking. I think you can do this in the controlDict you find in ..../OpenFOAM-2.1.x/etc/ .
Maybe you can put it in your own controlDict as well. I never tested this, but you can look here.
Mehrez likes this.
Bernhard is offline   Reply With Quote

Old   November 16, 2012, 11:31
Default
  #16
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Hi,

you can put a dimensionSet on your p vector that the dimension check is working!

Tobi
Mehrez likes this.
Tobi is offline   Reply With Quote

Old   November 16, 2012, 11:43
Default
  #17
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
Hi guys,

Again thanks for your precious help.

@ Bernhard : I have found the file " controlDict " , can you please precise me how to proceed.

@ Tobi : I've tried what you said but it doesn't work. it returns each time a different dimensions error ! I can't understand

Mehrez
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 11:50
Default
  #18
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
Post your dimenionSet and the output error.
Tobi is offline   Reply With Quote

Old   November 16, 2012, 11:56
Default
  #19
Member
 
Sami
Join Date: Nov 2012
Location: Cap Town, South Africa
Posts: 87
Rep Power: 13
Mehrez is on a distinguished road
Re Re [ 0 0 0 0 0 0 0 ] 0.01;
gradP gradP [ 0 0 0 0 0 0 0 ] (0.9 0.2 0);
U : dimensions [0 0 0 0 0 0 0];
p : dimensions [0 0 0 0 0 0 0];

In the blockMeshDict, I have commented the following line :
//convertToMeters 0.1;

and this is what I get :






Create time

Create mesh for time = 0

Reading transportProperties

Reading field p

Reading field U

Reading/calculating face flux field phi


Starting time loop

Time = 0.005

Courant Number mean: 0 max: 0


--> FOAM FATAL ERROR:
incompatible dimensions for operation
[U[0 -1 0 0 0 0 0] ] - [U[0 -2 0 0 0 0 0] ]

From function checkMethod(const fvMatrix<Type>&, const fvMatrix<Type>&)
in file /opt/openfoam211/src/finiteVolume/lnInclude/fvMatrix.C at line 1316.

FOAM aborting

#0 Foam::error:rintStack(Foam::Ostream&) in "/opt/openfoam211/platforms/linuxGccDPOpt/lib/libOpenFOAM.so"
#1 Foam::error::abort() in "/opt/openfoam211/platforms/linuxGccDPOpt/lib/libOpenFOAM.so"
#2
in "/opt/openfoam211/platforms/linuxGccDPOpt/bin/myIcoFoamB"
#3
in "/opt/openfoam211/platforms/linuxGccDPOpt/bin/myIcoFoamB"
#4
in "/opt/openfoam211/platforms/linuxGccDPOpt/bin/myIcoFoamB"
#5 __libc_start_main in "/lib/i386-linux-gnu/libc.so.6"
#6
in "/opt/openfoam211/platforms/linuxGccDPOpt/bin/myIcoFoamB"
Aborted (core dumped)
ubuntu@ubuntu-VirtualBox:~/OpenFOAM/mehrez-2.1.1/test/cavity$
Mehrez is offline   Reply With Quote

Old   November 16, 2012, 11:58
Default
  #20
Super Moderator
 
Tobi's Avatar
 
Tobias Holzmann
Join Date: Oct 2010
Location: Tussenhausen
Posts: 2,708
Blog Entries: 6
Rep Power: 51
Tobi has a spectacular aura aboutTobi has a spectacular aura aboutTobi has a spectacular aura about
Send a message via ICQ to Tobi Send a message via Skype™ to Tobi
What is that:

U : dimensions [0 0 0 0 0 0 0];
p : dimensions [0 0 0 0 0 0 0];


??? I dont understand the two lines.
Set

pGrad pGrad [0 -1 0 0 0 0 0]
Tobi is offline   Reply With Quote

Reply

Tags
navier-stokes solver, openfoam


Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are Off
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
suGWFoam: Richards equation solver for porous media flows liu OpenFOAM Announcements from Other Sources 4 February 10, 2021 15:14
Navier Stokes solver Siddharth Main CFD Forum 2 September 13, 2007 01:06
test prob for 2D unsteady navier stokes equation Shah Main CFD Forum 5 April 20, 2007 07:25
Navier Stokes Solver Khan Main CFD Forum 2 December 12, 2006 09:41
help: I am trying to solve Navier Stokes compressible and viscid flow Jose Choy Main CFD Forum 2 May 18, 2000 05:45


All times are GMT -4. The time now is 13:19.