Loop possible? [Blasius Equation with fsolve]

Hey ! So im pretty new to matlab, but I'm working my way into it. So now I have come across a problem, I'm not sure how to solve: I want to solve the blasius equation for a moving plate, thus different BC's than some might recognise (https://en.wikipedia.org/wiki/Blasius_boundary_layer) f''' + 1/2ff''=0, with BC's: f(0)=0, f'(0)=1 and f'(inf)=0. I know it's very easy to solve with the ODE (+shooting method) or rather the BVP command, but I need this example to apply it to a little more complex model, so I will need to do this with finite-differences(something similar to the keller box method). So since this will be boiling down to a problem formulation of a system of nonlinear equations, I figured why not use fsolve? I figured I would program a function containing all my BC's and equations, so I can just make a grid off "nodes" which will solve the equations at the corresponding node. However I seem to be stuck and or just stupid:

   %%Keller box Basius with fSolve
    eta0=0;
    etaInf=12;
    deltaEta=0.1; %stepsize
    N=(etaInf-eta0)/deltaEta; %number of nodes
    f=zeroes(N,1);g=ones(N,1);h=zeros(N,1); %freeing up space
    anfangswert=[0,1,0.3];
    %%call solver
    Sol = fsolve(blasiuskeller2(f,g,h,N,deltaEta), anfangswert);
    plot(eta0:edltaEta:etaInf, fval2);
  ---------------------------------------------------------------
    %%equations mit discretiszation 
    function fval = blasiuskeller2(f,g,h,N,deltaEta)
    % freeing up space
    fval1=zeros(N-1,1);fval2=zeros(N-1,1);fval3=zeros(N-1,1);
    % Define functions as F(X)=0
        for i=1:N-1
            fval1(i)=(f(i+1)-f(i))/deltaEta - (g(i+1)+g(i))/2;
            fval2(i)=(g(i+1)-g(i))/deltaEta - (h(i+1)+h(i))/2;
            fval3(i)=(h(i+1)-h(i))/deltaEta + (1/8)*(f(i+1)+f(i))*(h(i+1)+h(i));
        end;
    %BC's
    fval(:,1)=[f(1);g(1);fval2(1)]; %initial value
    for i=2:(N-1)
        fval(:,i)=[fval1(i-1);fval3(i-1);fval2(i)];
    end;
    fval(:,N)=[fval1(N-1);1-g(N);fval2(N-1)]; %Boundary value

Pretty sure this is possible, but not sure if I'm on the right track here. Any help would be greatly appreciated!

Answer 1

[Torsten]

I did not check your equations, but already your call to fsolve is incorrect:

anfangswert=horzcat(zeros(1,1:N),ones(1,1:N),0.3*ones(1,1:N));
sol = fsolve(@(x)blasiuskeller2(x(1:N),x(N+1:2*N),x(2*N+1:3*N),N,deltaEta), anfangswert);

Best wishes

Torsten.

[MaxPr]

Thanks that helped already. I do still get en error though: "Size inputs must be scalar."

Could it be that my way of cyling through the equations is off all together?

[Torsten]

Instead of

N=(etaInf-eta0)/deltaEta; %number of nodes

set

N=121;

explicitly.

Note that

number of nodes = number of intervals + 1

Best wishes

Torsten.

[MaxPr]

I found that zeros needs a scalar value. So

zeros(1,1:N)

doesn't work. So I tried something like this:

 [f{1,1:N}]=deal(zeros(1,1));
 [g{1,1:N}]=deal(ones(1,1));
 [h{1,1:N}]=deal(zeros(1,1));
 anfangswert=[f{1,1:N},g{1,1:N},h{1,1:N}];

However obviously matlab doesn't like the assignment to a non-cell array object.

Sry if I'm annoying, you helped me a lot so far!

[Torsten]

anfangswert=horzcat(zeros(1,N),ones(1,N),0.3*ones(1,N));

Best wishes

Torsten.

[MaxPr]

Thanks a lot!