Double Pendulum

[Ashish]

Hi Everyone! :)
This program traces the motion of the double pendulum to get beautiful patterns.
Press 'Space' for new settings.

Code: QB64: [Select]

1. 'Coded By Ashish on 4 March, 2018
2.  
3. _TITLE "Double Pendulum [Press Space for New Settings]"
4.  
5. SCREEN _NEWIMAGE(800, 600, 32)
6.  
7. TYPE vec2
8.     x AS SINGLE
9.     y AS SINGLE
10. END TYPE
11.  
12. TYPE pendlm
13.     pos AS vec2
14.     r AS SINGLE
15.     ang AS DOUBLE
16.     angInc AS DOUBLE
17.     angSize AS DOUBLE
18. END TYPE
19.  
20. DIM SHARED pendulum(1) AS pendlm
21. RANDOMIZE TIMER
22. change_settings:
23. pendulum(0).r = p5random(250, 350)
24. pendulum(0).pos.x = 400
25. pendulum(0).pos.y = 0
26. pendulum(0).angInc = p5random(0, _PI(2))
27. pendulum(0).angSize = p5random(_PI(.1), _PI(.6))
28.  
29. pendulum(1).angInc = p5random(0, _PI(2))
30. pendulum(1).r = p5random(100, 300)
31. pendulum(1).angSize = p5random(_PI(.1), _PI(.9))
32. px = 0: py = 0: a = 0: f = 0
33.  
34. tracer& = _NEWIMAGE(_WIDTH, _HEIGHT, 32)
35. DO
36.     CLS
37.     k& = _KEYHIT
38.     IF a < 255 THEN a = a + 1
39.     IF k& = ASC(" ") THEN _FREEIMAGE tracer&: GOTO change_settings
40.     tx1 = COS(pendulum(0).ang) * pendulum(0).r + pendulum(0).pos.x
41.     ty1 = SIN(pendulum(0).ang) * pendulum(0).r + pendulum(0).pos.y
42.  
43.     pendulum(0).ang = SIN(pendulum(0).angInc) * pendulum(0).angSize + _PI(.5)
44.     pendulum(1).pos.x = tx1: pendulum(1).pos.y = ty1
45.  
46.     tx2 = COS(pendulum(1).ang) * pendulum(1).r + pendulum(1).pos.x
47.     ty2 = SIN(pendulum(1).ang) * pendulum(1).r + pendulum(1).pos.y
48.  
49.     pendulum(1).ang = SIN(pendulum(1).angInc) * pendulum(0).angSize + pendulum(0).ang
50.  
51.     g = SIN(pendulum(0).angInc) * 128 + 128
52.     b = COS(pendulum(1).angInc) * 128 + 128
53.  
54.     _DEST tracer&
55.     COLOR _RGBA(0, g, b, a)
56.     IF f > 40 THEN
57.         thickLine px, py, tx2, ty2, 2
58.     END IF
59.     _DEST 0
60.     COLOR _RGB(255, 0, 0)
61.  
62.     LINE (pendulum(0).pos.x, pendulum(0).pos.y)-(tx1, ty1)
63.     CircleFill tx1, ty1, 30, _RGB(50, 50, 50)
64.     LINE (pendulum(1).pos.x, pendulum(1).pos.y)-(tx2, ty2)
65.     CircleFill tx2, ty2, 30, _RGB(50, 50, 50)
66.     _PUTIMAGE , tracer&
67.     _LIMIT 60
68.  
69.     _DISPLAY
70.     pendulum(0).angInc = pendulum(0).angInc + .03
71.     pendulum(1).angInc = pendulum(1).angInc + .05
72.     px = tx2: py = ty2
73.  
74.     f = f + 1
75. LOOP
76.  
77. 'By Fellippe
78. SUB thickLine (x1 AS SINGLE, y1 AS SINGLE, x2 AS SINGLE, y2 AS SINGLE, lineWeight%)
79.     DIM a AS SINGLE, x0 AS SINGLE, y0 AS SINGLE
80.     DIM prevDest AS LONG, prevColor AS _UNSIGNED LONG
81.     STATIC colorSample AS LONG
82.  
83.     IF colorSample = 0 THEN
84.         colorSample = _NEWIMAGE(1, 1, 32)
85.     END IF
86.  
87.     prevDest = _DEST
88.     prevColor = _DEFAULTCOLOR
89.     _DEST colorSample
90.     PSET (0, 0), prevColor
91.     _DEST prevDest
92.  
93.     a = _ATAN2(y2 - y1, x2 - x1)
94.     a = a + _PI / 2
95.     x0 = 0.5 * lineWeight% * COS(a)
96.     y0 = 0.5 * lineWeight% * SIN(a)
97.  
98.     _MAPTRIANGLE _SEAMLESS(0, 0)-(0, 0)-(0, 0), colorSample TO(x1 - x0, y1 - y0)-(x1 + x0, y1 + y0)-(x2 + x0, y2 + y0), , _SMOOTH
99.     _MAPTRIANGLE _SEAMLESS(0, 0)-(0, 0)-(0, 0), colorSample TO(x1 - x0, y1 - y0)-(x2 + x0, y2 + y0)-(x2 - x0, y2 - y0), , _SMOOTH
100. END SUB
101.  
102.  
103. SUB CircleFill (CX AS LONG, CY AS LONG, R AS LONG, C AS _UNSIGNED LONG)
104.     'This sub from here: http://www.[abandoned, outdated and now likely malicious qb64 dot net website - don’t go there]/forum/index.php?topic=1848.msg17254#msg17254
105.     DIM Radius AS LONG
106.     DIM RadiusError AS LONG
107.     DIM X AS LONG
108.     DIM Y AS LONG
109.  
110.     Radius = ABS(R)
111.     RadiusError = -Radius
112.     X = Radius
113.     Y = 0
114.  
115.     IF Radius = 0 THEN PSET (CX, CY), C: EXIT SUB
116.  
117.     ' Draw the middle span here so we don't draw it twice in the main loop,
118.     ' which would be a problem with blending turned on.
119.     LINE (CX - X, CY)-(CX + X, CY), C, BF
120.  
121.     WHILE X > Y
122.  
123.         RadiusError = RadiusError + Y * 2 + 1
124.  
125.         IF RadiusError >= 0 THEN
126.  
127.             IF X <> Y + 1 THEN
128.                 LINE (CX - Y, CY - X)-(CX + Y, CY - X), C, BF
129.                 LINE (CX - Y, CY + X)-(CX + Y, CY + X), C, BF
130.             END IF
131.  
132.             X = X - 1
133.             RadiusError = RadiusError - X * 2
134.  
135.         END IF
136.  
137.         Y = Y + 1
138.  
139.         LINE (CX - X, CY - Y)-(CX + X, CY - Y), C, BF
140.         LINE (CX - X, CY + Y)-(CX + X, CY + Y), C, BF
141.  
142.     WEND
143.  
144. END SUB
145.  
146. 'taken from p5js.bas
147. 'https://bit.ly/p5jsbas
148. FUNCTION p5random! (mn!, mx!)
149.     IF mn! > mx! THEN
150.         SWAP mn!, mx!
151.     END IF
152.     p5random! = RND * (mx! - mn!) + mn!
153. END FUNCTION
154.  
155.  

