Solving Heavy ODE System to Compute Round Jet Near Field

Clausius2
Science Advisor
Gold Member
Messages
1,433
Reaction score
7
In order to solve the near field description of a round jet, I want to work out the variables F(\eta), \rho(\eta) and Y(\eta) which represents the self similar stream function, density, and mass fraction respectively. The system obtained is:

\Big(\frac{F'}{\rho}\Big)''+\frac{F}{2}\Big(\frac{F'}{\rho}\Big)'=0 (1)

(\rho Y')'+\frac{F}{2}Y'=0 (2)

besides a function \rho=\rho(Y) which is known previously.

Boundary conditions are:

1) \eta\rightarrow+\infty; \frac{F'}{\rho}\rightarrow 0; Y \rightarrow 0;

2) \eta\rightarrow-\infty; \frac{F'}{\rho}\rightarrow 1; Y \rightarrow 1; F\rightarrow \eta;

The first question I have is how can I transform (1) into a system of three first order ordinary differential equations. I have done it yet before with Blasius type equations, but here the density makes it a bit difficult. The aim of my question is to compute both coupled equations with a Non Linear Shooting Method.

Thanks in advance.
 
Physics news on Phys.org
You can't do it analitically.It's like \int e^{\sin x} \ dx...Can"t find the sollution.I think there are numerical algorithms in which you can find the values you're looking for...

Daniel.
 
Dude, that's so interesting. I tell you what, I'm going to work on the easier one (that's a challenge for me):

f^{'''}+ff^{''}=0 with (0<\eta<+\infty)

f(0)=f^{'}(0)=0

f^{'}(+\infty)=1

I'll attempt to perfect a Runge-Kutta/ "Shooting Method". May take a few days but I'm a good programmer. Then I'll look at yours.
 
saltydog said:
Dude, that's so interesting. I tell you what, I'm going to work on the easier one (that's a challenge for me):

f^{'''}+ff^{''}=0 with (0<\eta<+\infty)

f(0)=f^{'}(0)=0

f^{'}(+\infty)=1

I'll attempt to perfect a Runge-Kutta/ "Shooting Method". May take a few days but I'm a good programmer. Then I'll look at yours.

Alright, maybe not a few days, but that's only because I spent all that time on Runge Kutta and those IDEs a few weeks ago. Anyway, the objective here is to find f^{''}(0) via "shooting method". I used simple interpolation to zero into the value of the second derivative at \eta=0. The value I got to 4 places is 0.4696. Attached is a plot of the LHS of the ODE in the range (0,20). The results were no larger than 10^-6. This gives me confidence that the numerical results are ok.
 

Attachments

  • Blasius1.JPG
    Blasius1.JPG
    5.8 KB · Views: 514
saltydog said:
Dude, that's so interesting. I tell you what, I'm going to work on the easier one (that's a challenge for me):

f^{'''}+\frac{ff^{''}}{2}=0 with (0<\eta<+\infty)

f(0)=f^{'}(0)=0

f^{'}(+\infty)=1

I'll attempt to perfect a Runge-Kutta/ "Shooting Method". May take a few days but I'm a good programmer. Then I'll look at yours.

That's a Blasius type problem. The fact is setting up a shooting method for this equation is a bit easier. I made it some time ago. But the problem with my equations is the coupling via density, and reducing it to five first order ODEs.
 
Clausius2 said:
In order to solve the near field description of a round jet, I want to work out the variables F(\eta), \rho(\eta) and Y(\eta) which represents the self similar stream function, density, and mass fraction respectively. The system obtained is:

\Big(\frac{F'}{\rho}\Big)''+\frac{F}{2}\Big(\frac{F'}{\rho}\Big)'=0 (1)

(\rho Y')'+\frac{F}{2}Y'=0 (2)

besides a function \rho=\rho(Y) which is known previously.

Clausius, can you kindly clear up something for me as I'm not familiar with the application. Are you saying that rho is a function of Y, that is \rho=g(Y(\eta)) and that this function is known? If that' so, then I can derive the following 5 equations:

F^{'}=vg(Y)

v^{'}=r

r^{'}+\frac{F}{2}r=0

Y^{'}=s

s^{'}+\frac{s^{2}g^{'}(Y)}{g(Y)}+\frac{Fs}{2g(Y)}=0

where g^{'}(Y) is the derivative with respect to Y.
 
saltydog said:
Clausius, can you kindly clear up something for me as I'm not familiar with the application. Are you saying that rho is a function of Y, that is \rho=g(Y(\eta)) and that this function is known? .

Yes it is known. Anyway, I have reduced it yet to a 5 first order ODEs. The aim of this is to solve each of them by a Runge-Kutta method, shooting from
-infinite, and converging to the shooting parameters via Newton-Raphson method.

I am having serious troubles at setting up the process, because Newton-Raphson formulae involves derivative respect to shooting parameters in the denominators.

Some web-link or reference to examples of this type of non linear problems solved by means of shooting method would be greatly appreciated.
 
Clausius2 said:
Yes it is known. Anyway, I have reduced it yet to a 5 first order ODEs. The aim of this is to solve each of them by a Runge-Kutta method, shooting from
-infinite, and converging to the shooting parameters via Newton-Raphson method.

I am having serious troubles at setting up the process, because Newton-Raphson formulae involves derivative respect to shooting parameters in the denominators.

Some web-link or reference to examples of this type of non linear problems solved by means of shooting method would be greatly appreciated.

Clausius, would you please report the five equations here so that I and others can take a look at them.

Thanks,
Salty
 
y_1=F
y_2=F'/\rho
y_3=(F'/\rho)'
y_4=Y
y_5=\rho Y

which yields:

y_1'=y_2\rho
y_2'=y_3
y_3'=-0.5y_1y_3
y_4'=y_5/\rho
y_5'=-y_1 y_5/(2\rho)
 
  • #10
Clausius, you said rho is known. My understanding from you is that it is some functional dependency on Y. Can you supply this functional dependence?

I'd like to begin working on the problem but can't until you define what rho is in terms of the function Y.
 
  • #11
saltydog said:
Clausius, you said rho is known. My understanding from you is that it is some functional dependency on Y. Can you supply this functional dependence?

I'd like to begin working on the problem but can't until you define what rho is in terms of the function Y.

Sorry.

\rho=\frac{1}{Y+\epsilon(1-Y)}

with

\epsilon<<1

I'm using Matlab. If you try it, please post me any incidence you have. Actually I am having problems with shooting parameters.
 
  • #12
Clausius2 said:
Sorry.

\rho=\frac{1}{Y+\epsilon(1-Y)}

with

\epsilon<<1

I'm using Matlab. If you try it, please post me any incidence you have. Actually I am having problems with shooting parameters.

Fantastic! Two equation in two unknowns. You know, I'm optimistic there are people in here that can do this one especially now that you've posted rho. I do intend to spend time with it but I'm kinda busy with work now. I'm patient; I spent 12 months on another one some time ago. I suppose you don't have that time though.
 
  • #13
saltydog said:
Fantastic! Two equation in two unknowns. You know, I'm optimistic there are people in here that can do this one especially now that you've posted rho. I do intend to spend time with it but I'm kinda busy with work now. I'm patient; I spent 12 months on another one some time ago. I suppose you don't have that time though.

Till Friday.
 
  • #14
Some qualitative progress:

Consider:

\rho=\frac{1}{Y+\epsilon(1-Y)}

And the boundary conditions:

1) \eta\rightarrow+\infty; \frac{F'}{\rho}\rightarrow 0; Y \rightarrow 0;

2) \eta\rightarrow-\infty; \frac{F'}{\rho}\rightarrow 1; Y \rightarrow 1; F\rightarrow \eta;

We can deduce the values of \rho at these boundaries:

\lim_{n\rightarrow+\infty}\rho=\frac{1}{\epsilon}

\lim_{n\rightarrow-\infty}\rho=1

Thus:

\lim_{n\rightarrow+\infty}F^{'}=0

\lim_{n\rightarrow-\infty}F^{'}=1

Am I interpreting the limits correctly?

If so, then I can begin to describe, qualitatively at least, the behavior of F and Y as per the attached plot (or some variant thereof). Note that Y is shown to be a sigmodal type function.

What do you think Clausius? Does the sigmoid function have any relevance in the application you're working on?

I plan to begin numerical work from this qualitative perspective using five first-order ODEs and the "shooting method" to search 3-D space for three of the five initial conditions, holding F(0) and Y(0) constant (at first). I'll attempt to set up a "quasi-random" search of the space to try and zero-in on a solution. May take a while . . .
 

Attachments

  • pf5.JPG
    pf5.JPG
    6.4 KB · Views: 449
  • #15
Clauius, I don't know where you at, but this is what I'm going with, "win-loose-or draw":

F^{'}=\frac{v}{y+\epsilon(1-y)}

Y^{'}=s

r^{'}=-\frac{F}{2}r

v^{'}=r

s^{'}=-s[Y+\epsilon(1-Y)][\frac{s(\epsilon-1)}{(y+\epsilon(1-y))^2}}+\frac{F}{2}]

with:

F(0)=a
F^{'}(0)=c
F^{''}(0)=k<0
Y(0)=b
Y^{'}(0)=d<0

and:

r(0)=k(b+\epsilon(1-b))-c(\epsilon d-d)
v(0)=c(b+\epsilon(1-b))
s(0)=d

Note that F^{''}(0)<0 and Y^{'}(0)<0 in keeping with the qualitative results stated earlier.

So, I'll keep F(0) and y(0) constant, (with small \epsilon), and search a 3-d space of F'(0), F''(0), and y'(0) using a Halton sequence (quasi-random): Do a Runge-Kutta at each point, check the boundary at some max point, find the closest match. If something pans-out, restrict the search around that point. The attached is the Halton sequence distribution.

You know, if anybody else has a better idea I'd be interested.
 

Attachments

  • halton.JPG
    halton.JPG
    20.5 KB · Views: 414
  • #16
Just a correction:

The expression for s^{'} above is incorrect. It should be:

s^{'}=\frac{\epsilon s^2-s^2-0.5 s v \rho}{\rho+0.5 F \rho^{2}}

I've coded the algorithm to search for solutions that I stated above. My initial results, although not successful, are encouraging as some numerical results seem to tending to the qualitative results stated earlier. I plan to continue working with it.

Clausisus, I hope I'm not distracting you in regards to your application. You know, some people work cross-word puzzles for fun. I seem to have chosen DE's for some reason . . .
 
  • #17
saltydog said:
Just a correction:

The expression for s^{'} above is incorrect. It should be:

s^{'}=\frac{\epsilon s^2-s^2-0.5 s v \rho}{\rho+0.5 F \rho^{2}}

I've coded the algorithm to search for solutions that I stated above. My initial results, although not successful, are encouraging as some numerical results seem to tending to the qualitative results stated earlier. I plan to continue working with it.

Clausisus, I hope I'm not distracting you in regards to your application. You know, some people work cross-word puzzles for fun. I seem to have chosen DE's for some reason . . .

THANKS for your help, saltydog.

I have solved the problem. If you want to know how to do about it, let me know. Hey! it seems this thread has been put upside down, now I am who is offering aid! :smile: .
 
  • #18
Clausius2 said:
THANKS for your help, saltydog.

I have solved the problem. If you want to know how to do about it, let me know. Hey! it seems this thread has been put upside down, now I am who is offering aid! :smile: .

Of course I'd like to know how you solved it. Please give me the details. My understanding is that the initial conditions F(0) and Y(0) needed to be known to solve it for a practical application. Is this not the case?

Thanks a bunch,
Salty
 
  • #19
saltydog said:
Of course I'd like to know how you solved it. Please give me the details. My understanding is that the initial conditions F(0) and Y(0) needed to be known to solve it for a practical application. Is this not the case?

Thanks a bunch,
Salty

What do you mean? I don't need to know the values of F and Y at \eta=0 in order to solve the problem. It is closed with the boundary conditions I posted.

Do you ever have implemented a code of shooting method? If so, let me know.
 
  • #20
Clausius2 said:
What do you mean? I don't need to know the values of F and Y at \eta=0 in order to solve the problem. It is closed with the boundary conditions I posted.

Do you ever have implemented a code of shooting method? If so, let me know.

Ok Clausius, I'm beginning to understand your 5 equations are nicer. However, I think you made a typo for y_5. It should be:

y_5=\rho Y^{'}

Right?

I do not understand how the problem is closed. I have the Mathematica code implemented to study five equations and is a simple matter to switch the derivatives for each one. I'll spend more time with it and your equations.
 
  • #21
yes, I made a mistake. You're right about y5.

The problem is closed. Go to my first thread and take a look at the boundary conditions. There are 5 boundary conditions for a system of equations of third order (the first one) and second order (the second one). You do not need any additional condition.
 
  • #22
Clausius2 said:
yes, I made a mistake. You're right about y5.

The problem is closed. Go to my first thread and take a look at the boundary conditions. There are 5 boundary conditions for a system of equations of third order (the first one) and second order (the second one). You do not need any additional condition.

Ok Clausius. I'm weak with boundary value problems like this (well, other stuf too but I digress). But I'm ok with that. This is my plan: Using my Mathematica program, I'm just going to search 5-space using a 5-D quasi-random distribution and attempt to zero-in to the "initial condition" which tends to the boundary conditions specified. I mean, I have a 2 GHz machine and I'll just let it run. Can you at least tell me if that is an appropriate means of attack?

Thanks,
Salty
 
  • #23
saltydog said:
Ok Clausius. I'm weak with boundary value problems like this (well, other stuf too but I digress). But I'm ok with that. This is my plan: Using my Mathematica program, I'm just going to search 5-space using a 5-D quasi-random distribution and attempt to zero-in to the "initial condition" which tends to the boundary conditions specified. I mean, I have a 2 GHz machine and I'll just let it run. Can you at least tell me if that is an appropriate means of attack?

Thanks,
Salty


Hummmmm . . . my initial results are looking encouraging . . . got to go though. Will get back to it tonight.
 
  • #24
Well, I have an answer but I don't like it:


F(\eta)=\left\{\begin{array}{cc}F_1(\eta),&\mbox{ if } \eta>0\\F_2(\eta),&\mbox{ if }\eta<0\end{array}\right

That is, it's discontinuous at the origin.

Same thing for Y(\eta). Back-substitution of the results into the ODEs produces low discrepancy from zero but I'm not impressed with it.

Anyway, if nobody wants to summarize a better solution, I'll write this one up in a few days.
 
  • #25
Not, the solution is not discountinous at all. Don't forget what is the origin of the problem. It cannot be a discontinuity anywhere.

Here is a plot of the solution F'/rho (solid line) and Y (dashed line). I hope you could see it.
 

Attachments

  • Capa1.jpeg
    Capa1.jpeg
    5.8 KB · Views: 495
  • #26
Clausius2 said:
Not, the solution is not discountinous at all. Don't forget what is the origin of the problem. It cannot be a discontinuity anywhere.

Here is a plot of the solution F'/rho (solid line) and Y (dashed line). I hope you could see it.

Thanks Clausius. What value of \epsilon did you use?

Also, what does the plot look like for \eta<0?
 
  • #27
Sorry, I missed your last post.

The value of \epsilon is 0.056. If you look at the vertical axis, there are also negative values of \eta computed.
 
Back
Top