The code speed is too slow

I need plot(z, W_res) but my code is incredibly slow and I can't do it.

How can I make the code faster so that it can calculate everything?

This is the expression I want to calculate:

My code:

D = 100;
e_f = 1.88e-18;
m = 9.11e-31;
k_b = 1.38e-23;
h_bar = 1.05e-34;
ksi = 1e-9;
t = 1;
A = 2*m/h_bar^2;
k_f = sqrt(e_f / (pi * k_b * 1.2));
N_0 = 10e28;
v_f = 1.57e6;
l = 1.5e-10;
tau = l / v_f;
gamma = h_bar / (2*pi * k_b * 1.2 * tau);
E = e_f / (pi * k_b * 1.2);
norm = 2*m / h_bar^2;
k0 = (ksi/h_bar)*sqrt(2.*m.*pi.*k_b.*1.2);
coef = 1;
w_n = @(n) t*(2*n+1);
n_max = 49;
cuA = m / (pi*h_bar*N_0*tau);
coef_K = cuA / (2*pi*ksi);
coef_w = (cuA * ksi * A) / (4*pi);
wt_sum = wt(1.2, 49);
k_pl = @(e) k0*sqrt( e + wt_sum); % k_+ off dim
s_w_func = @(e,z)( sin(k_pl(e).*z) .* sin(k_pl(e).*(z-D)) ) ./ (k_pl(e) .* sin(k_pl(e)*D));
S_w = @(z) imag(integral(@(e)s_w_func(e,z), 0, E,'ArrayValued',1 )); % off dim
K = @(e, z) k0*sqrt(e + 1i*wt_sum - coef_K * S_w(z));
inner_fun = @(e,z) 1/K(e,z) .* ( cos(integral(@(z_)K(e,z_), z-D, z)) ./ sin(integral(@(z_)K(e,z_), 0, D)) - cot(integral(@(z_)K(e,z_), 0, D)));
W_res = @(z)wt_sum - 1i*coef_w .* integral(@(e)inner_fun(e,z), 0, E,'ArrayValued',1); % off-dim, norm on pi*k_b*Tc
% Определение вектора z от 5 до 100 с шагом 0.1
z_values = 5:0.1:100;
% Инициализация вектора для хранения значений W_res
W_res_values = zeros(size(z_values));
% Вычисление значений W_res для каждого значения z
for i = 1:length(z_values)
    W_res_values(i) = W_res(z_values(i));
end
% Построение графика
figure;
plot(z_values, real(W_res_values), 'LineWidth', 2);
xlabel('z');
ylabel('Re(W\_res)');
title('График W\_res от z');
grid on;
function result = wt(t, n)
result = 0;
for i = 1:n
    result = result + t * (2 * i + 1);
end
end

[Alexander]

Could you format your code as code? Mark it and press on code button on the insert of the part of the upper line.

[Torsten]

Did you get senseful values for W_res_values with your code and only the speed is too slow ? Or are the results unsatisfactory, too ?

[Dmitry]

Torsten, I couldn't get any values because it took too long to count.

[Dmitry]

Alexander, sorry, What do you mean?

[Alexander]

Every thing okay now. My former view of you code was text. Now it is code.

[Torsten]

If the code does not produce a result, you should work on this first (e.g. by checking your constants and equations). What is striking are your extreme parameter values used in the setup.

[Dmitry]

Alexander, My code runs, but I can't get the result as it takes too long to count.

[Dmitry]

Torsten, All constants and equations are correct. Perhaps there is some other more optimal way to solve this problem?

[Torsten]

Even when I try to compute

sin(integral(@(z_)K(e,z_), 0, D))

for different values of e in between 0 and E in order to see whether it becomes 0 somewhere, integral does not finish.

Could you test it on your computer and report the result of

D = 100;
e_f = 1.88e-18;
m = 9.11e-31;
k_b = 1.38e-23;
h_bar = 1.05e-34;
ksi = 1e-9;
t = 1;
A = 2*m/h_bar^2;
k_f = sqrt(e_f / (pi * k_b * 1.2));
N_0 = 10e28;
v_f = 1.57e6;
l = 1.5e-10;
tau = l / v_f;
gamma = h_bar / (2*pi * k_b * 1.2 * tau);
E = e_f / (pi * k_b * 1.2);
norm = 2*m / h_bar^2;
k0 = (ksi/h_bar)*sqrt(2.*m.*pi.*k_b.*1.2);
coef = 1;
w_n = @(n) t*(2*n+1);
n_max = 49;
cuA = m / (pi*h_bar*N_0*tau);
coef_K = cuA / (2*pi*ksi);
coef_w = (cuA * ksi * A) / (4*pi);
wt_sum = wt(1.2, 49);
k_pl = @(e) k0*sqrt( e + wt_sum); % k_+ off dim
s_w_func = @(e,z)( sin(k_pl(e).*z) .* sin(k_pl(e).*(z-D)) ) ./ (k_pl(e) .* sin(k_pl(e)*D));
S_w = @(z) imag(integral(@(e)s_w_func(e,z), 0, E,'ArrayValued',1 )); % off dim
K = @(e, z) k0*sqrt(e + 1i*wt_sum - coef_K * S_w(z));
e = linspace(0,E,20)
result = arrayfun(@(e)sin(integral(@(z_)K(e,z_), 0, D)),e)

function result = wt(t, n)
result = 0;
for i = 1:n
    result = result + t * (2 * i + 1);
end
end

[Dmitry]

Torsten,

This code calculated:

[Torsten]

Do you see the problem ? There seems to be a singularity at e=0.