Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Heavy ODE System

  1. Mar 10, 2005 #1

    Clausius2

    User Avatar
    Science Advisor
    Gold Member

    In order to solve the near field description of a round jet, I want to work out the variables [tex] F(\eta) [/tex], [tex] \rho(\eta) [/tex] and [tex] Y(\eta) [/tex] which represents the self similar stream function, density, and mass fraction respectively. The system obtained is:

    [tex] \Big(\frac{F'}{\rho}\Big)''+\frac{F}{2}\Big(\frac{F'}{\rho}\Big)'=0 [/tex] (1)

    [tex] (\rho Y')'+\frac{F}{2}Y'=0 [/tex] (2)

    besides a function [tex]\rho=\rho(Y)[/tex] which is known previously.

    Boundary conditions are:

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

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

    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.
     
  2. jcsd
  3. Mar 11, 2005 #2

    dextercioby

    User Avatar
    Science Advisor
    Homework Helper

    You can't do it analitically.It's like [tex]\int e^{\sin x} \ dx[/tex]...Can"t find the sollution.I think there are numerical algorithms in which you can find the values you're looking for...

    Daniel.
     
  4. Mar 11, 2005 #3

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

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

    [tex] f^{'''}+ff^{''}=0[/tex] with [itex](0<\eta<+\infty)[/itex]

    [tex]f(0)=f^{'}(0)=0[/tex]

    [tex]f^{'}(+\infty)=1[/tex]

    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.
     
  5. Mar 11, 2005 #4

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

    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 [itex] f^{''}(0)[/itex] via "shooting method". I used simple interpolation to zero into the value of the second derivative at [itex]\eta=0[/itex]. 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.
     

    Attached Files:

  6. Mar 12, 2005 #5

    Clausius2

    User Avatar
    Science Advisor
    Gold Member

    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.
     
  7. Mar 12, 2005 #6

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

    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 [itex]\rho=g(Y(\eta))[/itex] and that this function is known? If that' so, then I can derive the following 5 equations:

    [tex]F^{'}=vg(Y)[/tex]

    [tex]v^{'}=r[/tex]

    [tex]r^{'}+\frac{F}{2}r=0[/tex]

    [tex]Y^{'}=s[/tex]

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

    where [itex]g^{'}(Y)[/itex] is the derivative with respect to Y.
     
  8. Mar 13, 2005 #7

    Clausius2

    User Avatar
    Science Advisor
    Gold Member

    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.
     
  9. Mar 13, 2005 #8

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

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

    Thanks,
    Salty
     
  10. Mar 14, 2005 #9

    Clausius2

    User Avatar
    Science Advisor
    Gold Member

    [tex] y_1=F[/tex]
    [tex] y_2=F'/\rho[/tex]
    [tex] y_3=(F'/\rho)'[/tex]
    [tex] y_4=Y[/tex]
    [tex] y_5=\rho Y[/tex]

    which yields:

    [tex] y_1'=y_2\rho[/tex]
    [tex] y_2'=y_3[/tex]
    [tex] y_3'=-0.5y_1y_3[/tex]
    [tex] y_4'=y_5/\rho[/tex]
    [tex] y_5'=-y_1 y_5/(2\rho)[/tex]
     
  11. Mar 14, 2005 #10

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

    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.
     
  12. Mar 15, 2005 #11

    Clausius2

    User Avatar
    Science Advisor
    Gold Member

    Sorry.

    [tex] \rho=\frac{1}{Y+\epsilon(1-Y)}[/tex]

    with

    [tex]\epsilon<<1[/tex]

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

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

    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.
     
  14. Mar 15, 2005 #13

    Clausius2

    User Avatar
    Science Advisor
    Gold Member

    Till Friday.
     
  15. Mar 15, 2005 #14

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

    Some qualitative progress:

    Consider:

    [tex]\rho=\frac{1}{Y+\epsilon(1-Y)}[/tex]

    And the boundary conditions:

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

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

    We can deduce the values of [itex]\rho[/itex] at these boundaries:

    [tex]\lim_{n\rightarrow+\infty}\rho=\frac{1}{\epsilon}[/tex]

    [tex]\lim_{n\rightarrow-\infty}\rho=1[/tex]

    Thus:

    [tex]\lim_{n\rightarrow+\infty}F^{'}=0[/tex]

    [tex]\lim_{n\rightarrow-\infty}F^{'}=1[/tex]

    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 . . .
     

    Attached Files:

    • pf5.JPG
      pf5.JPG
      File size:
      6.4 KB
      Views:
      60
  16. Mar 15, 2005 #15

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

    Clauius, I don't know where you at, but this is what I'm going with, "win-loose-or draw":

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

    [tex]Y^{'}=s[/tex]

    [tex]r^{'}=-\frac{F}{2}r[/tex]

    [tex]v^{'}=r[/tex]

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

    with:

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

    and:

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

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

    So, I'll keep F(0) and y(0) constant, (with small [itex]\epsilon[/itex]), 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.
     

    Attached Files:

  17. Mar 17, 2005 #16

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

    Just a correction:

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

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

    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 . . .
     
  18. Mar 17, 2005 #17

    Clausius2

    User Avatar
    Science Advisor
    Gold Member

    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! :rofl: .
     
  19. Mar 17, 2005 #18

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

    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
     
  20. Mar 18, 2005 #19

    Clausius2

    User Avatar
    Science Advisor
    Gold Member

    What do you mean? I don't need to know the values of F and Y at [tex]\eta=0[/tex] 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.
     
  21. Mar 19, 2005 #20

    saltydog

    User Avatar
    Science Advisor
    Homework Helper

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

    [tex] y_5=\rho Y^{'}[/tex]

    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.
     
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook

Have something to add?



Similar Discussions: Heavy ODE System
  1. System of ODEs (Replies: 2)

  2. Systems ODE (Replies: 1)

  3. System of ODE (Replies: 0)

  4. System of ODE (Replies: 0)

Loading...