problem with vpa and symsum
92 views (last 30 days)
Show older comments
george veropoulos
on 12 Nov 2024 at 14:58
Hi i make the fumtiom
function z= Escattheory(r,phi)
%UNTITLED Summary of this function goes here
% Detailed explanation goes here
syms kk
[f,N,Nc,a,ra,k0,Z0] = parameter();
E0=1;
Z1=k0.*r;
Z2=k0.*ra;
factor1=besselj(kk,Z2).*besselh(kk,2,Z1)./besselh(kk,2,Z2)
factor2=(-1j).^kk.*cos(kk.*phi)
factor=-E0.*factor1.*factor2.*en(kk);
%z=vpa(symsum(factor,k,0,inf),3)
f_all=symsum(factor,kk,0,300);
y=vpa(f_all,3)
end
when i call the function Escattheory(r,phi) with argument r=ra=1 and phi =0 the resultas is - 1.22 + 5.05e-5i
if i use the sum like this
function z= Escattheory_test(r,phi)
%UNTITLED Summary of this function goes here
% Detailed explanation goes here
Escattheory_test(r,phi)
function y=Escattheory_test(r,phi)
%UNTITLED Summary of this function goes here
% Detailed explanation goes here
[f,N,Nc,a,ra,k0,Z0] = parameter();
E0=1;
Z1=k0.*r;
Z2=k0.*ra;
kk=0:200
factor1=besselj(kk,Z2).*besselh(kk,2,Z1)./besselh(kk,2,Z2)
factor2=(-1j).^kk.*cos(kk.*phi)
factor=-E0.*factor1.*factor2.*e_n(kk);
%z=vpa(symsum(factor,k,0,inf),3)
y=sum(factor)
end
the reuslts is -1-j0 what happen ??
thank you
Accepted Answer
Torsten
on 13 Nov 2024 at 11:07
Moved: Torsten
on 13 Nov 2024 at 11:09
The numerical bessel functions are not able to compute your sum correctly (for higher values of k, they return NaN) (see below). That's the reason why you get different results.
r = 1;
phi = 0;
Escattheory(r,phi)
Escattheory_test(r,phi)
function z= Escattheory(r,phi)
%UNTITLED Summary of this function goes here
% Detailed explanation goes here
syms kk
[f,N,Nc,a,ra,k0,Z0] = parameter();
E0=1;
Z1=k0.*r;
Z2=k0.*ra;
factor1=besselj(kk,Z2).*besselh(kk,2,Z1)./besselh(kk,2,Z2);
factor2=(-1j).^kk.*cos(kk.*phi);
factor=-E0.*factor1.*factor2.*e_n(kk)
%z=vpa(symsum(factor,k,0,inf),3)
f_all=symsum(factor,kk,0,300);
z=vpa(f_all,3);
end
function y=Escattheory_test(r,phi)
%UNTITLED Summary of this function goes here
% Detailed explanation goes here
[f,N,Nc,a,ra,k0,Z0] = parameter();
E0=1;
Z1=k0.*r;
Z2=k0.*ra;
kk=0:300;
factor1=besselj(kk,Z2).*besselh(kk,2,Z1)./besselh(kk,2,Z2);
factor2=(-1j).^kk.*cos(kk.*phi);
factor=-E0.*factor1.*factor2.*en(kk)
%z=vpa(symsum(factor,k,0,inf),3)
y=sum(factor);
end
function [f,N,Nc,a,ra,k0,Z0] = parameter()
%UNTITLED Summary of this function goes here
c0=3e8;
Z0=120.*pi;
ra=1;
a=0.92;
N=30;
Nc=30;
f=300e6;
lambda=c0./f;
k0=2*pi./lambda;
end
function y=en(k)
if k==0
y=1;
else
y=2;
end
end
function y=e_n(k)
y = piecewise(k==0,1,k~=0,2);
end
13 Comments
Torsten
on 14 Nov 2024 at 12:37
Edited: Torsten
on 14 Nov 2024 at 13:11
As said, I don't know precisely why "piecewise" works with the symbolic solution while your if-construct gives a wrong answer.
My guess is that your function "e_n" is called with a symbolic (unspecified) variable (kk) that is not identified as "0" (which would give 1 for y), but as "else" (which always gives 2 for y).
But you might want to contact MATLAB support for a definite answer:
More Answers (0)
See Also
Products
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!