[_vince]

Nice work, really nice effect.  Remind me of this program:

Code: QB64: [Select]

1. defsng a-z
2. pi = 3.141593
3. r = 100
4. a1 = -pi/2
5. a2 = -pi/2
6. h=0.001
7. m=1000
8. g=1000
9. sw = 640
10. sh = 480
11.  
12. dim p1 as long
13. p1 = _newimage(sw,sh,12)
14. screen _newimage(sw,sh,12)
15.  
16. do
17.         a1p = (6/(m*r^2))*(2*m*v1-3*cos(a1-a2)*m*v2)/(16-9*(cos(a1-a2))^2)
18.         a1k1 = a1p
19.         a1k2 = a1p + h*a1k1/2
20.         a1k3 = a1p + h*a1k2/2
21.         a1k4 = a1p + h*a1k3
22.         a2p = (6/(m*r^2))*(8*m*v2-3*cos(a1-a2)*m*v1)/(16-9*(cos(a1-a2))^2)
23.         a2k1 = a2p
24.         a2k2 = a2p + h*a2k1/2
25.         a2k3 = a2p + h*a2k2/2
26.         a2k4 = a2p + h*a2k3
27.         na1=a1+h*(a1k1+2*a1k2+2*a1k3+a1k4)/6
28.         na2=a2+h*(a2k1+2*a2k2+2*a2k3+a2k4)/6
29.         v1p = -0.5*r^2*(a1p*a2p*sin(a1-a2)+3*g*sin(a1)/r)
30.         v1k1 = v1p
31.         v1k2 = v1p + h*v1k1/2
32.         v1k3 = v1p + h*v1k2/2
33.         v1k4 = v1p + h*v1k3
34.         v2p = -0.5*r^2*(-a1p*a2p*sin(a1-a2)+g*sin(a2)/r)
35.         v2k1 = v2p
36.         v2k2 = v2p + h*v2k1/2
37.         v2k3 = v2p + h*v2k2/2
38.         v2k4 = v2p + h*v2k3
39.         nv1=v1+h*(v1k1+2*v1k2+2*v1k3+v1k4)/6
40.         nv2=v2+h*(v2k1+2*v2k2+2*v2k3+v2k4)/6
41.         a1=na1
42.         a2=na2
43.         v1=nv1
44.         v2=nv2
45.  
46.         _dest p1
47.         pset (sw/2+r*cos(a1-pi/2)+r*cos(a2-pi/2), sh/2-r*sin(a1-pi/2)-r*sin(a2-pi/2)),12
48.         _dest 0
49.         _putimage, p1
50.  
51.         line (sw/2, sh/2)-step(r*cos(a1-pi/2), -r*sin(a1-pi/2))
52.         circle step(0,0),5
53.         line -step(r*cos(a2-pi/2), -r*sin(a2-pi/2))
54.         circle step(0,0),5
55.  
56.         _limit 200
57.         _display
58. loop until _keyhit=27
59. system
60.  

[bplus]

Nice V,

Ashish brings out the best in us. :)

[Ashish]

@v
Wow! Cool! Your version is very much accurate. :)