[Obad]

Hi, folks!


I am currently working on the simulation of the flow through a convergent-divergent nozzle at different back pressures as a part of a project at the university.
The program that I use for the programming is MATLAB.
I already read a lot of basic stuff about compressible flow (one dimensional flow, shock waves, etc.) and also managed to programm the flow behaviour (isentropic, inviscid) inside the nozzle using the analytical equations, including
the presence of normal shock waves.

Now I need to program the behaviour of under- and overexpanded nozzle flow, hence the reflections of the shock and expansion waves. Which is also the reason why I need some help.

I haven't had an appropriate lecture at the university yet, that covers of how to set up your own cfd simulation, so I'm pretty new to this stuff.
Nevertheless I already read about some basic cfd techniques. I read the book Fundamentals of CFD from Anderson.

So let's talk about my problem of how to simulate the under-, overexpanded jet behind the nozzle.

My idea of solving the flow field behind the nozzle is to solve the unsteady 2D Euler equations by using Maccormack's Technique (time-marching).
The mesh that I use is a simple rectangular one, since the flow field is symmetric about the axis of the nozzle one part of the mesh represents the area after the nozzle exit whereas the rest of
the mesh should represent the surroundings. See the uploaded picture to understand what I mean

So far I think I did nothing wrong, right? I figured out how to implement the mesh in my matlab program which is simply represented by a matrix and how to use the maccormack technique properly.

I think my main problem now are the boundary and initial conditions.
That's what I have tried:

The values that I want to calculate are the density, the x and y component of the velocity and the internal energy.
I initialized the nodes right at the nozzle exit with the values that I calculated in my previous program for calculating the values inside the nozzle, which are in fact the values that would exist when the nozzle
would be operated with an optimum back pressure. This values would stay constant during the calculation.
I initialized the boundary conditions of the surroundings as follows. At the left hand side of my mesh I initialized the first column of nodes with the atmospheric pressure and density, which would remain
fixed during the calculation and I set the velocity components in x and y direction to zero which I thought should change during the calculation. The nodes of the upper side of my mesh, hence the first row of my matrix
are initialized with the atmospheric values of the pressure and density and values of the x and y component of the velocity which are zero there, all the values at the upper side of the mesh are fixed values throughout
the calculation. For the lower side of the mesh, which represents the symmetry plane, I set the y component of the velocity to zero as well as the gradients of the density, the internal energy and the x component of the velocity.
Which means, that the values of the density, internal energy and the x component of the velocity become the values of the row above them, am I right?
The initial conditions that I used for the rest of the mesh are constant for every row. The values of every row correspond to the values that I set as a boundary condition at the first node at the left side of my mesh.

I also tried to describe the boundary and initial conditions in my sketch...


That is what I'm currently struggling with.

Could you please help me the set the boundary and initial conditions right. Or am I completely on the wrong track??

Cheers!

[Rami]

Over- and under- expanded jets have been studied a lot in the past, so you can find lots of papers on how to setup and solve such problems. Please note these problems are difficult to solve (even when inviscid isentropic flow is assumed) due to the resulting complex wave structure (e.g., oblique shocks, barrel shock, rarefaction waves and contact surfaces that are present, according to the pressure ratio). Your developed code should be able to cope with these (e.g., have shock capturing built into it) and you should have fine grid to resolve these features.

[Obad]

Hy Rami, thank you for your reply!

Yes there are a lot of papers on the internet, but most of them are about more advanced topics. I only want to simulate the diamond shapes, the accuracy is not so important.
I use the Maccormack technique with the equations in strong conservation form, which should be suitable for shock capturing. I already figured out which boundary conditions I have to use. I tested them by simulating my problem with OpenFoam.

Currently my problem is, that my calculation diverges. I checked if I implemented the Maccormack scheme correctly, by calculating the first iteration by hand. Everything should be ok with the implemented code. But still it diverges, I don't know why

Even if I apply my code to a simple problem like a supersonic flow through a pipe it diverges. I think that the velocity is the problem, it becomes either really high or really low (negativ values). The velocity then influences the calculation of my time step. I use the CFL criterion for the calculation of the time step.

[venkateshaero]

