how can I calculate my expression with double integration?

5 visualizaciones (últimos 30 días)
Dmitry
Dmitry el 25 de Feb. de 2024
Comentada: Dmitry el 25 de Feb. de 2024
how can I calculate my expression with double integration?
at startup, I get an error:
Operator '.*' is not supported for operands of type 'function_handle'.
Error in sw_plt2 (line 32)
s_w_func = ( sin(k_pl.*z) .* sin(k_pl.*(z-D)) ) ./ (k_pl .* sin(k_pl*D));
My code:
z = 5:0.1:100;
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);
w_z = zeros(size(z));
coef_K = cuA / (2*pi*ksi);
coef_w = (cuA * ksi * A) / (4*pi);
%syms e;
wt_sum = wt(1.2, 49);
k_pl = @(e) k0*sqrt( e + wt_sum); % k_+ off dim
s_w_func = ( sin(k_pl.*z) .* sin(k_pl.*(z-D)) ) ./ (k_pl .* sin(k_pl*D));
Operator '.*' is not supported for operands of type 'function_handle'.
S_w = imag( @(e) integral(s_w_func), 0, E ); % off dim
K = @(e, z) k0*sqrt(e + 1i*wt_sum - coef_K * S_w);
inner_fun = 1/K .* ( cos(@(z_)integral(K, z-D, z)) ./ sin(@(z_)integral(K, 0, D)) - cot(@(z_)integral(K, 0, D)));
W_res = wt_sum - 1i*coef_w .* (@(e) integral(inner_fun, 0, E)); % off-dim, norm on pi*k_b*Tc
figure;
plot(z, W_res);
xlabel('z');
ylabel('w_res(z)');
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
  2 comentarios
Dyuman Joshi
Dyuman Joshi el 25 de Feb. de 2024
s_w_func = ( sin(k_pl.*z) .* sin(k_pl.*(z-D)) ) ./ (k_pl .* sin(k_pl*D));
What is this line of code supposed to do?
Please share the expression in mathematical form so we can provide appropriate suggestions.
Dmitry
Dmitry el 25 de Feb. de 2024
I need plot(z, W_res)
As for "s_w_func", I want to integrate this by "e".
This is the expression I want to calculate:

Iniciar sesión para comentar.

Respuestas (2)

Torsten
Torsten el 25 de Feb. de 2024
Maybe you mean this, but I'm not sure:
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
function result = wt(t, n)
result = 0;
for i = 1:n
result = result + t * (2 * i + 1);
end
end
  5 comentarios
Torsten
Torsten el 25 de Feb. de 2024
Editada: Torsten el 25 de Feb. de 2024
First try for a single value:
W_res(1)
If this works, try for your array:
z = 5:0.1:100;
W_resnum = arrayfun(@(z) W_res(z),z)
But my guess is that the denominators in your functions will have zeros in the region of integration.
Dmitry
Dmitry el 25 de Feb. de 2024
How can I make these calculations faster?

Iniciar sesión para comentar.


John D'Errico
John D'Errico el 25 de Feb. de 2024
Editada: John D'Errico el 25 de Feb. de 2024
The important takeaway to use is, you CANNOT perform arithmetic operations between a pair of function handles. Maybe you think you should be able to, but you CANNOT. Yours is a common error we see by new users. Language syntax is irrefutable and inescapable.
HOWEVER, you can perform those same operations on the RESULT of a function handle. For example:
f1 = @(x) x + 1;
f2 = @(x) x.^2;
f3 = @(x) f1(x).*f2(x);
f3(5)
ans = 150
By way of comparison, see what happens when you try to multiply those function handles (I must do this at the end, as MATLAB in Answers will not allow me to do anything else after an error occurs.)
f1*f2
Operator '*' is not supported for operands of type 'function_handle'.
Do you see this is the same error you got? The solution is as I showed above. You need to operate on the RESULT of the function handles, not on the function handles themselves.

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by