, I need help! I do not get to have the graph of the solution of a differential equation.

• I have a differential equation of the first order. It describes the basic velocity profile of fluid flowing over an inclined plane.
• I used the Euler method with Matlab to solve the equation.
• My goal is to plot the graph of the solution of this equation. I want to have the velocity U on axis (ox) on the interval [0,18], in fonction of the height of the flow channel y, on axis (oy) in the interval [0,200] .
• No errors in the program are displayed, but the result I get is different from that of the article on which I work. The velocity profile should be parabolic while I get a straight line!
• I Attach the program that I do, the graph that I get and the graph that I must have.Could you please help me understand and solve the problem? I will be very grateful.

[Nadjah]

Sorry, I think I have not attached the graphs! Thank you.

Answer 1

[Roger Stafford]

I am guessing that you left out the parentheses in the denominator of F, which would presumably be this instead:

   F=@(x,s)(2-(x./100))./(1+sqrt(1+0.003*x.^2.*(1-exp(-x./25)).^2.*(200-x)./50));
                         ^                                                     ^

The values for U I get with this alteration do not agree with those in your 'profile vitesse.pdf' plot but the shape appears to be similar to that plot.

Also I find it curious that F, as you have defined it, does not depend on the second argument 's', so you are simply finding U(x) = the integral of F(t) from t = 0 to x where F does not depend on U. In other words, it is not really a differential equation but merely an integral which could have been solved using 'cumtrapz'.

[Nadjah]

Hello, thank you very much for your prompt response. You were right for the brackets! Finally it solved my problem. I found the same shape, however, the U values are different from those of the article, I don't understand why? interestingly, Using Maple, I find the same shape with the same values of U! I started to deal the problem using the cumtrap function by considering that U (x) is an integral of F (t), but I found a straight line too, finally, because the parentheses forgetten! thank you again.

[Roger Stafford]

Perhaps the discrepancy in U values is somehow related to the absence of U in the function F, which would make it a true differential equation. There must be some reason you wrote:

   F=@(x,s)(2-(x./100).....

as though there should be some such dependence.

[Nadjah]

Hello, I'm sure of the expression of F. The differential equation was: U '(x) + l ^ 2 * (U' (x)) ^ 2 = (L-x) / L, such that l = 0.41*x* [1-exp (-x / 25) ] and L = 200. I wrote again the equation as a full square for linearization, which gave the final expression of the equation U '(x) = F (x) ie U (x) = integral (F (x)).

[Roger Stafford]

I don't have time at the moment to check your formula in careful detail, but I have noticed one discrepancy, namely, that by my computations, the factor for which you give the value 0.003 ought to be (0.41)^2 = 0.1681 .

[Roger Stafford]

When I replaced 0.003 by 0.1681 in your function, the resulting curve looks much more like that shown in "profile vitesse.pdf", which undoubtedly means that the 0.003 is a mistake. Using 0.1681 leads to U = 16.2473 for x = 200, (but you will not get that much accuracy using an 'h' as large as 0.2.) The pdf figure shows U somewhere near 15.3 for this value of x. I don't know what accounts for this difference. Possible the integration method used in the pdf was somewhat inaccurate, or possibly the given 0.41 value in the definition of 'l' is not as accurate as that used in the pdf method.

[Nadjah]

Thank you very much, I found the error. Finally, F (x) = (2-(x./100))./(1+sqrt(1+0.003*x.^2.*(1-exp(-x./25)).^2.*(200-x))) ^ ie 0,003 = (4 * (0.41) ^ 2) / 200. Thus, the result is identical to that of the article. Your remark pushed me to revise the calculations. I didn't make this mistake when I used Maple, that's why I got the same results that article. I am new to matlab. This problem is my first application. I am trying to learn and apply in parallel on the problems of the article on which I work. Thank you again.

[Roger Stafford]

Not quite identical. There is the discrepancy between U being 16.2473 with matlab and somewhere near 15.3 for the article for x = 200.

[Nadjah]

Hello; The graph '' profile vitesse.pdf '' is the one I got with Maple and that I found (before to solve the problem with Matlab)) close to that of the article (herewith the graph of the article). Whilst if we correct the term in denominator: F=@(x,s)(2-(x./100))./(1+sqrt(1+0.003*x.^2.*(1-exp(-x./25)).^2.* (200-x)./50) ) which becomes: F=@(x,s)(2-(x./100))./(1+sqrt(1+0.003*x.^2.*(1-exp(-x./25)).^2.*(200-x))) you will notice that the result of Matlab is much better! Thank you.

Answer 2

[Andres Bayona]

Hi, I was wondering if it were possible that you shared the code with me. I would really appreciate it!