Hi,
I did in fluent same kind of problem..Its worked fine for me..If u need more info i can give for you. Intially i faced lot of issues.We have to set proper BC then only it will work, I am know how to do in MATLAB , If u want similar kind of program i can give that may be help you....

[Obad]

Originally Posted by venkateshaero

Hi,
I did in fluent same kind of problem..Its worked fine for me..If u need more info i can give for you. Intially i faced lot of issues.We have to set proper BC then only it will work, I am know how to do in MATLAB , If u want similar kind of program i can give that may be help you....

Hey venkateshaero!
Did you write your own code for the problem?
It would be nice if you could give me some indormation about your boundary conditions.

[Marouf Aldosari]

Hello Obad, I urgently need the Matlab code.

[Obad]

Hello,

I didn't manage to do it with matlab... Matlab is just too slow to handle such kind of calculations. Altough I only used matrices and barely for-loops Matlap just couldn't cope with it. However, I managed to do it with C++. I highly recommend you to do it in a "proper" programming language, since it's sO much faster.
It's already a while ago since I did that stuff, but if you have some questions just ask. I might be able to answer them.

Cheers!

[Suchit]

i am doing for sliding ring will you give me that MATLAB programme
on shahsuchit007@gmail.com
Thank you

[Numz]

Originally Posted by Obad

Hi, folks!


I am currently working on the simulation of the flow through a convergent-divergent nozzle at different back pressures as a part of a project at the university.
The program that I use for the programming is MATLAB.
I already read a lot of basic stuff about compressible flow (one dimensional flow, shock waves, etc.) and also managed to programm the flow behaviour (isentropic, inviscid) inside the nozzle using the analytical equations, including
the presence of normal shock waves.

Now I need to program the behaviour of under- and overexpanded nozzle flow, hence the reflections of the shock and expansion waves. Which is also the reason why I need some help.

I haven't had an appropriate lecture at the university yet, that covers of how to set up your own cfd simulation, so I'm pretty new to this stuff.
Nevertheless I already read about some basic cfd techniques. I read the book Fundamentals of CFD from Anderson.

So let's talk about my problem of how to simulate the under-, overexpanded jet behind the nozzle.

My idea of solving the flow field behind the nozzle is to solve the unsteady 2D Euler equations by using Maccormack's Technique (time-marching).
The mesh that I use is a simple rectangular one, since the flow field is symmetric about the axis of the nozzle one part of the mesh represents the area after the nozzle exit whereas the rest of
the mesh should represent the surroundings. See the uploaded picture to understand what I mean

So far I think I did nothing wrong, right? I figured out how to implement the mesh in my matlab program which is simply represented by a matrix and how to use the maccormack technique properly.

I think my main problem now are the boundary and initial conditions.
That's what I have tried:

The values that I want to calculate are the density, the x and y component of the velocity and the internal energy.
I initialized the nodes right at the nozzle exit with the values that I calculated in my previous program for calculating the values inside the nozzle, which are in fact the values that would exist when the nozzle
would be operated with an optimum back pressure. This values would stay constant during the calculation.
I initialized the boundary conditions of the surroundings as follows. At the left hand side of my mesh I initialized the first column of nodes with the atmospheric pressure and density, which would remain
fixed during the calculation and I set the velocity components in x and y direction to zero which I thought should change during the calculation. The nodes of the upper side of my mesh, hence the first row of my matrix
are initialized with the atmospheric values of the pressure and density and values of the x and y component of the velocity which are zero there, all the values at the upper side of the mesh are fixed values throughout
the calculation. For the lower side of the mesh, which represents the symmetry plane, I set the y component of the velocity to zero as well as the gradients of the density, the internal energy and the x component of the velocity.
Which means, that the values of the density, internal energy and the x component of the velocity become the values of the row above them, am I right?
The initial conditions that I used for the rest of the mesh are constant for every row. The values of every row correspond to the values that I set as a boundary condition at the first node at the left side of my mesh.

I also tried to describe the boundary and initial conditions in my sketch...


That is what I'm currently struggling with.

Could you please help me the set the boundary and initial conditions right. Or am I completely on the wrong track??

Cheers!

please give me your matlab code
nomi16425@gmail.com