Speed up numeric integral

24 次查看（过去 30 天）

显示更早的评论

Guan Hao 2024-7-24

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2139916-speed-up-numeric-integral

编辑： Guan Hao 2024-7-25

在 MATLAB Online 中打开

Hi,everyone.I got a problem of doing numeric integral on the matrix H (which is at the bottom part of my code).

The code works fine but I would like to speed it up? Am I doing something stupid to make it slow?Thanks for your help.

clear
L=0.1;
section=50;
a=L./2;
b=L./section;
v=3e8;
f1=1.2e9;
f2=4.8e9;
f3=7.5e9;
f4=10e9;
w1=(2.*pi.*f1);
w2=(2.*pi.*f2);
w3=(2.*pi.*f3);
w4=(2.*pi.*f4);
Z01=50;
Z02=75;
a0=(log(Z02./Z01))./(2..*L);
T=a;
Ub=a;
Lb=-a;
ub=a;
lb=-a;
k=5;
syms x y w m
%%% Revised Algorithm
l=0:(2.*L)./(2.*k):L  % Approximate sections dvided
% Approximate terms
sec=numel(l)-1;
Z1=[50 55 60 65 70 80 75];    % Set the impedance of the end of each sections
ac=10;                  % Choose the accuracy of dividing a section
for s=1:1:sec
    c=(l(s+1)-l(s))./ac; 
    app_cos=0.5./ac.*log(Z1(s+1)./Z1(s)).*sin(ac./2.*(w./v).*c).*cos(((ac+1)./2).*(w./v).*c)./sin(((c.*w)./(2.*v)));
    app_sin=0.5./ac.*log(Z1(s+1)./Z1(s)).*sin(ac./2.*(w./v).*c).*sin(((ac+1)./2).*(w./v).*c)./sin(((c.*w)./(2.*v)));
    amp=0.75;
    Ga(s)=(app_cos-(1i.*amp.*app_sin));
end
    Ta=(1-(abs(Ga(1:1:end)).^2));
%%% Turn closed loop terms into matrix presentation
for t=1:1:sec
    for r=1:1:sec
        E(t,r)=Ga(t).*Ga(r);
    end
end
E=triu(E,1);
for t=1:1:sec
    for r=1:1:sec
    
        D(t,r)=abs(E(r,t));
    end
end
C1=sum(D);
% Construct approximate loop terms 
for t=1:1:sec
    D1(t)=prod(Ta(1:t))./(1+(sum(C1(1:t))));
end
D=[1 D1(2:end)];
Dm=transpose(D).*D;
% Cmn
m=transpose(1:5);
n=1:5;
g=(v.^2.*w.^2.*cos((2.*w.*(- y + x))./v).*cos(20.*sym(pi).*m.*x).*cos(20.*sym(pi).*n.*y) + 10.*sym(pi).*m.*v.^3.*w.*sin(20.*sym(pi).*m.*x).*sin((2.*w.*(- y + x))./v).*cos(20.*sym(pi).*n.*y) - 10.*sym(pi).*n.*v.^3.*w.*sin(20.*sym(pi).*n.*y).*sin((2.*w.*(- y + x))./v).*cos(20.*sym(pi).*m.*x) + 100.*m.*n.*v.^4.*sym(pi).^2.*sin(20.*sym(pi).*m.*x).*sin(20.*sym(pi).*n.*y).*cos((2.*w.*(- y + x))./v))./(4.*(w.^2 - 100.*sym(pi).^2.*m.^2.*v.^2).*(w.^2 - 100.*sym(pi).^2.*n.^2.*v.^2));
for p=1:numel(l)-1
    x1=l(p);
    x2=l(p+1);
    for q=1:numel(l)-1
        y1=l(q);
        y2=l(p+1);
        Ht{p,q}=matlabFunction(((subs(g,[x y],[x2 y2])-subs(g,[x y],[x1 y1]))*Dm(p,q)));
    end
end
Hi=Ht{1};
for m=2:1:numel(Ht)
    Hi=Hi+vpa(Ht{m},3);
end
H=integral(Hi,w1,w2,"ArrayValued",true)