Numerical solution to a 3rd order D.E. - "use computer to solve"

[rudy]

Homework Statement

"Use a computer" to solve the third order differential equation by modifying the boundary conditions (see picture of problem).

Relevant Equations

f*(d^2f/dn^2) + 2*(d^3f/dn^3) = 0

f = 0 @ n = 0
f' = 0 @ n = 0
f' = 1 @ n = inf -or- f" = Constant @ n = 0 (see note)

Hello,

This problem comes from boundary layer theory in fluid mechanics, but we are studying it in heat transfer.

note: Since we are solving this numerically is has been suggested to replace the third boundary condition with f" = constant and then guess a constant. Then we are to check that using that constant f' converges at 1.

They even tell us the value of the constant (0.3321), however I have no clue how to check this answer. If f(n) is an undetermined function how can I solve this numerically?

If there is any information that would help let me know I will do my best to provide. I realize I haven't done much of an attempt at a solution, but I really don't know what to take for a first step! Any tips or suggestions appreciated.

-Rudy

 

[Chestermiller]

Do you have any idea where to start?

 

[rudy]

Hi Chester,

Starting with the first bullet-point suggested in the question I will guess f" = 1 at n = 0.

(They say later in the question that this is incorrect but they suggest it as a starting point so I'll go with it)

But that's only one point, I can't get f' using that. I have two more boundary conditions to work with, so I'm thinking I integrate f". That gives f' = [function of n] + C , where C = arbitrary constant and the [function of n] is the integral of f".

Using the second boundary condition I know that [function of n] + C = 0 at n = 0, but what is the function?!

Keep going, integrate f'. Then I know f = [another function of n] + C_2n + C = 0. Not very helpful either...

I think this is the wrong approach since the book says to solve it numerically. I am trying to solve analytically. I still don't know how to start this using a numerical approach. Hopefully that shows you what I've been doing thus far.

 

[Chestermiller]

Let $$y_1=f$$
$$y_2=f'$$
$$y_3=f''$$
Then:
$$\frac{dy_1}{d\eta}=y_2$$
$$\frac{dy_2}{d\eta}=y_3$$
$$\frac{dy_3}{d\eta}=-\frac{y_1y_3}{2}$$
3 simultaneous 1st order ODEs

 

[rudy]

That was very helpful, thank you! However, I appear to be hitting another snag.

I got an equation for y3 which will not ever approach 1 as n goes to infinity. I will attach a picture of my work, in the end y3 = 6 / (n^3 + C), which approaches zero no matter what C is. Don't I want this to approach 1?

I attached my work, I don't know what I am doing wrong... I also tried plugging y3 back into the equation for y2, but that also approaches zero no matter what the constant is.

Thank you for any help

 

[Chestermiller]

I don't understand what you are doing. You are supposed to be solving this numerically as an initial value problem. You guess an initial value of f'', and then you integrate to large eta (say eta = 10). If you guessed the correct initial value, y2 will approach 1. Do you know how to numerically integrate a set of simultaneous first order ODEs as an initial value problem?

 

[rudy]

Hi Chester-

Went back and reviewed my D.E. textbook on Euler's Method for solving systems. Sorry if it looked like I was looking for a free lunch, I really don't have much experience doing numerical techniques and I didn't recall this method. After reviewing the method I was able to verify that ~0.33 is the correct guess to make y2 converge on 1.

I included a screenshot of my code and output. In the output y2 overshoots 1, but it gets closer as I decrease the step size.

Thanks for your pointers and help setting up the system!

rudy

clear;
clc;
dn = 0.5;

%initial conditions

Guess = 0.330;
y1 = 0;
y2 = 0;
y3 = Guess;

fprintf("        y1        y2        y3        m        n        p\n\n");

for i = 1:15
    m = y2;
    n = y3;
    p = -0.5*y1*y3;
    M = [m n p];
    Y = [y1 y2 y3];
    P = [Y M];
    disp(P);
    y1 = y1 + m*dn;
    y2 = y2 + n*dn;
    y3 = y3 + p*dn;
end

        y1        y2        y3        m        n        p

         0         0    0.3300         0    0.3300         0

         0    0.1650    0.3300    0.1650    0.3300         0

    0.0825    0.3300    0.3300    0.3300    0.3300   -0.0136

    0.2475    0.4950    0.3232    0.4950    0.3232   -0.0400

    0.4950    0.6566    0.3032    0.6566    0.3032   -0.0750

    0.8233    0.8082    0.2657    0.8082    0.2657   -0.1094

    1.2274    0.9410    0.2110    0.9410    0.2110   -0.1295

    1.6979    1.0465    0.1463    1.0465    0.1463   -0.1242

    2.2212    1.1197    0.0842    1.1197    0.0842   -0.0935

    2.7810    1.1617    0.0374    1.1617    0.0374   -0.0520

    3.3619    1.1805    0.0114    1.1805    0.0114   -0.0192

    3.9521    1.1862    0.0018    1.1862    0.0018   -0.0036

    4.5452    1.1871    0.0000    1.1871    0.0000   -0.0000

    5.1387    1.1871   -0.0000    1.1871   -0.0000    0.0000

    5.7323    1.1871    0.0000    1.1871    0.0000   -0.0000

 

[Mark44]

@rudy, both your code and its output are text, so it would be better to post them directly in the input pane, rather than as png images.

For the code, use code tags, like so
[ code=matlab]<your code goes here>[ /code] -- don't include the spaces that I've added.

For the output, use code tags as well, but without the '=matlab' attribute.

 

[rudy]

Thanks Mark - reformatted

 

[Mark44]

Thanks Mark - reformatted

Much better -- thanks!