Problem Using ode45 Function

Dear MATLAB users,

I am having this problem while using ode45:

Error using *

Incorrect dimensions for matrix multiplication. Check that the number of columns in the first matrix matches the number of rows in the second matrix.

To operate on each element of the matrix individually, use TIMES (.*) for elementwise multiplication.

Error in EDmulti_rectang_attach_auto_7Eylul>zortit (line 221)

dq=-inv(A)*B*q+inv(A)*F*sin(30*t);

Error in odearguments (line 92)

f0 = ode(t0,y0,args{:}); % ODE15I sets args{1} to yp0.

Error in ode45 (line 107)

odearguments(odeIsFuncHandle,odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);

Error in EDmulti_rectang_attach_auto_7Eylul (line 216)

[t,q]=ode45(@zortit,[0 3], [zeros(1,2*Nx*Ny)])

Related documentation

My original code is like this:

clear all;
global A B F
E=0.7e11;
roplhac=2630;
vu=0.35;
a=0.59;
b=0.4;
psi=b/a;
h=0.002;
ropl=roplhac*h;
roek1hac=7428.57;
roek2hac=7500;
c1=0.040;
d1=0.060;
he1=0.040;
roek1=0.72/(c1*d1);
c2=0.040;
d2=0.060;
he2=0.040;
roek2=0.72/(c2*d2);
rortek=(roek1+roek2)/2;
yogor=rortek/ropl;
yogor1=roek1/ropl;
yogor2=roek2/ropl;
x01=a/2-c1/2;
y01=b/2-d1/2;
ksi1=x01/a;
zeta1=y01/b;
gama1=c1/a;
delta1=d1/b;
x02=a/4-c2/2;
y02=b/4-d2/2;
ksi2=x02/a;
zeta2=y02/b;
gama2=c2/a;
delta2=d2/b;
ksid=[ksi1 ksi2];
zetad=[zeta1 zeta2];
gamad=[gama1 gama2];
deltad=[delta1 delta2];
D=E*h^3/(12*(1-vu^2));
DK=sqrt(D/(ropl*a^4))/(2*pi);
Nx=5; Ny=5;
KAT1=zeros(Nx,Ny,Nx,Ny);
KAT2=zeros(Nx,Ny,Nx,Ny);
durum1=0;
durum2=0;
durum3=0;
durum4=0;
for p=1:2
    for i=1:Nx
        for j=1:Ny
            for r=1:Nx
                for s=1:Ny
                    if p==1
                    if i~=r & j~=s
                        KAT1(i,j,r,s)=KAT1(i,j,r,s)+fKAT1(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                        elseif i==r&j~=s
                            KAT1(i,j,r,s)=KAT1(i,j,r,s)+fKAT2(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                            elseif i~=r&j==s
                                KAT1(i,j,r,s)=KAT1(i,j,r,s)+fKAT3(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                                elseif i==r&j==s
                                    KAT1(i,j,r,s)=KAT1(i,j,r,s)+fKAT4(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                    end
                    elseif p==2
                    if i~=r & j~=s
                        KAT2(i,j,r,s)=KAT2(i,j,r,s)+fKAT1(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                        %durum1=durum1+1;
                        elseif i==r&j~=s
                            KAT2(i,j,r,s)=KAT2(i,j,r,s)+fKAT2(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                            %durum2=durum2+1;
                            elseif i~=r&j==s
                                KAT2(i,j,r,s)=KAT2(i,j,r,s)+fKAT3(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                                %durum3=durum3+1;
                                elseif i==r&j==s
                                    KAT2(i,j,r,s)=KAT2(i,j,r,s)+fKAT4(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                                    %durum4=durum4+1;
                    end
                    end
                m=Ny*(r-1)+s;
                n=Ny*(i-1)+j;
                IEK1(n,m)=KAT1(i,j,r,s);
                IEK2(n,m)=KAT2(i,j,r,s);
               end
            end
        end
    end
end
IEKS=4*(yogor1)*IEK1+4*(yogor2)*IEK2;
omegam=zeros(Nx*Ny,1);
for rr=1:Nx
    for ss=1:Ny
        m=Ny*(rr-1)+ss;
          omegam(m)=pi^2*(rr^2+(ss/psi)^2); 
    end
end
sort(omegam);
OMEGA2=zeros(Nx*Ny,Nx*Ny);
for zz=1:Nx*Ny
    OMEGA2(zz,zz)=omegam(zz)^2;
end
[VV,DD]=eig(OMEGA2,(eye(Nx*Ny,Nx*Ny)+IEKS));
dddd=diag(DD);
eeee=sqrt(dddd);
ffff=sort(eeee);
frekanslar_Hz=DK*(ffff);
nx=30;ny=20;
deltax=a/nx;
deltay=b/ny;
MODSHAPE=zeros((nx+1),(ny+1));
k=1;
%for k=1:Nx*Ny
    for jj=1:(ny+1)
        y(jj)=(jj-1)*deltay;
        for kk=1:(nx+1)
            x(kk)=(kk-1)*deltax;
            for i=1:Nx
                for j=1:Ny
                    n=Ny*(i-1)+j;
                    MODSHAPE(kk,jj)=MODSHAPE(kk,jj)+sin(i*pi*x(kk)/a)*sin(j*pi*y(jj)/b)*VV(n,1);
                end
            end
        end
    end
AA=eye(Nx*Ny,Nx*Ny)+IEKS;
A=[eye(Nx*Ny,Nx*Ny) zeros(Nx*Ny,Nx*Ny);
  zeros(Nx*Ny,Nx*Ny) AA];
B=[zeros(Nx*Ny,Nx*Ny) eye(Nx*Ny,Nx*Ny);
   OMEGA2             zeros(Nx*Ny,Nx*Ny)];
ff=zeros(Nx*Ny,1);
F0=100;
for r=1:Nx
    for s=1:Ny
        m=Ny*(r-1)+s;
        ff(m)=(4*F0/(a*b))*sin(r*pi*(1/4))*sin(s*pi*(1/4));
    end
end
F=[zeros(1,Nx*Ny) ff']';
[t,q]=ode45(@zortit,[0 3], [zeros(1,2*Nx*Ny)])

Error using *
Incorrect dimensions for matrix multiplication. Check that the number of columns in the first matrix matches the number of rows in the second matrix. To operate on each element of the matrix
individually, use TIMES (.*) for elementwise multiplication.

Error in solution>zortit (line 143)
dq=-inv(A)*B*q+inv(A)*F*sin(30*t);

Error in odearguments (line 92)
f0 = ode(t0,y0,args{:}); % ODE15I sets args{1} to yp0.

Error in ode45 (line 107)
odearguments(odeIsFuncHandle,odeTreatAsMFile, solver_name, ode, tspan, y0, options, varargin);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [dq]=zortit(q,t)
global A B F
            dq=-inv(A)*B*q+inv(A)*F*sin(30*t);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [K1] = fKAT1(i,j,r,s,ksi,zeta,gama,delta)
K1= (1/((i^2-r^2)*(j^2-s^2)*pi^2))*( r*sin(i*pi*(ksi+gama))*cos(r*pi*(ksi+gama))-...
i*cos(i*pi*(ksi+gama))*sin(r*pi*(ksi+gama))-...
 r*sin(i*pi*ksi)*cos(r*pi*ksi) + i*cos(i*pi*ksi)*sin(r*pi*ksi))*...
( s*sin(j*pi*(zeta+delta))*cos(s*pi*(zeta+delta)) - ...
j*cos(j*pi*(zeta+delta))*sin(s*pi*(zeta+delta)) - ...
 s*sin(j*pi*zeta)*cos(s*pi*zeta) + ...
 j*cos(j*pi*zeta)*sin(s*pi*zeta) );
end
%%%
function [K2] = fKAT2(i,j,r,s,ksi,zeta,gama,delta)
K2=(1/((j^2-s^2)*pi))*...
(gama/2-(1/(4*r*pi))*(sin(2*r*pi*(ksi+gama))-sin(2*r*pi*ksi)))*...
( s*sin(j*pi*(zeta+delta))*cos(s*pi*(zeta+delta)) - ...
j*cos(j*pi*(zeta+delta))*sin(s*pi*(zeta+delta)) - ...
s*sin(j*pi*zeta)*cos(s*pi*zeta) + ...
j*cos(j*pi*zeta)*sin(s*pi*zeta) );
end
%%%
function [K3] = fKAT3(i,j,r,s,ksi,zeta,gama,delta)
K3=(1/((i^2-r^2)*pi))*...
(delta/2-(1/(4*s*pi))*(sin(2*s*pi*(zeta+delta))-sin(2*s*pi*zeta)))*...
( r*sin(i*pi*(ksi+gama))*cos(r*pi*(ksi+gama))-...
i*cos(i*pi*(ksi+gama))*sin(r*pi*(ksi+gama))-...
r*sin(i*pi*ksi)*cos(r*pi*ksi) + i*cos(i*pi*ksi)*sin(r*pi*ksi)  );
end
%%%
function [K4] = fKAT4(i,j,r,s,ksi,zeta,gama,delta)
K4=(gama/2-(1/(4*r*pi))*(sin(2*r*pi*(ksi+gama))-sin(2*r*pi*ksi)))*...
(delta/2-(1/(4*s*pi))*(sin(2*s*pi*(zeta+delta))-sin(2*s*pi*zeta)));
end

Do you have suggestion? Where is my fault about ode45?

Thank you very much in advance.

[Walter Roberson]

function [dq]=zortit(q,t)
global A B F
            dq=-inv(A)*B*q+inv(A)*F*sin(30*t);
end

Relative to the function, A, B, and F are all constants. So you do not need to keep recalculating inv(A)*B and inv(A)*F : you can calculate those ahead of time, something like

global AB AF
AB = inv(A)*B;
AF = inv(A)*F;

....
function [dq]=zortit(q,t)
global AB AF
            dq = -AB*q + AF*sin(30*t);
end

But then:

Using inv() is seldom desirable for efficiency and numeric stability reasons. You should instead be using something like

AB = A\B;
AF = A\F;

Then there is the fact that using global is not recommended; see paramfun and write code something like

[t, q] = ode45(@(t,q)zortit(t,q,AB,AF),[0 3], [zeros(1,2*Nx*Ny)])
function [dq]=zortit(t, q, AB, AF)
  dq = -AB*q + AF*sin(30*t);
end

These are improvements to the code, not something that would "repair" the code.

Answer 1

[Cris LaPierre]

The error I see is that you have not defined the inputs to your ode function in the correct order. Namely, t must be the first input.

function [dq]=zortit(t,q)

At least that gets rid of the error:

global A B F
E=0.7e11;
roplhac=2630;
vu=0.35;
a=0.59;
b=0.4;
psi=b/a;
h=0.002;
ropl=roplhac*h;
roek1hac=7428.57;
roek2hac=7500;
c1=0.040;
d1=0.060;
he1=0.040;
roek1=0.72/(c1*d1);
c2=0.040;
d2=0.060;
he2=0.040;
roek2=0.72/(c2*d2);
rortek=(roek1+roek2)/2;
yogor=rortek/ropl;
yogor1=roek1/ropl;
yogor2=roek2/ropl;
x01=a/2-c1/2;
y01=b/2-d1/2;
ksi1=x01/a;
zeta1=y01/b;
gama1=c1/a;
delta1=d1/b;
x02=a/4-c2/2;
y02=b/4-d2/2;
ksi2=x02/a;
zeta2=y02/b;
gama2=c2/a;
delta2=d2/b;
ksid=[ksi1 ksi2];
zetad=[zeta1 zeta2];
gamad=[gama1 gama2];
deltad=[delta1 delta2];
D=E*h^3/(12*(1-vu^2));
DK=sqrt(D/(ropl*a^4))/(2*pi);
Nx=5; Ny=5;
KAT1=zeros(Nx,Ny,Nx,Ny);
KAT2=zeros(Nx,Ny,Nx,Ny);
durum1=0;
durum2=0;
durum3=0;
durum4=0;
for p=1:2
    for i=1:Nx
        for j=1:Ny
            for r=1:Nx
                for s=1:Ny
                    if p==1
                        if i~=r & j~=s
                            KAT1(i,j,r,s)=KAT1(i,j,r,s)+fKAT1(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                        elseif i==r&j~=s
                            KAT1(i,j,r,s)=KAT1(i,j,r,s)+fKAT2(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                        elseif i~=r&j==s
                            KAT1(i,j,r,s)=KAT1(i,j,r,s)+fKAT3(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                        elseif i==r&j==s
                            KAT1(i,j,r,s)=KAT1(i,j,r,s)+fKAT4(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                        end
                    elseif p==2
                        if i~=r & j~=s
                            KAT2(i,j,r,s)=KAT2(i,j,r,s)+fKAT1(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                            %durum1=durum1+1;
                        elseif i==r&j~=s
                            KAT2(i,j,r,s)=KAT2(i,j,r,s)+fKAT2(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                            %durum2=durum2+1;
                        elseif i~=r&j==s
                            KAT2(i,j,r,s)=KAT2(i,j,r,s)+fKAT3(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                            %durum3=durum3+1;
                        elseif i==r&j==s
                            KAT2(i,j,r,s)=KAT2(i,j,r,s)+fKAT4(i,j,r,s,ksid(p),zetad(p),gamad(p),deltad(p));
                            %durum4=durum4+1;
                        end
                    end
                    m=Ny*(r-1)+s;
                    n=Ny*(i-1)+j;
                    IEK1(n,m)=KAT1(i,j,r,s);
                    IEK2(n,m)=KAT2(i,j,r,s);
                end
            end
        end
    end
end
IEKS=4*(yogor1)*IEK1+4*(yogor2)*IEK2;
omegam=zeros(Nx*Ny,1);
for rr=1:Nx
    for ss=1:Ny
        m=Ny*(rr-1)+ss;
        omegam(m)=pi^2*(rr^2+(ss/psi)^2);
    end
end
sort(omegam);
OMEGA2=zeros(Nx*Ny,Nx*Ny);
for zz=1:Nx*Ny
    OMEGA2(zz,zz)=omegam(zz)^2;
end
[VV,DD]=eig(OMEGA2,(eye(Nx*Ny,Nx*Ny)+IEKS));
dddd=diag(DD);
eeee=sqrt(dddd);
ffff=sort(eeee);
frekanslar_Hz=DK*(ffff);
nx=30;ny=20;
deltax=a/nx;
deltay=b/ny;
MODSHAPE=zeros((nx+1),(ny+1));
k=1;
%for k=1:Nx*Ny
for jj=1:(ny+1)
    y(jj)=(jj-1)*deltay;
    for kk=1:(nx+1)
        x(kk)=(kk-1)*deltax;
        for i=1:Nx
            for j=1:Ny
                n=Ny*(i-1)+j;
                MODSHAPE(kk,jj)=MODSHAPE(kk,jj)+sin(i*pi*x(kk)/a)*sin(j*pi*y(jj)/b)*VV(n,1);
            end
        end
    end
end
AA=eye(Nx*Ny,Nx*Ny)+IEKS;
A=[eye(Nx*Ny,Nx*Ny) zeros(Nx*Ny,Nx*Ny);
    zeros(Nx*Ny,Nx*Ny) AA];
B=[zeros(Nx*Ny,Nx*Ny) eye(Nx*Ny,Nx*Ny);
    OMEGA2             zeros(Nx*Ny,Nx*Ny)];
ff=zeros(Nx*Ny,1);
F0=100;
for r=1:Nx
    for s=1:Ny
        m=Ny*(r-1)+s;
        ff(m)=(4*F0/(a*b))*sin(r*pi*(1/4))*sin(s*pi*(1/4));
    end
end
F=[zeros(1,Nx*Ny) ff']';
[t,q]=ode45(@zortit,[0 3], [zeros(1,2*Nx*Ny)])

t = 2497×1

0 0.0002 0.0005 0.0007 0.0009 0.0012 0.0014 0.0016 0.0019 0.0022

q = 2497×50

0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 -0.0000 0.0000 0.0001 0.0000 -0.0001 0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 -0.0001 0.0002 0.0004 0.0001 -0.0005 0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 -0.0001 0.0003 0.0008 0.0003 -0.0012 0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 -0.0003 0.0006 0.0015 0.0005 -0.0021 0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 -0.0004 0.0009 0.0023 0.0008 -0.0034 0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 -0.0007 0.0013 0.0033 0.0012 -0.0049 0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 -0.0010 0.0017 0.0045 0.0016 -0.0067 0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 -0.0013 0.0022 0.0059 0.0022 -0.0090 0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 -0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 -0.0000 -0.0000 -0.0000 0.0000 0.0000 0.0000 0.0000 0.0000 -0.0000 -0.0020 0.0029 0.0084 0.0033 -0.0128

plot(t,q)

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function [dq]=zortit(t,q)
global A B F
dq=-inv(A)*B*q+inv(A)*F*sin(30*t);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [K1] = fKAT1(i,j,r,s,ksi,zeta,gama,delta)
K1= (1/((i^2-r^2)*(j^2-s^2)*pi^2))*( r*sin(i*pi*(ksi+gama))*cos(r*pi*(ksi+gama))-...
    i*cos(i*pi*(ksi+gama))*sin(r*pi*(ksi+gama))-...
    r*sin(i*pi*ksi)*cos(r*pi*ksi) + i*cos(i*pi*ksi)*sin(r*pi*ksi))*...
    ( s*sin(j*pi*(zeta+delta))*cos(s*pi*(zeta+delta)) - ...
    j*cos(j*pi*(zeta+delta))*sin(s*pi*(zeta+delta)) - ...
    s*sin(j*pi*zeta)*cos(s*pi*zeta) + ...
    j*cos(j*pi*zeta)*sin(s*pi*zeta) );
end
%%%
function [K2] = fKAT2(i,j,r,s,ksi,zeta,gama,delta)
K2=(1/((j^2-s^2)*pi))*...
    (gama/2-(1/(4*r*pi))*(sin(2*r*pi*(ksi+gama))-sin(2*r*pi*ksi)))*...
    ( s*sin(j*pi*(zeta+delta))*cos(s*pi*(zeta+delta)) - ...
    j*cos(j*pi*(zeta+delta))*sin(s*pi*(zeta+delta)) - ...
    s*sin(j*pi*zeta)*cos(s*pi*zeta) + ...
    j*cos(j*pi*zeta)*sin(s*pi*zeta) );
end
%%%
function [K3] = fKAT3(i,j,r,s,ksi,zeta,gama,delta)
K3=(1/((i^2-r^2)*pi))*...
    (delta/2-(1/(4*s*pi))*(sin(2*s*pi*(zeta+delta))-sin(2*s*pi*zeta)))*...
    ( r*sin(i*pi*(ksi+gama))*cos(r*pi*(ksi+gama))-...
    i*cos(i*pi*(ksi+gama))*sin(r*pi*(ksi+gama))-...
    r*sin(i*pi*ksi)*cos(r*pi*ksi) + i*cos(i*pi*ksi)*sin(r*pi*ksi)  );
end
%%%
function [K4] = fKAT4(i,j,r,s,ksi,zeta,gama,delta)
K4=(gama/2-(1/(4*r*pi))*(sin(2*r*pi*(ksi+gama))-sin(2*r*pi*ksi)))*...
    (delta/2-(1/(4*s*pi))*(sin(2*s*pi*(zeta+delta))-sin(2*s*pi*zeta)));
end