getting error as Matrix dimensions must agree." can someone help me resolve this

getting error for my code as matrix dimension must agree "have tried to solve coupled wave equation using ode 45 command " error exaclty is as

Error using .*

Matrix dimensions must agree.

Error in odefun (line 58)

dE_omega_dz(1:length(E_omega)/2) =

fft(1i*conj(E1_t).*E2_t.*exp(1i*Dk.*z)

...

Error in odearguments (line 87)

f0 = feval(ode,t0,y0,args{:}); %

ODE15I sets args{1} to yp0.

Error in ode45 (line 113)

[neq, tspan, ntspan, next, t0,

tfinal, tdir, y0, f0, odeArgs,

odeFcn, ...

Error in odefun1 (line 54)

[z, E_omega] =

ode45(@odefun,zspan,E_omega);

code is as- "FUNCTION"

function dE_omega_dz = odefun(z, E_omega,~,~)
dE_omega_dz=zeros(length(E_omega),1);

z=z*10^4;
display(z);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
lambda1=800*10^-9;
c=(3*10^8)
lambda2=400*10^-9;
o=1;
y= 28.78076 %22.4431   %22.39002159      %20:0.25:50  22.39002159     22.443 
ne2=1.5687;
no2=1.6934;
no1=1.6614;
r22=2.1*10^-12                   %electro-optic coefficient in m/v
j=1;
L=1
t=1
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
for Ez=0:500*10^3:85*10^5
noE= no1*(1-0.5*no1^2*r22*Ez)
NEE= ((((sin(y)).^2)/((ne2)^2))+(((cos(y))^2)/((no2)^2)))^-0.5
deltak=-(((4*3.14*(NEE-noE))/ lambda1)) 
Dk(j)=(deltak);
% Dk=deltak;
V=(Ez*4)/10^6
V1(t)=V;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
l=800*10^-9;        % lambda
c=30*10^8;
pi=3.1415926535;
omega=2*pi*c/l;
n_2=6.6508*10^-20;
L_NL=6.8286e-18*1.5; 
LGVM=0.6*10^-3;
I0=0.4*10^11;
% you would have to split the fields back in two:
E1_omega = E_omega(1:end/2);
E2_omega = E_omega(end/2+1:end);

% go back to time space to calculate the nonlinear part:
E1_t = ifft(E1_omega);
E2_t = ifft(E2_omega);
figure(786)
x=-300*10^-15:1*10^-15:300*10^-15;
plot(x,fftshift(E1_t.^2));
% pause(0.2)
N=max(size(E1_t));
to=120e-15/1.655; % initial pulse widthin second
dt=1/120e-15;
dw=1/N/dt*2*pi;
%dw=2*pi*c/l;
 w=(-1*N/2:1:N/2-1)*dw;

% and calculate the derivatives:
dE_omega_dz(1:length(E_omega)/2) = fft(1i*conj(E1_t).*E2_t.*exp(1i*Dk.*z) ...
    + 1i*2*pi*n_2*I0*L_NL/l*(abs(E1_t.^2 + ...
    2*abs(E2_t.^2)).*E1_t));
dE_omega_dz(length(E_omega)/2+1:length(E_omega)) = 1i*w.*L_NL/LGVM * E2_t + fft(1i*E1_t.*E1_t.*exp(-1i*Dk.*z)....
    + 1i*4*pi*n_2*I0*L_NL/l*(2*abs(E1_t.^2 + ...
    abs(E2_t.^2)).*E1_t));


j=j+1;
t=t+1;
%k=k+1;
L=L+1;
o=o+1;
end
end

Calling function is as -

%first electric field envelop in fourier space
clc
clear all
close all
Po=0.4*10^11;
% Dk=-60;
Ao=sqrt(Po);
C=0.5;
pi=3.1415926535;
x=-300*10^-15:1*10^-15:300*10^-15;  
c=3e8;
n=1.655;
l=800*10^-9; 
w1=2*pi*c/l;
z=30*10^-3;
beta=-1.5*10^-29;       
i=sqrt(-1);
t=120*10^-15;
Wo=120*10^-15/1.655; %initial pulse width in second
u=Ao*exp(1i*(beta*z-w1*t)).*exp(-((1+i*(-C))/2)*(x./Wo).^2)
 figure(1)
 
plot(x,abs(u).^2,'b');
title('Input Pulse'); xlabel('Time'); ylabel('Amplitude');
hold on
grid on
E11_omega=fft(fftshift(u));
disp(E11_omega)

%second electric field envelop in fourier space
P1=0.4*10^11;
A1=sqrt(P1);
% C=5;
u1=A1*(exp(1i*(beta*z-2*w1*t)).*exp(-((1+i*(-C))/2)*(x./(Wo)).^2));
% figure(2)
% plot(x,abs(u1).^2,'b');
% title('Input Pulse'); xlabel('Time'); ylabel('Amplitude');
% hold on
% grid on
E22_omega=fft(fftshift(u1));
disp(E22_omega)

%  patching of the E1 and E2 electric fields with:
E_omega=[E11_omega,E22_omega]
disp(E_omega)
% figure(3)
% plot(abs(E_omega),'r');


%a=[1:0.1:3]
zspan=[0 30*10^-7];
% Z1=zspan*10;
E_omega_0=[0 0.0625]
[z, E_omega] = ode45(@odefun,zspan,E_omega);
E11_omega = E_omega(:, 1:end/2);
E22_omega = E_omega(:, end/2+1:end);
E1_t = ifft(E11_omega, [], 2);
E2_t = ifft(E22_omega, [], 2);

figure(101)
subplot(2, 2, 1)
[X, Y] = meshgrid(1:size(E11_omega, 2), z);
mesh(X, Y, abs(fftshift(E11_omega, 2)).^2);
subplot(2, 2, 2)
[X, Y] = meshgrid(1:size(E11_omega, 2), z);
mesh(X, Y,abs(fftshift(E22_omega, 2)).^2);
subplot(2, 2, 3)
[X, Y] = meshgrid(x, z);
mesh(X, Y,abs(fftshift(E1_t, 2)).^2);
subplot(2, 2, 4)
[X, Y] = meshgrid(x, z);
mesh(X, Y,abs(fftshift(E2_t, 2)).^2);

% qq(x,:)=E11_omega(end,:);
figure(1)
 plot(x,fftshift(E1_t(end,:)).^2);
% pause(0.1)
% figure(4)
%  plot(z,E_omega,'linewidth',1)
%     xlabel('z')
%    ylabel('phase (radian)')

Answer 1

[Walter Roberson]

You have

dE_omega_dz(1:length(E_omega)/2) = fft(1i*conj(E1_t).*E2_t.*exp(1i*Dk.*z) ...
    + 1i*2*pi*n_2*I0*L_NL/l*(abs(E1_t.^2 + ...
    2*abs(E2_t.^2)).*E1_t));

This is in a loop in which you are doing

Dk(j)=(deltak);

where you are incrementing j for each iteration of the for Ez loop. Therefore after the first iteration of for Ez, Dk is not a scalar, and the 601 x 1 expression that you get from 1i*conj(E1_t).*E2_t is being .* by a 1 * j variable. In R2016a and earlier that is an error. In R2016b and later, it is permitted and would mean the same thing as if you had used

bsxfun(@times, 1i*conj(E1_t).*E2_t, exp(1i*Dk.*z)

producing a 601 x j result -- a size that is increasing in each iteration of the for Ez loop. But the output size dE_omega_dz(1:length(E_omega)/2) is fixed, leading to a problem of trying to output the wrong amount of information.

Perhaps you should be indexing Dk(j) in the statement.

Answer 2

[Cris LaPierre]

The error is the variables conj(E1_t), E2_t, Dk, and z in the expression

dE_omega_dz(1:length(E_omega)/2) = fft(1i*conj(E1_t).*E2_t.*exp(1i*Dk.*z)

must all be the same size, but at least one of them is not.

A quick debug of your code shows that Dk and z are both scalars (1 x 1) in the offending line, while your other variables are 601 x 1

[Walter Roberson]

Dk is only a scalar on the first iteration of the for Ez loop.

[madhan ravi]

https://www.tutorialspoint.com/matlab/matlab_vector_scalar_multiplication.htm

[Cris LaPierre]

Also, running in 18b gives a slightly different error message

Unable to perform assignment because the left and right sides have a different number of elements.

Error in Untitled>odefun (line 128)
dE_omega_dz(1:length(E_omega)/2) = fft(1i*conj(E1_t).*E2_t.*exp(1i*Dk.*z) ...

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

Error in ode45 (line 115)
  odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);

This is the second timestep (z=0.03). The error suggests there is a dimension mismatch between the left and right sides of the '='. The left side is 601 x 1. The right side is now 601 x 2. Dk is no longer scalar but a 1 x 2. I think you want to do what Walter suggested and index Dk.

[Walter Roberson]

Yup, this is as I explained in my Answer. R2016b or later you get size mismatch because of Dk being nonscalar.

[Cris LaPierre]

@asim asrar, you will get a similar error in line 131 because you use Dk in that calculation as well.

I no longer get this error when I run your code after making those changes.

[asim asrar]

@Cris LaPierre &@ Walter Roberson for your valuable suggestion , ill amend the code accordingly