Skip to content

Commit 217447c

Browse files
committed
Example7 Added
1 parent 073b054 commit 217447c

File tree

9 files changed

+210
-30
lines changed

9 files changed

+210
-30
lines changed

Animator7.m

Lines changed: 121 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,121 @@
1+
function Animator7(Th, tt,Param)
2+
% Author: Mansour Torabi
3+
% Email: smtoraabi@ymail.com
4+
%%
5+
try
6+
filename = 'Pic/Anim7.gif';
7+
Hf = figure;
8+
set(Hf,'color',[1 1 1]);
9+
10+
th0 = Th(1,1); ths = Th(1,2); dx = Th(1,3);
11+
R0 = Param(1); r0 = Param(2); l0 = Param(3);
12+
13+
RotFcn = @(X,th)[cos(th), -sin(th); sin(th), cos(th)]*X;
14+
15+
%% Main Param
16+
17+
R = 3.5; r = r0*R/R0; l = l0*R/R0;
18+
R_od = R-r;
19+
20+
xRo = 5; yRo = 7;
21+
22+
%% SemiCircle and Disc
23+
24+
aa = linspace(-pi,0,100);
25+
xR = R*cos(aa) + xRo;
26+
yR = R*sin(aa) + yRo;
27+
plot(xR, yR, 'LineWidth',2, 'Color',[0,0.6,0]);
28+
hold on;
29+
30+
P_disc = [xRo, yRo] + R_od*[sin(th0), -cos(th0)];
31+
Hod = plot(P_disc(1), P_disc(2), 'ro');
32+
33+
aa = linspace(-pi,pi,100);
34+
xr = r*cos(aa) + P_disc(1);
35+
yr = r*sin(aa) + P_disc(2);
36+
37+
Hd = plot(xr, yr, 'LineWidth',5, 'Color','r');
38+
39+
% Line inside disc
40+
xld = [0, 0]; yld = [0, -r];
41+
XXld = [xld;yld];
42+
alpha = R_od * th0 / r;
43+
Xld = RotFcn(XXld, alpha) + repmat([P_disc(1); P_disc(2)],1,length(xld));
44+
45+
Hld = plot(Xld(1,:), Xld(2,:), 'r','LineWidth',2);
46+
%% Spring and Mass
47+
Ls = l + dx;
48+
49+
a = 0.1; t0 = linspace(0,Ls,1e3);
50+
x = a*sin(100/Ls*t0); y = - t0;
51+
XX = [x;y];
52+
Xrs = RotFcn(XX, ths) + repmat([P_disc(1); P_disc(2)],1,length(t0));
53+
Hs = line([P_disc(1), Xrs(1,:)],[P_disc(2), Xrs(2,:)],'linewidth',1.2);
54+
55+
Pm = R_od*[sin(th0), -cos(th0)] + Ls*[sin(ths), -cos(ths)] + [xRo, yRo];
56+
Hm = plot(Pm(1), Pm(2), 'O', 'markersize', 12, ...
57+
'linewidth', 2, 'color', 'r', 'markerfacecolor', 'r');
58+
%% Axis
59+
60+
box off; axis equal
61+
set(gca,'xlim',[0 10],'ylim',[0 9],'xtick',[], 'ytick', [])
62+
63+
%% Text
64+
Txt = sprintf('Time: %0.2f sec', tt(1));
65+
Htxt = text(3.5, 9, Txt);
66+
set(Htxt, 'fontsize', 16, 'fontweight', 'bold');
67+
68+
frame = getframe(Hf); im = frame2im(frame); [imind,cm] = rgb2ind(im,256);
69+
imwrite(imind,cm,filename,'gif', 'Loopcount',inf, 'DelayTime',0.05);
70+
71+
for i = 1:size(Th,1)
72+
pause(0.00001)
73+
th0 = Th(i,1); ths = Th(i,2); dx = Th(i,3);
74+
%% SemiCircle and Disc
75+
76+
P_disc = [xRo, yRo] + R_od*[sin(th0), -cos(th0)];
77+
set(Hod, 'xdata',P_disc(1), 'ydata', P_disc(2))
78+
79+
aa = linspace(-pi,pi,100);
80+
xr = r*cos(aa) + P_disc(1);
81+
yr = r*sin(aa) + P_disc(2);
82+
83+
set(Hd, 'xdata',xr, 'ydata', yr);
84+
85+
% Line inside disc
86+
xld = [0, 0]; yld = [0, -r];
87+
XXld = [xld;yld];
88+
alpha = R_od * th0 / r;
89+
Xld = RotFcn(XXld, alpha) + repmat([P_disc(1); P_disc(2)],1,length(xld));
90+
set(Hld, 'xdata', Xld(1,:), 'ydata', Xld(2,:));
91+
92+
%% Spring and Mass
93+
Ls = l + dx;
94+
95+
a = 0.1; t0 = linspace(0,Ls,1e3);
96+
x = a*sin(100/Ls*t0); y = - t0;
97+
XX = [x;y];
98+
Xrs = RotFcn(XX, ths) + repmat([P_disc(1); P_disc(2)],1,length(t0));
99+
100+
set(Hs, 'xdata', Xrs(1,:), 'ydata', Xrs(2,:));
101+
102+
Pm = R_od*[sin(th0), -cos(th0)] + Ls*[sin(ths), -cos(ths)] + [xRo, yRo];
103+
set(Hm,'xdata', Pm(1), 'ydata', Pm(2));
104+
105+
%%
106+
107+
Txt = sprintf('Time: %0.2f sec', tt(i));
108+
set(Htxt, 'string', Txt);
109+
110+
111+
frame = getframe(Hf); im = frame2im(frame); [imind,cm] = rgb2ind(im,256);
112+
imwrite(imind,cm,filename,'gif', 'WriteMode','append', 'DelayTime',0.05);
113+
114+
end
115+
catch
116+
fprintf('You Forced to stop Animation!!\n\n')
117+
end
118+
119+
120+
121+

