# Untitled

Jan 30th, 2021
1,645
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. clc;
2. clear;
3. xdel(winsid());
4.
5. // Initializing variables
6. h = 0.05;
7. x = 3;
8. y = 1;
9. z = 2;
10. xt = 0.2;
11. yt = 0.1;
12. zt = 0.05;
13.
14. // Start values equals zero
15. for a = 1:2000
16.  x_st = 0;
17.  y_st = 0;
18.  z_st = 0;
19.  xt_st = 0;
20.  yt_st = 0;
21.  zt_st = 0;
22. end
23.
24. // myval equals zero in all dimensions
25. for i = 1:2000
26.  for j = 1:3
27.     myvalue = 0;
28.  end;
29. end
30.
31. // Using RK4 method to solve equations for 2000 points
32. for i = 1:2000
33.  myvalue(i,1) = y + z;
34.  myvalue(i,2) = -x + 0.5*y;
35.  myvalue(i,3) = x^2 - z;
36.
37.  kx1 = h*(y + z);
38.  ky1 = h*(-x + 0.5*y);
39.  kz1 = h*(x^2 - z);
40.
41.  kx2 = h*((y+ky1)/2 + (z+kz1)/2);
42.  ky2 = h*((-x+kx1)/2 + 0.5*((y+ky1)/2));
43.  kz2 = h*(((x+kx1)^2)/2 - (z+kz1)/2);
44.
45.  kx3 = h*((y+ky2)/2 + (z+kz2)/2);
46.  ky3 = h*((-x+kx2)/2 + 0.5*((y+ky2)/2));
47.  kz3 = h*(((x+kx2)^2)/2 - (z+kz2)/2);
48.
49.  kx4 = h*((y+ky3) + (z+kz3));
50.  ky4 = h*((-x+kx3) + 0.5*(y+ky3));
51.  kz4 = h*((x+kx3)^2 - (z+kz3));
52.
53.  xt = (kx1+2*kx2+2*kx3+kx4)/6 + xt;
54.  yt = (ky1+2*ky2+2*ky3+ky4)/6 + yt;
55.  zt = (kz1+2*kz2+2*kz3+kz4)/6 + zt;
56.
57.  xt_st(i) = xt;
58.  yt_st(i) = yt;
59.  zt_st(i) = zt;
60.
61.  x_st(i) = x;
62.  y_st(i) = y;
63.  z_st(i) = z;
64.
65.  x = x + h;
66.  y = y + h;
67.  z = z + h;
68.
69.  if i < 51
70.      mprintf("x(t) = %f, y(t) = %f, z(t) = %f\n", xt_st(i), yt_st(i), zt_st(i));
71.  end;
72.
73. end;
74.
75. plot(x_st, xt_st, 'g*-');
76. plot(y_st, yt_st, 'b*-');
77. plot(z_st, zt_st, 'r*-');
78. h1 = legend(['X(t)'; 'Y(t)'; 'Z(t)']);
79. ylabel("t");
80.
81. figure();
82. param3d(xt_st, yt_st, zt_st);
83. h1=legend(['X(t) = f(Y(t), Z(t))']);
84.
85. //graphs after compression
86. minValue = min(xt_st);
87. maxValue = max(yt_st);
88. for i = 1:2000
89.  xt_st(i) = ((xt_st(i)-minValue) * 2)/(maxValue-minValue) - 1;
90.  yt_st(i) = ((yt_st(i)-minValue) * 2)/(maxValue-minValue) - 1;
91.  zt_st(i) = ((zt_st(i)-minValue) * 2)/(maxValue-minValue) - 1;
92. end
93. figure();
94.
95. plot(x_st, xt_st,'g*-');
96. plot(y_st, yt_st,'b*-');
97. plot(z_st, zt_st,'r*-');
98. h1=legend(['X(t)'; 'Y(t)'; 'Z(t)']);
99. ylabel("t");
100.
101. figure();
102. param3d(xt_st, yt_st, zt_st);