How to integrate with 3 variables?

1 次查看（过去 30 天）
ylcnt
ylcnt 2017-10-21
回答： Walter Roberson 2017-10-24
Hi,
I have to take the integration of the complex function. When i run, it doesn't give any result. Could anybody please help me to achieve it? And i need a general view to calculate integration of complex functions. Thanks in advance.
Here is my code :
function Fig1_pc_
clc
global z; global kppx; global kppy; global u;
alphas=0.03;
N=3;
lambda=1.55e-6;
k=2*pi/lambda;
w=-3;
eta=1e-3;
vsc=1e-6;
eps=(vsc^3)/(eta^4);
XT=1e-4;
L=50;
A_T=1.863e-2; A_S=1.9e-4; A_TS=9.41e-3;
mux=1; muy=2;
j=1;
I0L1=0;
for SNRx=0:2:20
SNR=10^(SNRx/10);
SNRn(j)=SNRx;
for l1=1:1:N
alphasl1=alphas/sqrt(l1);
alpha_l1=1/(2*k*alphasl1*alphasl1);
A_l1=(((-1)^(l1-1))*factorial(N)/(factorial(l1)*factorial(N-l1)*N));
for l2=1:1:N
alphasl2=alphas/sqrt(l2);
alpha_l2=1/(2*k*alphasl2*alphasl2);
A_l2=(((-1)^(l2-1))*factorial(N)/(factorial(l2)*factorial(N-l2)*N));
for ld1=1:1:N
alphasld1=alphas/sqrt(ld1);
alpha_ld1=1/(2*k*alphasld1*alphasld1);
A_ld1=(((-1)^(ld1-1))*factorial(N)/(factorial(ld1)*factorial(N-ld1)*N));
for ld2=1:1:N
alphasld2=alphas/sqrt(ld2);
alpha_ld2=1/(2*k*alphasld2*alphasld2);
A_ld2=(((-1)^(ld2-1))*factorial(N)/(factorial(ld2)*factorial(N-ld2)*N));
D2=A_ld1*A_ld2*(1/(0.5+alpha_ld1*1i*L))*(1/(0.5+alpha_ld2*1i*L))
D2A=A_ld1*A_ld2*(1/(0.5+alpha_ld1*1i*L))*(1/(0.5-alpha_ld2*1i*L))
end
end
NR=@(z,kppx,kppy) ((-(A_l1*A_l2*(1/(0.5+alpha_l1*1i*L))*(1/(0.5+alpha_l2*1i*L))*exp(-0.5i*(L-z).*(kppx.*kppx+kppy.*kppy).*((0.5+alpha_l1*1i*z)/(0.5+alpha_l1*1i*L)+(0.5+alpha_l2*1i*z)/(0.5+alpha_l2*1i*L))))/D2+...
(A_l1*A_l2*(1/(0.5+alpha_l1*1i*L))*(1/(0.5-alpha_l2*1i*L))*exp(-0.5i*(L-z).*(kppx.*kppx+kppy.*kppy).*((0.5+alpha_l1*1i*z)/(0.5+alpha_l1*1i*L)-(0.5-alpha_l2*1i*z)/(0.5-alpha_l2*1i*L))))/D2A).*...
(0.388*(1e-8)*mux*muy*XT*(eps^(-1/3))*(w^(-2))*((mux*mux*kppx.*kppx+muy*muy*kppy.*kppy).^(-11/6)).*...
(1+2.35*(eta^(2/3))*((mux*mux*kppx.*kppx+muy*muy*kppy.*kppy).^(1/3))).*...
(w*w*exp(-A_T*(8.284*(eta^(4/3))*((mux*mux*kppx.*kppx+muy*muy*kppy.*kppy).^(2/3))))+...
exp(-A_S*(8.284*(eta^(4/3))*((mux*mux*kppx.*kppx+muy*muy*kppy.*kppy).^(2/3))))-...
2*w*exp(-A_TS*(8.284*(eta^(4/3))*((mux*mux*kppx.*kppx+muy*muy*kppy.*kppy).^(2/3)))))));
m2=4*pi*real(integral3(NR1,0,L,0,inf,0,inf))
BERx = @(u) ((0.5./(sqrt(m2)*u*sqrt(2*pi))).*exp(-((log(u)+0.5*m2).^2)/(2*m2)).*erfc(SNR*u/sqrt(8)));
BER1=I0L1+integral(BERx,0,inf);
end
end
trns1(j)=BER1
j=j+1;
end
hold on
plot(SNRn,trns1,'--k','LineWidth',2);
hold off
xlabel( '\it\varsigma\rm','FontSize',24,'FontSize',24,'FontName','Times New Roman'); set(gca,'FontSize',24,'FontName','Times New Roman');
ylabel('Transmittance, \tau','FontSize',24,'FontName','Times New Roman');set(gca,'FontSize',24,'FontName','Times New Roman');
  2 个评论
Walter Roberson
Walter Roberson 2017-10-21
You define a function handle NR which you do not use, but you apply integral3 to NR1 which you have not defined
ylcnt
ylcnt 2017-10-22
There is a typo error. The above question is valid when NR is used.

请先登录，再进行评论。

请先登录，再回答此问题。

回答（1 个）

Walter Roberson
Walter Roberson 2017-10-24
If you have the symbolic toolbox, then you can reduce NR by one variable by integrating symbolically from Z = 0 to L.
Unfortunately, the real part of the resulting function goes to infinity as either kppx or kppy go to 0, making the fine details near 0 to be quite important. The adaptive numeric integral routines take a long long time. It is a fairly messy complex function.

类别

在 Help Center 和 File Exchange 中查找有关 Numeric Solvers 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by