Coupled ODE eigenvalue problem with Chebfun

Hello Everyone.

I am trying to solve an eigenvalue problem with Chebfun but I getting an error. My code is attached below.

dom = [0,1];
opr = chebop(@(y,u,b,s) [-10*4^2*(1+y*0)^(-1.495)*5^2.*s./(50^2*0.001)...
    +4^2*u+(1.495*0^2.*u./(1+y*0).^2)-(1.495^2*0^2*u./(1+y*0).^2)...
    +1j*4^2*2*1.495*0*(1+y*0)^(-1.495)*50^2.*u./(4*2*50*0.001*(1+y*0))...
    +4^2*(1+y*0).^(-2*1.495)*5^4.*u./(50^2*0.001^2)...
    -2*1.495*0*diff(u,1)./(1+y*0)-1j*4^3*(1+y*0)^(-1.495)*5^2.*diff(u,1)/(4*2*50*0.001)...
    +1j*4*2*(1+y*0)^(-1.495)*5^2.*diff(u,1)/(2*50*0.001)...
    +1j*2*2*(1+y*0)^(-1.495)*5^2.*diff(u,1)/(2*50*0.001)-diff(u,2);...
    0.001*b+1j*5*2*u+0.001*diff(b,2)/4^2;-u./(1+y*0)+...
    1e-5*(-4^2.*s+1j*(1+1.495)*0.*s./(1+y*0)+diff(s,2))], dom);
oprlam = chebop(@(y,u,b,s) [(-1j*4^2*(1+y*0)^(-2*1.495)*5^4.*b/(5*2*50^2*0.001^2)...
    -4^2*u/0.001);b;s], dom);
opr.bc = @(y,u,b,s) [u(0); b(0); s(0);u(1);b(1);s(1)];

[V, D] = eigs(opr, oprlam);
D = diag(D);

The error that I am getting is "Dimensions of matrices being concatenated are notconsistent". I am not sure how to get rid of the error. Please help me here. Thank you.

[Hiro Yoshino]

what is chebop like?

[GAURAV]

@Hiro hello. It forms opr(var) = lambda*oprlam (var) type of equation. An example is given here: https://www.chebfun.org/examples/ode-eig/WaveDecay.html

Also, there is this code:

% Assign the equation to two chebops N and B such that N(u) = lambda*B (u).
dom = [0 1];
N = chebop(@(x,u1,u10,u2,u3,u4,u5,u6,u7,u8,u9) [diff(u1)-1/(1-x*0.79)*u2-u7-u9;diff(u2)-u3;diff(u3)-u4;diff(u4)-u2-u3-u6-u7-u9;diff(u5)-u6;diff(u6)-u4-u5;diff(u7)-u8;diff(u8)-u2-u7-u8-u9;diff(u9)-u10;diff(u10)-u2-u7-u9-u10], dom);
B = chebop(@(x,u1,u10,u2,u3,u4,u5,u6,u7,u8,u9) [0.*u1;0.*u2;0.*u3;u1;0.*u5;0.*u6;0.*u7;1.*u7;0.*u9;1.*u9], dom);

% Assign boundary conditions to the chebop.
N.bc = @(x,u1,u10,u2,u3,u4,u5,u6,u7,u8,u9) [u2(0); u3(0); u5(0); u8(0); u10(0); u2(1); u3(1); u5(1); u8(1); u10(1)];

%% Setup preferences for solving the problem.
% Create a CHEBOPPREF object for passing preferences.

options = cheboppref();
options.discretization = 'values';

%% Solve the eigenvalue problem.
[V, D] = eigs(N, B);
D = diag(D);

[Torsten]

You will have to write your first equation as two first-order differential equations and split your equations in real and imaginary parts such that you have 6 boundary conditions for 6 equations. At the moment, you have 6 boundary conditions for 2 equations - and that's a setting the solver doesn't like.

[GAURAV]

@Torsten Hello. Thank you for replying. There are three equations. Is it possible that instead of splitting into real and imaginary I split the equations into 6 first order equation and then try to solve?

[Torsten]

I don't know if CHEBFUN can handle complex-valued ODEs. You will have to look into the documentation or test with a simple example.

[GAURAV]

@Torsten Thank you. I'll try some experiments.