EVAL6.m

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -48,7 +48,7 @@
4848
set(H,'interpreter','latex','fontsize',18,'location','SouthWest');
4949

5050
hx = xlabel('Time (sec)'); set(hx, 'fontsize', 18);
51-
hy = ylabel('Angles (rad)- Lenghth (m)'); set(hy, 'fontsize', 18);
51+
hy = ylabel('Angles (rad)'); set(hy, 'fontsize', 18);
5252
set(gca, 'fontsize', 18);
5353
saveas(gcf, 'Pic/Ex6.png')
5454
%%

EVAL7.m

Lines changed: 22 additions & 27 deletions
Original file line numberDiff line numberDiff line change
@@ -4,54 +4,49 @@
44
%%
55
clc, clear
66

7-
syms th0 th x Dth0 Dth Dx
8-
syms k M m l R r g
7+
syms th0 ths x Dth0 Dths Dx
8+
syms R r M J m k l g
99

1010
%% Kinetic and Potential Energy
1111

12-
v1x = l1*Dth1*cos(th1) ;
13-
v1y = -l1*Dth1*sin(th1);
12+
VoM = (R-r)*[cos(th0), sin(th0)];
13+
Wd = -(R-r)*Dth0/r;
1414

15-
v2x = l2*Dth2*cos(th2) ;
16-
v2y = -l2*Dth2*sin(th2);
15+
Vm = (R-r)*Dth0*[cos(th0), sin(th0)] + (l+x)*Dths*[cos(ths), sin(ths)] + Dx*[sin(ths), -cos(ths)];
1716

18-
v1t = v1x^2 + v1y^2;
19-
v2t = v2x^2 + v2y^2;
17+
yM = -(R-r)*cos(th0);
18+
ym = yM - (l+x)*cos(ths);
2019

21-
T = 1/2*m1*v1t + 1/2*m2*v2t;
20+
T = 1/2*M*(VoM*VoM.') + 1/2*m*(Vm*Vm.') + 1/2*J*Wd^2;
2221

23-
dXX = l0 + l2*sin(th2) - l1*sin(th1);
24-
dYY = l1*cos(th1) - l2*cos(th2);
25-
dx = (dXX^2 + dYY^2)^0.5 - l3;
26-
27-
V1 = -m1*g*(l1*cos(th1)) + 1/2*k*dx^2;
28-
V2 = -m2*g*(l2*cos(th2));
29-
V = V1 + V2;
22+
V = M*g*yM + m*g*ym + 1/2*k*x^2;
3023

3124
L = T - V;
3225
%%
33-
q = [th1, th2];
34-
Dq = [Dth1, Dth2];
26+
q = [th0, ths, x];
27+
Dq = [Dth0, Dths, Dx];
3528
tt = linspace(0, 20, 500);
3629
Eq = LagrangeDynamicEqDeriver(L, q, Dq);
37-
l0n = 2; l1n = 1; l2n = 1.5; l3n = 2;
38-
[SS, xx] = DynamicEqSolver(Eq, q, Dq, [k m1 m2 l0 l1 l2 l3 g],...
39-
[20,1,3, l0n, l1n, l2n, l3n, 9.81], tt, [pi/6, pi/2.5, 0, 0]);
30+
R0 = 5; r0 = 1; l0 = 2;
31+
[SS, xx] = DynamicEqSolver(Eq, q, Dq, [R r M J m k l g],...
32+
[R0, r0, 1, 2, 3, 30, l0, 9.81], tt, [pi/3, pi/2, 0, 0, 0, 0]);
4033
%%
4134
figure;
4235
plot(tt,xx(:,1),'color',[0,0.6,0],'linewidth',2);
4336
hold on;
4437
plot(tt,xx(:,2),'b','linewidth',2);
38+
plot(tt,xx(:,3),'r','linewidth',2);
4539

