How to integrate with 3 variables?

Question

ylcnt 2017-10-21

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/362542-how-to-integrate-with-3-variables

回答： 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 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 2017-10-22

There is a typo error. The above question is valid when NR is used.

请先登录，再进行评论。

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

Answer 1

Walter Roberson 2017-10-24

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/362542-how-to-integrate-with-3-variables#answer_287544

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.