Is there a way to improve the code speed?

Question

Kai Wang 2023-10-28

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2039821-is-there-a-way-to-improve-the-code-speed

评论： Kai Wang 2023-11-29

APLL.slx

在 MATLAB Online 中打开

Hi,

Here below is the code I made for plotting the output frequency over time of a classic Phase Lock Loop module.

The code is horrible in speed, much slower than the simulink model which I also attached here.

My team leader coded similar thing in Julia but in a quite different way, while the speed is ultra fast. Within less than 3 seconds can he get the plot. I have realized there is something like VectorContinuousCallback in the Julia code for solving ODEs and many other special things.

Can we improve the code speed in matlab although without such functions?

clear;close all;clc;
tr(1)=5e-6;tref=25e-9;tdd=250e-12;ip=60e-6;in=60e-6;c0=40e-12;c1=360e-12;c2=20e-12;r1=14000;r2=1000;f0=14.16e9;kvco=300e6;n=360;
m=2500;ic=zeros(1,m);tc=sym('tc',[1 m]);p=sym('p',[1 m]);q=sym('q',[1 m]);s=sym('s',[1 m]);phy=sym('phy',[1 m]);syms f(t) g(t) h(t) vc0(t) vc1(t) vc2(t) t x a b;
ode=[r1*c1*diff(vc1,t)+vc1==vc0, r2*c2*diff(vc2,t)+vc2==vc0, c0*diff(vc0,t)+c1*diff(vc1,t)+c2*diff(vc2,t)==0];conds=[vc0(0)==0,vc1(0)==0,vc2(0)==0];f(t)=dsolve(ode,conds).vc0;g(t)=dsolve(ode,conds).vc1;h(t)=dsolve(ode,conds).vc2;
t0=vpasolve(mod(0,1)+int((f0+kvco*h)/n,0,x)==1,x,[0,inf]);
if t0==tr(1)
    ic(1)=ip-in;tc(1)=tr(1);p(1)=f(tc(1));q(1)=g(tc(1));s(1)=h(tc(1));phy(1)=0;
elseif t0>tr(1)
    ic(1)=ip;tc(1)=tr(1);p(1)=f(tc(1));q(1)=g(tc(1));s(1)=h(tc(1));phy(1)=int((f0+kvco*h)/n,0,tc(1));
elseif t0<tr(1)
    ic(1)=-in;tc(1)=t0;p(1)=f(tc(1));q(1)=g(tc(1));s(1)=h(tc(1));phy(1)=0;
end
fplot(f0+kvco*h,[0,double(tc(1))]);hold on
for i=2:m
    if ic(i-1)==-in
        ic(i)=ip-in;
        if tc(i-1)-tr(1)>=0
            tc(i)=tr(1)+tref*(floor((tc(i-1)-tr(1))/tref)+1);
        else
            tc(i)=tr(1);
        end
        ode=[r1*c1*diff(vc1,t)+vc1==vc0, r2*c2*diff(vc2,t)+vc2==vc0, c0*diff(vc0,t)+c1*diff(vc1,t)+c2*diff(vc2,t)==ic(i-1)];conds=[vc0(0)==p(i-1),vc1(0)==q(i-1),vc2(0)==s(i-1)];f(t)=dsolve(ode,conds).vc0;g(t)=dsolve(ode,conds).vc1;h(t)=dsolve(ode,conds).vc2;
        p(i)=f(tc(i)-tc(i-1));q(i)=g(tc(i)-tc(i-1));s(i)=h(tc(i)-tc(i-1));phy(i)=phy(i-1)+vpaintegral((f0+kvco*h)/n,0,tc(i)-tc(i-1));
    elseif ic(i-1)==ip
        ic(i)=ip-in;
        ode=[r1*c1*diff(vc1,t)+vc1==vc0, r2*c2*diff(vc2,t)+vc2==vc0, c0*diff(vc0,t)+c1*diff(vc1,t)+c2*diff(vc2,t)==ic(i-1)];conds=[vc0(0)==p(i-1),vc1(0)==q(i-1),vc2(0)==s(i-1)];f(t)=dsolve(ode,conds).vc0;g(t)=dsolve(ode,conds).vc1;h(t)=dsolve(ode,conds).vc2;
        tc(i)=tc(i-1)+vpasolve(mod(phy(i-1),1)+int((f0+kvco*h)/n,0,x)==1,x,[0,inf]);
        p(i)=f(tc(i)-tc(i-1));q(i)=g(tc(i)-tc(i-1));s(i)=h(tc(i)-tc(i-1));phy(i)=0;
    elseif ic(i-1)==ip-in
        if tc(i-1)-tr(1)>=0
            a=tr(1)+tref*(floor((tc(i-1)-tr(1))/tref)+1);
        else
            a=tr(1);
        end
        ode=[r1*c1*diff(vc1,t)+vc1==vc0, r2*c2*diff(vc2,t)+vc2==vc0, c0*diff(vc0,t)+c1*diff(vc1,t)+c2*diff(vc2,t)==ic(i-1)];conds=[vc0(0)==p(i-1),vc1(0)==q(i-1),vc2(0)==s(i-1)];f(t)=dsolve(ode,conds).vc0;g(t)=dsolve(ode,conds).vc1;h(t)=dsolve(ode,conds).vc2;
        b=tc(i-1)+vpasolve(mod(phy(i-1),1)+int((f0+kvco*h)/n,0,x)==1,x,[0,inf]);
        if b==a
            ic(i)=ip-in;tc(i)=a;p(i)=f(tc(i)-tc(i-1));q(i)=g(tc(i)-tc(i-1));s(i)=h(tc(i)-tc(i-1));phy(i)=0;
        elseif b>a
            ic(i)=ip;tc(i)=a;p(i)=f(tc(i)-tc(i-1));q(i)=g(tc(i)-tc(i-1));s(i)=h(tc(i)-tc(i-1));phy(i)=phy(i-1)+vpaintegral((f0+kvco*h)/n,0,tc(i)-tc(i-1));
        else
            ic(i)=-in;tc(i)=b;p(i)=f(tc(i)-tc(i-1));q(i)=g(tc(i)-tc(i-1));s(i)=h(tc(i)-tc(i-1));phy(i)=0;
        end
    end
    fplot(f0+kvco*h(t-tc(i-1)),[double(tc(i-1)),double(tc(i))])
end

3 个评论
显示 1更早的评论隐藏 1更早的评论

Dyuman Joshi 2023-11-26

编辑：Dyuman Joshi 2023-11-26

It's not clear to me what you are trying to do and it's hard to follow the code. So I will ask some questions -

Any particular reason why you have defined tc, phy, p, q, and s as a symbolic variable and not numerical like ic?

The plots are added to the same figure, so they don't seem to be the bottle neck.

Have you tried using profiler on your code?

Kai Wang 2023-11-29

Thank you for your questions.

I defined those parameters as symbolic variables because I assumed they would provide better accuracy in vpasolve. Maybe it is unnecessary.

I did not know such thing as profiler. Thank you for reminding, let me check with that.

请先登录，再进行评论。

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

Answer 1

Sulaymon Eshkabilov 2023-10-29

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2039821-is-there-a-way-to-improve-the-code-speed#answer_1342781

Without digging much into your code, I'd suggest to put:

fplot(f0+kvco*h(t-tc(i-1)),[double(tc(i-1)),double(tc(i))]) outside of the loop that speeds up the simulation of your code signifcantly. Because you are plotting calculated data for over 2500 times <-- m=2500 which requires lots of time.

All the best.