46-
S1 = sprintf('$ \\theta_1$');
47-
S2 = sprintf('$ \\theta_2$');
48-
H = legend(S1, S2);
40+
S1 = sprintf('$ \\theta_0$');
41+
S2 = sprintf('$ \\theta_s$');
42+
S3 = sprintf('$x_s$');
43+
H = legend(S1, S2, S3);
4944
set(H,'interpreter','latex','fontsize',18,'location','SouthWest');
5045

5146
hx = xlabel('Time (sec)'); set(hx, 'fontsize', 18);
52-
hy = ylabel('Angles (rad)- Lenghth (m)'); set(hy, 'fontsize', 18);
47+
hy = ylabel('Angles (rad) - Length (m)'); set(hy, 'fontsize', 18);
5348
set(gca, 'fontsize', 18);
54-
saveas(gcf, 'Pic/Ex6.png')
49+
saveas(gcf, 'Pic/Ex7.png')
5550
%%
56-
Animator6(xx(:,1:2), tt, [l0n, l1n, l2n])
51+
Animator7(xx(:,1:3), tt, [R0, r0, l0])
5752

Pic/Anim6.gif

-1.47 MB
Loading

Pic/Anim7.gif

2.16 MB
Loading

Pic/Ex6.png

-4.02 KB
Loading

Pic/Ex7.png

36.8 KB
Loading

Pic/Ex7A.png

43.5 KB
Loading

README.md

Lines changed: 66 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -3,14 +3,15 @@
33

44
Using the above library, one can derive differential equations for any dynamic systems and solve response of the system for a given conditions.
55

6-
Functinality of the library has been Illustrated by the following examples:
6+
Functionality of the library has been illustrated by the following examples:
77

88
1. Double Pendulum
99
2. Spring Pendulum
1010
3. Pendulum with Spring-loaded support
1111
4. Double Pendulum with free support
1212
5. Double Spring Pendulum
13-
6. Couple Pendulum
13+
6. Coupled Pendulum
14+
7. Spring Pendulum with Rolling base inside a semicircle
1415

1516
## Example 1: Double Pendulum
1617

@@ -383,5 +384,68 @@ Slider Position, Pendulum Anlges:
383384
</tr>
384385
</table>
385386

387+
## Example 7: Spring Pendulum with Rolling base inside a semicircle
388+
389+
**Problem Definition:**
390+
<table style="width:100%">
391+
<tr>
392+
<td width="10%"> </td>
393+
<td width="80%">
394+
<img src="../master/Pic/Ex7A.png" />
395+
</td>
396+
<td width="10%"> </td>
397+
</tr>
398+
</table>
399+
400+
**How to solve:**
401+
402+
Just run the ```EVAL7.m``` to **derive equations** and solve intial condition problem:
403+
404+
### Code Usage:
405+
``` MATLAB
406+
syms th0 ths x Dth0 Dths Dx
407+
syms R r M J m k l g
408+
409+
%% Kinetic and Potential Energy
410+
411+
VoM = (R-r)*[cos(th0), sin(th0)];
412+
Wd = -(R-r)*Dth0/r;
413+
414+
Vm = (R-r)*Dth0*[cos(th0), sin(th0)] + (l+x)*Dths*[cos(ths), sin(ths)] + Dx*[sin(ths), -cos(ths)];
415+
416+
yM = -(R-r)*cos(th0);
417+
ym = yM - (l+x)*cos(ths);
418+
419+
T = 1/2*M*(VoM*VoM.') + 1/2*m*(Vm*Vm.') + 1/2*J*Wd^2;
420+
421+
V = M*g*yM + m*g*ym + 1/2*k*x^2;
422+
423+
L = T - V;
424+
%%
425+
q = [th0, ths, x];
426+
Dq = [Dth0, Dths, Dx];
427+
tt = linspace(0, 20, 500);
428+
Eq = LagrangeDynamicEqDeriver(L, q, Dq);
429+
R0 = 5; r0 = 1; l0 = 2;
430+
[SS, xx] = DynamicEqSolver(Eq, q, Dq, [R r M J m k l g],...
431+
[R0, r0, 1, 2, 3, 30, l0, 9.81], tt, [pi/3, pi/2, 0, 0, 0, 0]);
432+
433+
```
434+
Slider Position, Pendulum Anlges:
435+
<table style="width:100%">
436+
<tr>
437+
<th>Angles of spring length:</th>
438+
<th>Animated Response:</th>
439+
</tr>
440+
<tr>
441+
<td width="50%">
442+
<img src="../master/Pic/Ex7.png" />
443+
</td>
444+
<td width="50%">
445+
<img src="../master/Pic/Anim7.gif" />
446+
</td>
447+
</tr>
448+
</table>
449+
386450
# Contact
387451
Email: smtoraabi@ymail.com

0 commit comments

Comments
 (0)