how can improve my code performanece, it uses a syms vector, i tried to delete y=t= sym(zeros(1, m + 3)); and it's faster but the solution of the function its wrong
1 vue (au cours des 30 derniers jours)
Afficher commentaires plus anciens
function [fi, fi1, Ji, Ji1] = funcion(L,r1,s1,t1,C0,C1,C2,C3,FS,b,D0,D1,D2,AstG,S,R,DgS,ra,n)
%definimos las condiciones iniciales
m = 40;
y = sym(zeros(1, m + 3));
t = sym(zeros(1, m + 3));
y(1)=C0;
y(2)=C1;
t(1)=C2;
t(2)=C3;
for k = 3:m
y(k)=((1+FS)*(k+1)*t(k+2)+(2*R*k+2*R)*t(k+1)+(R^2*(k-1)+2*R^2)*t(k)-(2*R*k*(k+1)+2*R*(k+1))*y(k+2)+(-R^2*(k-1)*k+D0*S-b*S*n-2*R^2*k)*y(k+1)+(D1*S-2*b*S*n*R)*y(k)+(D2*S-b*S*n*R^2)*y(k-1))/((1+DgS)*(k+1)*(k+2));
t(k)=(-(4*S*R*k*(k+1)+4*S*R*(k+1))*t(k+2)+(1+FS-b*S*ra-6*S*R^2*(k-1)*k-12*S*R^2*k)*t(k+1)+(-4*S*R^3*(k-2)*(k-1)+2*R-4*b*S*ra*R-12*S*R^3*(k-1))*t(k)+(-S*R^4*(k-3)*(k-2)+R^2-6*b*S*ra*R^2-4*S*R^4*(k-2))*t(k-1)-(4*b*S*ra*R^3)*t(k-2)-(b*S*ra*R^4)*t(k-3)-(1+FS)*(k+1)*y(k+2)-2*R*k*y(k+1)-R^2*(k-1)*y(k))/(S*(k+1)*(k+2));
end
j = 0:m+2;
fi = sum(y);
fi1 = y * j';
Ji = sum(t);
Ji1 = t * j';
end
0 commentaires
Réponses (2)
Walter Roberson
le 1 Juin 2023
t(k)=(-(4*S*R*k*(k+1)+4*S*R*(k+1))*t(k+2)+(1+FS-b*S*ra-6*S*R^2*(k-1)*k-12*S*R^2*k)*t(k+1)+(-4*S*R^3*(k-2)*(k-1)+2*R-4*b*S*ra*R-12*S*R^3*(k-1))*t(k)+(-S*R^4*(k-3)*(k-2)+R^2-6*b*S*ra*R^2-4*S*R^4*(k-2))*t(k-1)-(4*b*S*ra*R^3)*t(k-2)-(b*S*ra*R^4)*t(k-3)-(1+FS)*(k+1)*y(k+2)-2*R*k*y(k+1)-R^2*(k-1)*y(k))/(S*(k+1)*(k+2));
% ^^^^^^
I pointed to a place in the code where you use t(k-3) when k starts from 3.
6 commentaires
Walter Roberson
le 5 Juin 2023
According to my tests, a pure numeric solution (no syms at all) is pretty much 100 times faster, and is accurate to one part in
Side note: my tests with lower m values showed that at least for these random inputs, the calculations are divergent, increasing in value moderately rapidly with larger m values.
To do further timing and accuracy tests, we will need your inputs.
format long g
rng(655321)
values = num2cell(randn(1, 19) * 10);
[L,r1,s1,t1,C0,C1,C2,C3,FS,b,D0,D1,D2,AstG,S,R,DgS,ra,n] = values{:};
tic
[OUT1, OUT2, OUT3, OUT4] = funcion(L,r1,s1,t1,C0,C1,C2,C3,FS,b,D0,D1,D2,AstG,S,R,DgS,ra,n);
OUT1n = double(OUT1)
OUT2n = double(OUT2)
OUT3n = double(OUT3)
OUT4n = double(OUT4)
toc
tic
[OUT1_num, OUT2_num, OUT3_num, OUT4_num] = funcion_num(L,r1,s1,t1,C0,C1,C2,C3,FS,b,D0,D1,D2,AstG,S,R,DgS,ra,n);
OUT1n_num = double(OUT1_num)
OUT2n_num = double(OUT2_num)
OUT3n_num = double(OUT3_num)
OUT4n_num = double(OUT4_num)
toc
OUT1n - OUT1n_num, ans/OUT1n
OUT2n - OUT2n_num, ans/OUT2n
OUT3n - OUT3n_num, ans/OUT3n
OUT4n - OUT4n_num, ans/OUT4n
values_sym = num2cell(sym(values));
[L,r1,s1,t1,C0,C1,C2,C3,FS,b,D0,D1,D2,AstG,S,R,DgS,ra,n] = values{:};
tic
[OUT1s, OUT2s, OUT3s, OUT4s] = funcion(L,r1,s1,t1,C0,C1,C2,C3,FS,b,D0,D1,D2,AstG,S,R,DgS,ra,n);
OUT1ns = double(OUT1s)
OUT2ns = double(OUT2s)
OUT3ns = double(OUT3s)
OUT4ns = double(OUT4s)
toc
OUT1n - OUT1ns
OUT2n - OUT2ns
OUT3n - OUT3ns
OUT4n - OUT4ns
function [fi, fi1, Ji, Ji1] = funcion(L,r1,s1,t1,C0,C1,C2,C3,FS,b,D0,D1,D2,AstG,S,R,DgS,ra,n)
% syms x ; % Definimos las variables simbólicas
%definimos las condiciones iniciales
m = 100;
% mi=m+3; %# de terminos
% % Definimos una variable "y" como un vector de longitud m n función de la variable simbólica L.
y = sym(zeros(1,m));
t = sym(zeros(1,m));
y(1)=C0;
y(2)=C1;
t(1)=C2;
t(2)=C3;
y(3)=((1+FS)*t(2)+(2*R)*t(1)-(2*R*(1))*y(2)+(D0*S-b*S*n)*y(1)+L*r1/AstG)/((1+DgS)*2);
t(3)=(-(4*S*R*(1))*t(2)+(1+FS-b*S*ra)*t(1)-(1+FS)*(1)*y(2))/(S*2);
y(4)=((1+FS)*(2)*t(3)+(2*R*1+2*R)*t(2)+(2*R^2)*t(1)-(2*R*1*(2)+2*R*(2))*y(3)+(D0*S-b*S*n-2*R^2*1)*y(2)+(D1*S-2*b*S*n*R)*y(1)+L^2*s1/AstG)/((1+DgS)*6);
t(4)=(-(4*S*R*1*(2)+4*S*R*(2))*t(3)+(1+FS-b*S*ra-12*S*R^2*1)*t(2)+(2*R-4*b*S*ra*R)*t(1)-(1+FS)*(2)*y(3)-2*R*1*y(2))/(S*6);
y(5)=((1+FS)*(3)*t(4)+(2*R*2+2*R)*t(3)+(R^2*(1)+2*R^2)*t(2)-(2*R*2*(3)+2*R*(3))*y(4)+(-R^2*(1)*2+D0*S-b*S*n-2*R^2*2)*y(3)+(D1*S-2*b*S*n*R)*y(2)+(D2*S-b*S*n*R^2)*y(1)+(t1*L^3)/AstG)/((1+DgS)*12);
t(5)=(-(4*S*R*2*(3)+4*S*R*(3))*t(4)+(1+FS-b*S*ra-6*S*R^2*(1)*2-12*S*R^2*2)*t(3)+(2*R-4*b*S*ra*R-12*S*R^3*(1))*t(2)+(R^2-6*b*S*ra*R^2)*t(1)-(1+FS)*(3)*y(4)-2*R*2*y(3)-R^2*(1)*y(2))/(S*12);
y(6)=((1+FS)*(4)*t(5)+(2*R*3+2*R)*t(4)+(R^2*(2)+2*R^2)*t(3)-(2*R*3*(4)+2*R*(4))*y(5)+(-R^2*(2)*3+D0*S-b*S*n-2*R^2*3)*y(4)+(D1*S-2*b*S*n*R)*y(3)+(D2*S-b*S*n*R^2)*y(2))/((1+DgS)*20);
t(6)=(-(4*S*R*3*(4)+4*S*R*(4))*t(5)+(1+FS-b*S*ra-6*S*R^2*(2)*3-12*S*R^2*3)*t(4)+(-4*S*R^3*(1)*(2)+2*R-4*b*S*ra*R-12*S*R^3*(2))*t(3)+(R^2-6*b*S*ra*R^2-4*S*R^4*(1))*t(2)-(4*b*S*ra*R^3)*t(1)-(1+FS)*(4)*y(5)-2*R*3*y(4)-R^2*(2)*y(3))/(S*20);
for k = 4:m
y(k+3)=((1+FS)*(k+1)*t(k+2)+(2*R*k+2*R)*t(k+1)+(R^2*(k-1)+2*R^2)*t(k)-(2*R*k*(k+1)+2*R*(k+1))*y(k+2)+(-R^2*(k-1)*k+D0*S-b*S*n-2*R^2*k)*y(k+1)+(D1*S-2*b*S*n*R)*y(k)+(D2*S-b*S*n*R^2)*y(k-1))/((1+DgS)*(k+1)*(k+2));
t(k+3)=(-(4*S*R*k*(k+1)+4*S*R*(k+1))*t(k+2)+(1+FS-b*S*ra-6*S*R^2*(k-1)*k-12*S*R^2*k)*t(k+1)+(-4*S*R^3*(k-2)*(k-1)+2*R-4*b*S*ra*R-12*S*R^3*(k-1))*t(k)+(-S*R^4*(k-3)*(k-2)+R^2-6*b*S*ra*R^2-4*S*R^4*(k-2))*t(k-1)-(4*b*S*ra*R^3)*t(k-2)-(b*S*ra*R^4)*t(k-3)-(1+FS)*(k+1)*y(k+2)-2*R*k*y(k+1)-R^2*(k-1)*y(k))/(S*(k+1)*(k+2));
end
j = 0:m+2;
fi = sum(y);
fi1 = y * j';
Ji = sum(t);
Ji1 = t * j';
end
function [fi, fi1, Ji, Ji1] = funcion_num(L,r1,s1,t1,C0,C1,C2,C3,FS,b,D0,D1,D2,AstG,S,R,DgS,ra,n)
% syms x ; % Definimos las variables simbólicas
%definimos las condiciones iniciales
m = 100;
% mi=m+3; %# de terminos
% % Definimos una variable "y" como un vector de longitud m n función de la variable simbólica L.
y = (zeros(1,m));
t = (zeros(1,m));
y(1)=C0;
y(2)=C1;
t(1)=C2;
t(2)=C3;
y(3)=((1+FS)*t(2)+(2*R)*t(1)-(2*R*(1))*y(2)+(D0*S-b*S*n)*y(1)+L*r1/AstG)/((1+DgS)*2);
t(3)=(-(4*S*R*(1))*t(2)+(1+FS-b*S*ra)*t(1)-(1+FS)*(1)*y(2))/(S*2);
y(4)=((1+FS)*(2)*t(3)+(2*R*1+2*R)*t(2)+(2*R^2)*t(1)-(2*R*1*(2)+2*R*(2))*y(3)+(D0*S-b*S*n-2*R^2*1)*y(2)+(D1*S-2*b*S*n*R)*y(1)+L^2*s1/AstG)/((1+DgS)*6);
t(4)=(-(4*S*R*1*(2)+4*S*R*(2))*t(3)+(1+FS-b*S*ra-12*S*R^2*1)*t(2)+(2*R-4*b*S*ra*R)*t(1)-(1+FS)*(2)*y(3)-2*R*1*y(2))/(S*6);
y(5)=((1+FS)*(3)*t(4)+(2*R*2+2*R)*t(3)+(R^2*(1)+2*R^2)*t(2)-(2*R*2*(3)+2*R*(3))*y(4)+(-R^2*(1)*2+D0*S-b*S*n-2*R^2*2)*y(3)+(D1*S-2*b*S*n*R)*y(2)+(D2*S-b*S*n*R^2)*y(1)+(t1*L^3)/AstG)/((1+DgS)*12);
t(5)=(-(4*S*R*2*(3)+4*S*R*(3))*t(4)+(1+FS-b*S*ra-6*S*R^2*(1)*2-12*S*R^2*2)*t(3)+(2*R-4*b*S*ra*R-12*S*R^3*(1))*t(2)+(R^2-6*b*S*ra*R^2)*t(1)-(1+FS)*(3)*y(4)-2*R*2*y(3)-R^2*(1)*y(2))/(S*12);
y(6)=((1+FS)*(4)*t(5)+(2*R*3+2*R)*t(4)+(R^2*(2)+2*R^2)*t(3)-(2*R*3*(4)+2*R*(4))*y(5)+(-R^2*(2)*3+D0*S-b*S*n-2*R^2*3)*y(4)+(D1*S-2*b*S*n*R)*y(3)+(D2*S-b*S*n*R^2)*y(2))/((1+DgS)*20);
t(6)=(-(4*S*R*3*(4)+4*S*R*(4))*t(5)+(1+FS-b*S*ra-6*S*R^2*(2)*3-12*S*R^2*3)*t(4)+(-4*S*R^3*(1)*(2)+2*R-4*b*S*ra*R-12*S*R^3*(2))*t(3)+(R^2-6*b*S*ra*R^2-4*S*R^4*(1))*t(2)-(4*b*S*ra*R^3)*t(1)-(1+FS)*(4)*y(5)-2*R*3*y(4)-R^2*(2)*y(3))/(S*20);
for k = 4:m
y(k+3)=((1+FS)*(k+1)*t(k+2)+(2*R*k+2*R)*t(k+1)+(R^2*(k-1)+2*R^2)*t(k)-(2*R*k*(k+1)+2*R*(k+1))*y(k+2)+(-R^2*(k-1)*k+D0*S-b*S*n-2*R^2*k)*y(k+1)+(D1*S-2*b*S*n*R)*y(k)+(D2*S-b*S*n*R^2)*y(k-1))/((1+DgS)*(k+1)*(k+2));
t(k+3)=(-(4*S*R*k*(k+1)+4*S*R*(k+1))*t(k+2)+(1+FS-b*S*ra-6*S*R^2*(k-1)*k-12*S*R^2*k)*t(k+1)+(-4*S*R^3*(k-2)*(k-1)+2*R-4*b*S*ra*R-12*S*R^3*(k-1))*t(k)+(-S*R^4*(k-3)*(k-2)+R^2-6*b*S*ra*R^2-4*S*R^4*(k-2))*t(k-1)-(4*b*S*ra*R^3)*t(k-2)-(b*S*ra*R^4)*t(k-3)-(1+FS)*(k+1)*y(k+2)-2*R*k*y(k+1)-R^2*(k-1)*y(k))/(S*(k+1)*(k+2));
end
j = 0:m+2;
fi = sum(y);
fi1 = y * j';
Ji = sum(t);
Ji1 = t * j';
end
Voir également
Catégories
En savoir plus sur Number Theory dans Help Center et File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!