sir how to use p loop in this equation 34
so for p llop we have to use eq 34
this is program
a(1)=1;
b(1)=1;
c(1)=1;
k2=1;
T=0.1;
r=2.8
k1=0.3;
for i=1:20
for p=1:i-1
a(i+1)=(1/(i+1))*(a(i)-b(i))*k2;
S=S+a(p)*c(i-p);
Z=Z+a(p)*b(i-p);
b(i+1)=(1/(i+1))*(r*a(i)-T*S-T*b(i));
c(i+1)=(1/(i+1))*(T*Z-k1*c(i));
end
end

 Réponse acceptée

Torsten
Torsten le 8 Mar 2022
Modifié(e) : Torsten le 8 Mar 2022
There are some mistakes in your article.
The code is based on
% Set parameters
X0 = 0.0;
Y0 = 1.0;
Z0 = 0.0;
N = 120;
SIGMA = 10;
R = 28;
B = 8/3;
% Initialize coefficient arrays
a = zeros(N+1,1);
b = zeros(N+1,1);
c = zeros(N+1,1);
a(1) = X0;
b(1) = Y0;
c(1) = Z0;
% Compute Taylor coefficients
for k = 1:N
SB = 0.0;
SC = 0.0;
for i= 0:k-1
SB = SB + a(i+1)*c(k-i);
SC = SC + a(i+1)*b(k-i);
end
a(k+1) = (SIGMA*(b(k) - a(k)))/k;
b(k+1) = ((R*a(k) - b(k)) - SB)/k ;
c(k+1) = (-B*c(k) + SC)/k ;
end
% Evaluate Taylor series at T
T = 0:0.01:0.2;
X = zeros(size(T));
Y = zeros(size(T));
Z = zeros(size(T));
for i = 1:numel(T)
tau = T(i);
x = a(N+1);
y = b(N+1);
z = c(N+1);
for k = N:-1:1
x = x*tau + a(k);
y = y*tau + b(k);
z = z*tau + c(k);
end
X(i) = x;
Y(i) = y;
Z(i) = z;
end
figure(1)
plot(T,X,'Color','red')
hold on
figure(2)
plot(T,Y,'Color','red')
hold on
figure(3)
plot(T,Z,'Color','red')
hold on
% Compare with ODE solution
fun = @(t,x)[SIGMA*(x(2)-x(1));x(1)*(R-x(3))-x(2);x(1)*x(2)-B*x(3)];
x0 = [X0;Y0;Z0];
tspan = T;
options = odeset ("RelTol", 1e-8, "AbsTol", 1e-8);
[TT,XX] = ode15s(fun,tspan,x0,options);
figure(1)
plot(TT,XX(:,1),'Color','blue')
figure(2)
plot(TT,XX(:,2),'Color','blue')
figure(3)
plot(TT,XX(:,3),'Color','blue')

14 commentaires

