
clear all
syms s1 s2 zz lamda alpha c a w;%variance and parameter
C=10^-14;
Z=50:10:300;
k=2*pi/lamda;
rho=(0.545*c^2*k^2*zz)^(-3/5);
f1=k^2/(zz^2)*s1*s2;
f2=exp(j*k/(2*zz)*(s1^2-s2^2))*exp(-1/rho^2*(s1^2+s2^2))*besseli(0,2*s1*s2/rho^2);
%上面f1,f2为被积函数中与源无关部分,下面为sgb
sgb=besselj(0,alpha*s1)*besselj(0,alpha*s2)*exp(-s1^8/w^8)*exp(-s2^8/w^8);
f=f1*f2*sgb;%integrand
f=subs(f,{zz,lamda,alpha,c,w},{Z,1550*10^-9,600,C,0.05});%这是一个矢量;
disp(f);
intValue=ones(size(Z));
for i=1:length(Z)
disp(i);
ff=inline(f(i));
intValue(i)=quad2dg(ff,0,0.1,0,0.1,10^-20);
end;
plot(Z',abs(intValue))
syms s1 s2 zz lamda alpha c a w;%variance and parameter
C=10^-14;
Z=50:10:300;
k=2*pi/lamda;
rho=(0.545*c^2*k^2*zz)^(-3/5);
f1=k^2/(zz^2)*s1*s2;
f2=exp(j*k/(2*zz)*(s1^2-s2^2))*exp(-1/rho^2*(s1^2+s2^2))*besseli(0,2*s1*s2/rho^2);
%上面f1,f2为被积函数中与源无关部分,下面为sgb
sgb=besselj(0,alpha*s1)*besselj(0,alpha*s2)*exp(-s1^8/w^8)*exp(-s2^8/w^8);
f=f1*f2*sgb;%integrand
f=subs(f,{zz,lamda,alpha,c,w},{Z,1550*10^-9,600,C,0.05});%这是一个矢量;
disp(f);
intValue=ones(size(Z));
for i=1:length(Z)
disp(i);
ff=inline(f(i));
intValue(i)=quad2dg(ff,0,0.1,0,0.1,10^-20);
end;
plot(Z',abs(intValue))
总是显示有奇点,请教大侠,如何处理?问题1,如何确定奇点的位置?
问题2 为了完成积分确定了位置然后怎样处理呢?
注:quad2dg可以改为dblquad,谢谢!请高手费心!