shiv gaur
shiv gaur le 8 Mar 2022
dear sir your approach is appreciable you are very close to my problem I have attach the paper of lorenz using power series using equation 32 to 34 very nearly to taylor so sir you are req to fig 5
Torsten
Torsten le 8 Mar 2022
The code is the correct implementation of a power-series solution to the Lorenz equations.
So now it's up to you to do with the code what you like to do.
shiv gaur
shiv gaur le 8 Mar 2022
but why 3d plot is not butterfly and t vs Y is not matching lorenz
shiv gaur
shiv gaur le 8 Mar 2022
thanks for answer pl apply to this to my paper in later part
Torsten
Torsten le 8 Mar 2022
Modifié(e) : Torsten le 8 Mar 2022
Here is your butterfly:
% Set parameters
Tstart = 0.0;
Tend = 1000.0;
Nt = 20000;
dT = (Tend-Tstart)/Nt;
X0 = 0.0;
Y0 = 1.0;
Z0 = 0.0;
N = 20;
SIGMA = 10;
R = 28;
B = 8/3;
%
% Initialize coefficient arrays
T = zeros(Nt+1,1);
X = zeros(Nt+1,1);
Y = zeros(Nt+1,1);
Z = zeros(Nt+1,1);
a = zeros(N+1,1);
b = zeros(N+1,1);
c = zeros(N+1,1);
T(1) = 0.0;
X(1) = X0;
Y(1) = Y0;
Z(1) = Z0;
for j = 2:Nt+1
% Compute Taylor coefficients
a(1) = X(j-1);
b(1) = Y(j-1);
c(1) = Z(j-1);
for k = 1:N
SB = 0.0;
SC = 0.0;
for i= 0:k-1
SB = SB + a(i+1)*c(k-i);
SC = SC + a(i+1)*b(k-i);
end
a(k+1) = (SIGMA*(b(k) - a(k)))/k;
b(k+1) = ((R*a(k) - b(k)) - SB)/k ;
c(k+1) = (-B*c(k) + SC)/k ;
end
% Evaluate Taylor series at actual T
x = a(N+1);
y = b(N+1);
z = c(N+1);
for k = N:-1:1
x = x*dT + a(k);
y = y*dT + b(k);
z = z*dT + c(k);
end
%x = a(1);
%y = b(1);
%z = c(1);
%for k = 2:N+1
% x = x + a(k)*dT^(k-1);
% y = y + b(k)*dT^(k-1);
% z = z + c(k)*dT^(k-1);
%end
% Prepare for T = T + dT
T(j) = T(j-1) + dT;
X(j) = x;
Y(j) = y;
Z(j) = z;
end
figure(1)
plot(T,X,'Color','red')
hold on
figure(2)
plot(T,Y,'Color','red')
hold on
figure(3)
plot(T,Z,'Color','red')
hold on
figure(4)
plot3(X,Y,Z,'Color','red')
hold on
% Compare with ODE solution
fun = @(t,x)[SIGMA*(x(2)-x(1));x(1)*(R-x(3))-x(2);x(1)*x(2)-B*x(3)];
x0 = [X0;Y0;Z0];
tspan = [Tstart,Tend];
options = odeset ("RelTol", 1e-7, "AbsTol", 1e-7, "InitialStep", 1e-9);
[TT,XX] = ode45(fun,tspan,x0,options);
figure(1)
plot(TT,XX(:,1),'Color','blue')
figure(2)
plot(TT,XX(:,2),'Color','blue')
figure(3)
plot(TT,XX(:,3),'Color','blue')
figure(4)
plot3(XX(:,1),XX(:,2),XX(:,3),'Color','blue')
Torsten
Torsten le 8 Mar 2022
thanks for answer pl apply to this to my paper in later part
As I mentionned, there are errors in your paper.
You'd better use the one I gave you the link for.
shiv gaur
shiv gaur le 8 Mar 2022
sir how I use u=summation (a(i)*tau^i) so from calculating the taylor coeff is same as cauchy product
which you have calculate and finally my equation is
u=u+ a(i)*tau^i;
v=v+b(i)*tau^i;
w=w+c(i)*tau^i;
so plor3d(u,v,w) looks like butterfly as in fig5
shiv gaur
shiv gaur le 8 Mar 2022
this is last help
Torsten
Torsten le 8 Mar 2022
Modifié(e) : Torsten le 8 Mar 2022
sir how I use u=summation (a(i)*tau^i) so from calculating the taylor coeff is same as cauchy product
which you have calculate and finally my equation is
u=u+ a(i)*tau^i;
v=v+b(i)*tau^i;
w=w+c(i)*tau^i;
so plor3d(u,v,w) looks like butterfly as in fig5
% Evaluate Taylor series at actual T
x = a(N+1);
y = b(N+1);
z = c(N+1);
for k = N:-1:1
x = x*dT + a(k);
y = y*dT + b(k);
z = z*dT + c(k);
end
is the part of the code where you do
x = x + a(i)*tau^i
y = y + b(i)*tau^i
z = z + b(i)*tau^i
It's done with Horner's scheme for polynomials because it uses less multiplications.
You might want to compare with
x = a(1);
y = b(1);
z = c(1);
for k = 2:N+1
x = x + a(k)*dT^(k-1);
y = y + b(k)*dT^(k-1);
z = z + c(k)*dT^(k-1);
end
Both methods should give more or less the same result, but yours is less efficient.
shiv gaur
shiv gaur le 8 Mar 2022
in this result is not mat
Torsten
Torsten le 8 Mar 2022
Pardon ?
shiv gaur
shiv gaur le 8 Mar 2022
in this result is not matching so pl need ful
shiv gaur
shiv gaur le 8 Mar 2022
shiv gaur
shiv gaur le 8 Mar 2022
thanks very much

Connectez-vous pour commenter.

Plus de réponses (0)

Catégories

En savoir plus sur MATLAB dans Centre d'aide et File Exchange

Produits

Version

R2021b

Tags

Community Treasure Hunt

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

Start Hunting!

Translated by