Animate a pendulum

From Rosetta Code
Revision as of 12:35, 15 December 2010 by Oenone (talk | contribs) (add Ada)
Capture of the Oz version.
Task
Animate a pendulum
You are encouraged to solve this task according to the task description, using any language you may know.

One good way of making an animation is by simulating a physical system and illustrating the variables in that system using a dynamically changing graphical display. The classic such physical system is a simple gravity pendulum.

For this task, create a simple physical model of a pendulum and animate it.

Ada

This does not use a GUI, it simply animates the pendulum and prints out the positions. If you want, you can replace the output method with graphical update methods.

X and Y are relative positions of the pendulum to the anchor.

pendulum.ads: <lang Ada>generic

  type Float_Type is digits <>;
  Gravitation : Float_Type;

package Pendulums is

  type Pendulum is private;
  function New_Pendulum (Length : Float_Type;
                         Theta0 : Float_Type) return Pendulum;
  function Get_X (From : Pendulum) return Float_Type;
  function Get_Y (From : Pendulum) return Float_Type;
  procedure Update_Pendulum (Item : in out Pendulum; Time : in Duration);

private

  type Pendulum is record
     Length   : Float_Type;
     Theta    : Float_Type;
     X        : Float_Type;
     Y        : Float_Type;
     Velocity : Float_Type;
  end record;

end Pendulums;</lang>

pendulum.adb: <lang Ada>generic

  type Float_Type is digits <>;
  Gravitation : Float_Type;

package Pendulums is

  type Pendulum is private;
  function New_Pendulum (Length : Float_Type;
                         Theta0 : Float_Type) return Pendulum;
  function Get_X (From : Pendulum) return Float_Type;
  function Get_Y (From : Pendulum) return Float_Type;
  procedure Update_Pendulum (Item : in out Pendulum; Time : in Duration);

private

  type Pendulum is record
     Length   : Float_Type;
     Theta    : Float_Type;
     X        : Float_Type;
     Y        : Float_Type;
     Velocity : Float_Type;
  end record;

end Pendulums;</lang>

example main.adb: <lang Ada>with Ada.Text_IO; with Ada.Calendar; with Pendulums;

procedure Main is

  package Float_Pendulum is new Pendulums (Float, -9.81);
  use Float_Pendulum;
  use type Ada.Calendar.Time;
  My_Pendulum : Pendulum := New_Pendulum (10.0, 30.0);
  Now, Before : Ada.Calendar.Time;

begin

  Before := Ada.Calendar.Clock;
  loop
     Now := Ada.Calendar.Clock;
     Update_Pendulum (My_Pendulum, Now - Before);
     -- output positions relative to origin
     -- replace with graphical output if wanted
     Ada.Text_IO.Put_Line (" X: " & Float'Image (Get_X (My_Pendulum)) &
                           " Y: " & Float'Image (Get_Y (My_Pendulum)));
     Before := Now;
     Delay 0.1;
  end loop;

end Main;</lang>

C

Library: SDL
Library: Cairo

<lang c>#include <stdio.h>

  1. include <stdlib.h>
  2. include <stdbool.h>
  3. include <math.h>
  4. include <cairo.h>
  5. include <SDL.h>
  1. define WIDTH 320.0
  2. define HEIGHT 240.0
  1. ifndef M_PI
  2. define M_PI 3.14159265358979323846
  3. endif

const double L = HEIGHT*150.0/200.0; // length double theta = M_PI/4.0; double dtheta = 0.0; int homeX = WIDTH / 2; int homeY = HEIGHT*25.0/200.0;

static uint32_t update_anim(uint32_t iv, void *p) {

 double scaling = 30.0 / (L*L); // this seems better for this impl.
 
 double firstDDtheta = -sin(theta) * scaling;
 double midDtheta = dtheta + firstDDtheta;
 double midtheta = theta + (dtheta + midDtheta)/2.0;
 double midDDtheta = -sin(midtheta) * scaling;
 midDtheta = dtheta + (firstDDtheta + midDDtheta)/2.0;
 midtheta = theta + (dtheta + midDtheta)/2.0;
 midDDtheta = -sin(midtheta) * scaling;
 double lastDtheta = midDtheta + midDDtheta;
 double lasttheta = midtheta + (midDtheta + lastDtheta)/2.0;
 double lastDDtheta = -sin(lasttheta) * scaling;
 lastDtheta = midDtheta + (midDDtheta + lastDDtheta)/2.0;
 lasttheta = midtheta + (midDtheta + lastDtheta)/2.0;
 dtheta = lastDtheta;
 theta = lasttheta;
 return iv;

}

void show_pendulum(SDL_Surface *s) {

 cairo_surface_t *surf = 
   cairo_image_surface_create_for_data( s->pixels,

CAIRO_FORMAT_RGB24, s->w, s->h, s->pitch);

 cairo_t *ct = cairo_create(surf);
 double x = homeX + L*sin(theta);
 double y = homeY + L*cos(theta);
 cairo_set_source_rgba(ct, 0, 0, 0, 1);
 cairo_new_path(ct);
 cairo_rectangle(ct, 0, 0, WIDTH, HEIGHT);
 cairo_fill(ct);
 cairo_set_source_rgba(ct, 1, 1, 1, 1);
 
 // rod
 cairo_new_path(ct);
 cairo_move_to(ct, homeX, homeY);
 cairo_line_to(ct, x, y);
 cairo_stroke(ct);
 cairo_set_source_rgba(ct, 1, 1, 0, 1);
 
 // bob
 cairo_new_path(ct);
 cairo_move_to(ct, x, y);
 cairo_arc(ct, x, y, 20.0, 0, 2*M_PI);
 cairo_fill(ct);
 cairo_set_source_rgba(ct, 0.7, 0.7, 0.7, 1);
 // plate
 cairo_new_path(ct);
 cairo_move_to(ct, 0, homeY);
 cairo_line_to(ct, WIDTH, homeY);
 cairo_stroke(ct);
 // pivot
 cairo_new_path(ct);
 cairo_move_to(ct, homeX, homeY);
 cairo_arc(ct, homeX, homeY, 5.0, 0, 2*M_PI);
 cairo_fill(ct);
 cairo_surface_destroy(surf);
 cairo_destroy(ct);

}

int main() {

 SDL_Surface *scr, *t;
 SDL_Event event[1];
 bool quit = false;
 if ( SDL_Init(SDL_INIT_VIDEO | SDL_INIT_TIMER) >= 0 ) {
   atexit(SDL_Quit);
   if ( SDL_SetVideoMode(WIDTH, HEIGHT, 32, SDL_DOUBLEBUF) != NULL ) {
     
     scr = SDL_GetVideoSurface();
     SDL_AddTimer(30, update_anim, NULL);
     event->type = SDL_VIDEOEXPOSE;
     SDL_PushEvent(event);
     while(SDL_WaitEvent(event) && !quit) {

switch(event->type) { case SDL_VIDEOEXPOSE: while(SDL_LockSurface(scr) != 0) SDL_Delay(1);

show_pendulum(scr);

SDL_UnlockSurface(scr); SDL_Flip(scr); event->type = SDL_VIDEOEXPOSE; SDL_PushEvent(event); break; case SDL_KEYDOWN: if (event->key.keysym.sym == SDLK_q) quit = true; break; }

     }
     SDL_FreeSurface(scr);
   }
 }
 return 0;

}</lang>

C#

Library: Windows Forms

<lang csharp> using System; using System.Drawing; using System.Windows.Forms;

class CSharpPendulum {

   Form _form;
   Timer _timer;
   
   double _angle = Math.PI / 2, 
          _angleAccel, 
          _angleVelocity = 0, 
          _dt = 0.1;
   
   int _length = 50;
   [STAThread]
   static void Main()
   {
       var p = new CSharpPendulum();
   }
   public CSharpPendulum()
   {
       _form = new Form() { Text = "Pendulum", Width = 200, Height = 200 };
       _timer = new Timer() { Interval = 30 };
       _timer.Tick += delegate(object sender, EventArgs e)
       {
           int anchorX = (_form.Width / 2) - 12,
               anchorY = _form.Height / 4,
               ballX = anchorX + (int)(Math.Sin(_angle) * _length),
               ballY = anchorY + (int)(Math.Cos(_angle) * _length);
           _angleAccel = -9.81 / _length * Math.Sin(_angle);
           _angleVelocity += _angleAccel * _dt;
           _angle += _angleVelocity * _dt;
         
           Bitmap dblBuffer = new Bitmap(_form.Width, _form.Height);
           Graphics g = Graphics.FromImage(dblBuffer);
           Graphics f = Graphics.FromHwnd(_form.Handle);
           g.DrawLine(Pens.Black, new Point(anchorX, anchorY), new Point(ballX, ballY));
           g.FillEllipse(Brushes.Black, anchorX - 3, anchorY - 4, 7, 7);
           g.FillEllipse(Brushes.DarkGoldenrod, ballX - 7, ballY - 7, 14, 14);
           
           f.Clear(Color.White);
           f.DrawImage(dblBuffer, new Point(0, 0));    
       };
       _timer.Start();
       Application.Run(_form);
   }     

} </lang>

Clojure

Clojure solution using an atom and a separate rendering thread

Library: Swing
Library: AWT

<lang clojure> (ns pendulum

 (:import
   (javax.swing JFrame)
   (java.awt Canvas Graphics Color)))

(def length 200) (def width (* 2 (+ 50 length))) (def height (* 3 (/ length 2))) (def dt 0.1) (def g 9.812) (def k (- (/ g length))) (def anchor-x (/ width 2)) (def anchor-y (/ height 8)) (def angle (atom (/ (Math/PI) 2)))

(defn draw [#^Canvas canvas angle]

 (let [buffer  (.getBufferStrategy canvas)
       g       (.getDrawGraphics buffer)
       ball-x (+ anchor-x (* (Math/sin angle) length))
       ball-y (+ anchor-y (* (Math/cos angle) length))]
   (try      
     (doto g
       (.setColor Color/BLACK)
       (.fillRect 0 0 width height)
       (.setColor Color/RED)
       (.drawLine anchor-x anchor-y ball-x ball-y)
       (.setColor Color/YELLOW)
       (.fillOval (- anchor-x 3) (- anchor-y 4) 7 7)
       (.fillOval (- ball-x 7) (- ball-y 7) 14 14))      
     (finally (.dispose g)))
   (if-not (.contentsLost buffer)
     (.show buffer)) ))

(defn start-renderer [canvas]

 (->>
   (fn [] (draw canvas @angle) (recur))
   (new Thread)
   (.start)))

(defn -main [& args]

 (let [frame  (JFrame. "Pendulum")
       canvas (Canvas.)]

   (doto frame
     (.setSize width height)      
     (.setDefaultCloseOperation JFrame/EXIT_ON_CLOSE)
     (.setResizable false)
     (.add canvas)
     (.setVisible true))

   (doto canvas
     (.createBufferStrategy 2)      
     (.setVisible true)
     (.requestFocus))

   (start-renderer canvas)

   (loop [v 0]      
     (swap! angle #(+ % (* v dt)))
     (Thread/sleep 15)
     (recur (+ v (* k (Math/sin @angle) dt)))) ))

(-main) </lang>

E

Works with: E-on-Java

(Uses Java Swing for GUI. The animation logic is independent, however.)

The angle of a pendulum with length and acceleration due to gravity with all its mass at the end and no friction/air resistance has an acceleration at any given moment of

This simulation uses this formula directly, updating the velocity from the acceleration and the position from the velocity; inaccuracy results from the finite timestep.

The event flow works like this: The clock object created by the simulation steps the simulation on the specified in the interval. The simulation writes its output to angle, which is a Lamport slot which can notify of updates. The whenever set up by makeDisplayComponent listens for updates and triggers redrawing as long as interest has been expressed, which is done whenever the component actually redraws, which happens only if the component's window is still on screen. When the window is closed, additionally, the simulation itself is stopped and the application allowed to exit. (This logic is more general than necessary; it is designed to be suitable for a larger application as well.)

<lang e>#!/usr/bin/env rune pragma.syntax("0.9")

def pi := (-1.0).acos() def makeEPainter := <unsafe:com.zooko.tray.makeEPainter> def makeLamportSlot := <import:org.erights.e.elib.slot.makeLamportSlot> def whenever := <import:org.erights.e.elib.slot.whenever> def colors := <import:java.awt.makeColor>

  1. --------------------------------------------------------------
  2. --- Definitions

def makePendulumSim(length_m :float64,

                   gravity_mps2 :float64,
                   initialAngle_rad :float64,
                   timestep_ms :int) {
 var velocity := 0
 def &angle := makeLamportSlot(initialAngle_rad)
 def k := -gravity_mps2/length_m
 def timestep_s := timestep_ms / 1000
 def clock := timer.every(timestep_ms, fn _ {
   def acceleration := k * angle.sin()
   velocity += acceleration * timestep_s
   angle    += velocity     * timestep_s
 })
 return [clock, &angle]

}

def makeDisplayComponent(&angle) {

 def c
 def updater := whenever([&angle], fn { c.repaint() })
 
 bind c := makeEPainter(def paintCallback {
   to paintComponent(g) {
     try {
       def originX := c.getWidth() // 2
       def originY := c.getHeight() // 2
       def pendRadius := (originX.min(originY) * 0.95).round()
       def ballRadius := (originX.min(originY) * 0.04).round()
       def ballX := (originX + angle.sin() * pendRadius).round()
       def ballY := (originY + angle.cos() * pendRadius).round()
       g.setColor(colors.getWhite())
       g.fillRect(0, 0, c.getWidth(), c.getHeight())
       g.setColor(colors.getBlack())
       
       g.fillOval(originX - 2, originY - 2, 4, 4)
       g.drawLine(originX, originY, ballX, ballY)
       g.fillOval(ballX - ballRadius, ballY - ballRadius, ballRadius * 2, ballRadius * 2)
     
       updater[] # provoke interest provided that we did get drawn (window not closed)
     } catch p {
       stderr.println(`In paint callback: $p${p.eStack()}`)
     }
   }
 })
 
 c.setPreferredSize(<awt:makeDimension>(300, 300))
 return c

}

  1. --------------------------------------------------------------
  2. --- Application setup

def [clock, &angle] := makePendulumSim(1, 9.80665, pi*99/100, 10)

  1. Initialize AWT, move to AWT event thread

when (currentVat.morphInto("awt")) -> {

 # Create the window
 def frame := <unsafe:javax.swing.makeJFrame>("Pendulum")
 frame.setContentPane(def display := makeDisplayComponent(&angle))
 frame.addWindowListener(def mainWindowListener {
   to windowClosing(_) {
     clock.stop()
     interp.continueAtTop()
   }
   match _ {}
 })
 frame.setLocation(50, 50)
 frame.pack()
 # Start and become visible
 frame.show()
 clock.start()

}

interp.blockAtTop()</lang>

Factor

Approximation of the pendulum for small swings : theta = theta0 * cos(omega0 * t) <lang factor>USING: accessors alarms arrays calendar colors.constants kernel locals math math.constants math.functions math.rectangles math.vectors opengl sequences system ui ui.gadgets ui.render ; IN: pendulum

CONSTANT: g 9.81 CONSTANT: l 20 CONSTANT: theta0 0.5

current-time ( -- time ) nano-count -9 10^ * ;
T0 ( -- T0 ) 2 pi l g / sqrt * * ;
omega0 ( -- omega0 ) 2 pi * T0 / ;
theta ( -- theta ) current-time omega0 * cos theta0 * ;
relative-xy ( theta l -- xy )
   [ [ sin ] [ cos ] bi ]
   [ [ * ] curry ] bi* bi@ 2array ;
theta-to-xy ( origin theta l -- xy ) relative-xy v+ ;

TUPLE: pendulum-gadget < gadget alarm ;

O ( gadget -- origin ) rect-bounds [ drop ] [ first 2 / ] bi* 0 2array ;
window-l ( gadget -- l ) rect-bounds [ drop ] [ second ] bi* ;
gadget-xy ( gadget -- xy ) [ O ] [ drop theta ] [ window-l ] tri theta-to-xy ;

M: pendulum-gadget draw-gadget*

   COLOR: black gl-color
   [ O ] [ gadget-xy ] bi gl-line ;

M:: pendulum-gadget graft* ( gadget -- )

   [ gadget relayout-1 ]
   20 milliseconds every gadget (>>alarm) ;

M: pendulum-gadget ungraft* alarm>> cancel-alarm ;

<pendulum-gadget> ( -- gadget )
   pendulum-gadget new 
   { 500 500 } >>pref-dim ;
pendulum-main ( -- )
   [ <pendulum-gadget> "pendulum" open-window ] with-ui ;

MAIN: pendulum-main </lang>

F#

A nice application of F#'s support for units of measure. <lang fsharp>open System open System.Drawing open System.Windows.Forms

// define units of measurement [<Measure>] type m; // metres [<Measure>] type s; // seconds

// a pendulum is represented as a record of physical quantities type Pendulum =

{ length   : float<m>
  gravity  : float<m/s^2>
  velocity : float<m/s>
  angle    : float
}

// calculate the next state of a pendulum let next pendulum deltaT : Pendulum =

 let k = -pendulum.gravity / pendulum.length
 let acceleration = k * Math.Sin pendulum.angle * 1.0<m> 
 let newVelocity = pendulum.velocity + acceleration * deltaT
 let newAngle = pendulum.angle + newVelocity * deltaT / 1.0<m>
 { pendulum with velocity = newVelocity; angle = newAngle }

// paint a pendulum (using hard-coded screen coordinates) let paint pendulum (gr: System.Drawing.Graphics) =

 let homeX = 160
 let homeY = 50
 let length = 140.0
 // draw plate
 gr.DrawLine( new Pen(Brushes.Gray, width=2.0f), 0, homeY, 320, homeY )
 // draw pivot
 gr.FillEllipse( Brushes.Gray,           homeX-5, homeY-5, 10, 10 )
 gr.DrawEllipse( new Pen(Brushes.Black), homeX-5, homeY-5, 10, 10 )
 // draw the pendulum itself
 let x = homeX + int( length * Math.Sin pendulum.angle )
 let y = homeY + int( length * Math.Cos pendulum.angle )
 // draw rod
 gr.DrawLine( new Pen(Brushes.Black, width=3.0f), homeX, homeY, x, y )
 // draw bob
 gr.FillEllipse( Brushes.Yellow,         x-15, y-15, 30, 30 )
 gr.DrawEllipse( new Pen(Brushes.Black), x-15, y-15, 30, 30 )

// defines an operator "-?" that calculates the time from t2 to t1 // where t2 is optional let (-?) (t1: DateTime) (t2: DateTime option) : float =

 match t2 with
 | None   -> 0.0 // only one timepoint given -> difference is 0
 | Some t -> (t1 - t).TotalSeconds * 1.0

// our main window is double-buffered form that reacts to paint events type PendulumForm() as self =

 inherit Form(Width=325, Height=240, Text="Pendulum")
 let mutable pendulum = { length   = 1.0<m>;
                          gravity  = 9.81<m/s^2>
                          velocity = 0.0<m/s>
                          angle    = Math.PI / 2.0
                        }
 let mutable lastPaintedAt = None
 let updateFreq = 0.01
 do self.DoubleBuffered <- true
    self.Paint.Add( fun args ->
      let now = DateTime.Now
      let deltaT = now -? lastPaintedAt |> min 0.01 
      lastPaintedAt <- Some now
      pendulum <- next pendulum deltaT
      let gr = args.Graphics
      gr.Clear( Color.LightGray )
      paint pendulum gr
      // initiate a new paint event after a while (non-blocking)
      async { do! Async.Sleep( int( 1000.0 * updateFreq / 1.0 ) )
              self.Invalidate()
           }
      |> Async.Start 
    )

[<STAThread>] Application.Run( new PendulumForm( Visible=true ) )</lang>

Haskell

Using

Library: HGL

from HackageDB

<lang haskell>import Graphics.HGL.Draw.Monad (Graphic, ) import Graphics.HGL.Draw.Picture import Graphics.HGL.Utils import Graphics.HGL.Window import Graphics.HGL.Run

import Control.Exception (bracket, ) import Control.Arrow

toInt = fromIntegral.round

pendulum = runGraphics $

 bracket
   (openWindowEx "Pendulum animation task" Nothing (600,400) DoubleBuffered (Just 30))
   closeWindow
   (\w -> mapM_ ((\ g -> setGraphic w g >> getWindowTick w).

(\ (x, y) -> overGraphic (line (300, 0) (x, y)) (ellipse (x - 12, y + 12) (x + 12, y - 12)) )) pts)

where
   dt = 1/30
   t = - pi/4
   l = 1
   g = 9.812 
   nextAVT (a,v,t) = (a', v', t + v' * dt) where

a' = - (g / l) * sin t v' = v + a' * dt

   pts = map (\(_,t,_) -> (toInt.(300+).(300*).cos &&& toInt. (300*).sin) (pi/2+0.6*t) )

$ iterate nextAVT (- (g / l) * sin t, t, 0)</lang> Use (interpreter ghci):

*Main> pendulum

HicEst

DIFFEQ and the callback procedure pendulum numerically integrate the pendulum equation. The display window can be resized during the run, but for window width not equal to 2*height the pendulum rod becomes a rubber band instead: <lang HicEst>REAL :: msec=10, Lrod=1, dBob=0.03, g=9.81, Theta(2), dTheta(2) BobMargins = ALIAS(ls, rs, ts, bs) ! box margins to draw the bob


Theta = (1, 0)  ! initial angle and velocity start_t = TIME()

DO i = 1, 1E100  ! "forever"

  end_t = TIME()     ! to integrate in real-time sections:
  DIFFEQ(Callback="pendulum", T=end_t, Y=Theta, DY=dTheta, T0=start_t)
  xBob = (SIN(Theta(1)) + 1) / 2
  yBob = COS(Theta(1)) - dBob
  ! create or clear window and draw pendulum bob at (xBob, yBob):
  WINDOW(WIN=wh, LeftSpace=0, RightSpace=0, TopSpace=0, BottomSpace=0, Up=999)
  BobMargins = (xBob-dBob, 1-xBob-dBob, yBob-dBob, 1-yBob-dBob)
  WINDOW(WIN=wh, LeftSpace=ls, RightSpace=rs, TopSpace=ts, BottomSpace=bs)
  WRITE(WIN=wh, DeCoRation='EL=4, BC=4') ! flooded red ellipse as bob
  ! draw the rod hanging from the center of the window:
  WINDOW(WIN=wh, LeftSpace=0.5, TopSpace=0, RightSpace=rs+dBob)
  WRITE(WIN=wh, DeCoRation='LI=0 0; 1 1, FC=4.02') ! red pendulum rod
  SYSTEM(WAIT=msec)
  start_t = end_t

ENDDO

END

SUBROUTINE pendulum  ! Theta" = - (g/Lrod) * SIN(Theta)

 dTheta(1) = Theta(2)              ! Theta' = Theta(2)  substitution
 dTheta(2) = -g/Lrod*SIN(Theta(1)) ! Theta" = Theta(2)' = -g/Lrod*SIN(Theta(1))

END</lang>

J

<lang j>require 'gl2 trig' coinsert 'jgl2'

DT =: %30 NB. seconds ANGLE=: 0.25p1 NB. radians L =: 1 NB. metres G =: 9.80665 NB. ms_2 VEL =: 0 NB. ms_1

PEND=: noun define pc pend;pn "Pendulum"; xywh 0 0 320 200;cc isi isigraph rightmove bottommove; pas 0 0;pcenter; rem form end; )

pend_run =: verb def ' wd PEND,;pshow;timer ,":DT * 1000 ' pend_close =: verb def ' wd timer 0; pclose ' pend_isi_paint=: verb def ' drawPendulum ANGLE '

sys_timer_z_=: verb define

 recalcAngle 
 wd 'psel pend; setinvalid isi'

)

recalcAngle=: verb define

 accel=. - (G % L) * sin ANGLE
 VEL  =: VEL + accel * DT
 ANGLE=: ANGLE + VEL * DT

)

drawPendulum=: verb define

 width=. {. glqwh
 ps=. (-: width) , 40
 pe=. ps + 280 <.@* (cos , sin) 0.5p1 + y    NB. adjust orientation
 glrgb 91 91 91
 glbrush
 gllines ps , pe
 glellipse (,~ ps - -:) 40 15
 glellipse (,~ pe - -:) 20 20
 glrect 0 0 ,width, 40

)

pend_run NB. run animation</lang>

Java

Library: Swing
Library: AWT

<lang java>import java.awt.*; import javax.swing.*;

public class Pendulum extends JPanel implements Runnable {

   private double angle = Math.PI / 2;
   private int length;
   public Pendulum(int length) {
       this.length = length;
       setDoubleBuffered(true);
   }
   @Override
   public void paint(Graphics g) {
       g.setColor(Color.WHITE);
       g.fillRect(0, 0, getWidth(), getHeight());
       g.setColor(Color.BLACK);
       int anchorX = getWidth() / 2, anchorY = getHeight() / 4;
       int ballX = anchorX + (int) (Math.sin(angle) * length);
       int ballY = anchorY + (int) (Math.cos(angle) * length);
       g.drawLine(anchorX, anchorY, ballX, ballY);
       g.fillOval(anchorX - 3, anchorY - 4, 7, 7);
       g.fillOval(ballX - 7, ballY - 7, 14, 14);
   }
   public void run() {
       double angleAccel, angleVelocity = 0, dt = 0.1;
       while (true) {
           angleAccel = -9.81 / length * Math.sin(angle);
           angleVelocity += angleAccel * dt;
           angle += angleVelocity * dt;
           repaint();
           try { Thread.sleep(15); } catch (InterruptedException ex) {}
       }
   }
   @Override
   public Dimension getPreferredSize() {
       return new Dimension(2 * length + 50, length / 2 * 3);
   }
   public static void main(String[] args) {
       JFrame f = new JFrame("Pendulum");
       Pendulum p = new Pendulum(200);
       f.add(p);
       f.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
       f.pack();
       f.setVisible(true);
       new Thread(p).start();
   }

}</lang>

JavaScript + <canvas>

Translation of: E

(plus gratuitous motion blur)

<lang javascript><html><head>

 <title>Pendulum</title>

</head><body style="background: gray;">

<canvas id="canvas" width="600" height="600">

Sorry, your browser does not support the <canvas> used to display the pendulum animation.

</canvas> <script>

 function PendulumSim(length_m, gravity_mps2, initialAngle_rad, timestep_ms, callback) {
   var velocity = 0;
   var angle = initialAngle_rad;
   var k = -gravity_mps2/length_m;
   var timestep_s = timestep_ms / 1000;
   return setInterval(function () {
     var acceleration = k * Math.sin(angle);
     velocity += acceleration * timestep_s;
     angle    += velocity     * timestep_s;
     callback(angle);
   }, timestep_ms);
 }
 
 var canvas = document.getElementById('canvas');
 var context = canvas.getContext('2d');
 var prev=0;
 var sim = PendulumSim(1, 9.80665, Math.PI*99/100, 10, function (angle) {
   var rPend = Math.min(canvas.width, canvas.height) * 0.47;
   var rBall = Math.min(canvas.width, canvas.height) * 0.02;
   var rBar = Math.min(canvas.width, canvas.height) * 0.005;
   var ballX = Math.sin(angle) * rPend;
   var ballY = Math.cos(angle) * rPend;
   context.fillStyle = "rgba(255,255,255,0.51)";
   context.globalCompositeOperation = "destination-out";
   context.fillRect(0, 0, canvas.width, canvas.height);
   
   context.fillStyle = "yellow";
   context.strokeStyle = "rgba(0,0,0,"+Math.max(0,1-Math.abs(prev-angle)*10)+")";
   context.globalCompositeOperation = "source-over";
   context.save();
     context.translate(canvas.width/2, canvas.height/2);
     context.rotate(angle);
     
     context.beginPath();
     context.rect(-rBar, -rBar, rBar*2, rPend+rBar*2);
     context.fill();
     context.stroke();
     
     context.beginPath();
     context.arc(0, rPend, rBall, 0, Math.PI*2, false);
     context.fill();
     context.stroke();
   context.restore();
   prev=angle;
 });

</script>

</body></html></lang>

Liberty BASIC

<lang lb>nomainwin

   WindowWidth = 400
   WindowHeight = 300
   open "Pendulum" for graphics_nsb_nf as #main
   #main, "down;fill white; flush"
   #main, "color black"
   #main, "trapclose [quit.main]"
   Angle = asn(1)
   DeltaT = 0.1
   PendLength = 150
   FixX = int(WindowWidth / 2)
   FixY = 40
   timer 30, [swing]
   wait

[swing]

   #main, "cls"
   #main, "discard"
   PlumbobX = FixX + int(sin(Angle) * PendLength)
   PlumbobY = FixY + int(cos(Angle) * PendLength)
   AngAccel = -9.81 / PendLength * sin(Angle)
   AngVelocity = AngVelocity + AngAccel * DeltaT
   Angle = Angle + AngVelocity * DeltaT
   #main, "backcolor black"
   #main, "place ";FixX;" ";FixY
   #main, "circlefilled 3"
   #main, "line ";FixX;" ";FixY;" ";PlumbobX;" ";PlumbobY
   #main, "backcolor red"
   #main, "circlefilled 10"
   wait

[quit.main]

   close #main
   end</lang>

Works with: UCB Logo

<lang logo>make "angle 45 make "L 1 make "bob 10

to draw.pendulum

 clearscreen
 seth :angle+180		; down on screen is 180
 forward :L*100-:bob
 penup
 forward :bob
 pendown
 arc 360 :bob

end

make "G 9.80665 make "dt 1/30 make "acc 0 make "vel 0

to step.pendulum

 make "acc  -:G / :L * sin :angle
 make "vel   :vel   + :acc * :dt
 make "angle :angle + :vel * :dt
 wait :dt*60
 draw.pendulum

end

hideturtle until [key?] [step.pendulum]</lang>

MATLAB

pendulum.m <lang MATLAB>%This is a numerical simulation of a pendulum with a massless pivot arm.

%% User Defined Parameters %Define external parameters g = -9.8; deltaTime = 1/50; %Decreasing this will increase simulation accuracy endTime = 16;

%Define pendulum rodPivotPoint = [2 2]; %rectangular coordinates rodLength = 1; mass = 1; %of the bob radius = .2; %of the bob theta = 45; %degrees, defines initial position of the bob velocity = [0 0]; %cylindrical coordinates; first entry is radial velocity,

                 %second entry is angular velocity

%% Simulation assert(radius < rodLength,'Pendulum bob radius must be less than the length of the rod.');

position = rodPivotPoint - (rodLength*[-sind(theta) cosd(theta)]); %in rectangular coordinates

%Generate graphics, render pendulum figure; axesHandle = gca; xlim(axesHandle, [(rodPivotPoint(1) - rodLength - radius) (rodPivotPoint(1) + rodLength + radius)] ); ylim(axesHandle, [(rodPivotPoint(2) - rodLength - radius) (rodPivotPoint(2) + rodLength + radius)] );

rectHandle = rectangle('Position',[(position - radius/2) radius radius],...

   'Curvature',[1,1],'FaceColor','g'); %Pendulum bob

hold on plot(rodPivotPoint(1),rodPivotPoint(2),'^'); %pendulum pivot lineHandle = line([rodPivotPoint(1) position(1)],...

   [rodPivotPoint(2) position(2)]); %pendulum rod

hold off

%Run simulation, all calculations are performed in cylindrical coordinates for time = (deltaTime:deltaTime:endTime)

   drawnow; %Forces MATLAB to render the pendulum
   
   %Find total force
   gravitationalForceCylindrical = [mass*g*cosd(theta) mass*g*sind(theta)];
   
   %This code is just incase you want to add more forces,e.g friction
   totalForce = gravitationalForceCylindrical; 
   
   %If the rod isn't massless or is a spring, etc., modify this line
   %accordingly
   rodForce = [-totalForce(1) 0]; %cylindrical coordinates
   
   totalForce = totalForce + rodForce;
   
   acceleration = totalForce / mass; %F = ma
   velocity = velocity + acceleration * deltaTime;
   rodLength = rodLength + velocity(1) * deltaTime;
   theta = theta + velocity(2) * deltaTime;
   
   position = rodPivotPoint - (rodLength*[-sind(theta) cosd(theta)]);
   
   %Update figure with new position info
   set(rectHandle,'Position',[(position - radius/2) radius radius]);
   set(lineHandle,'XData',[rodPivotPoint(1) position(1)],'YData',...
       [rodPivotPoint(2) position(2)]);

end</lang>

Oz

Inspired by the E and Ruby versions.

<lang oz>declare

 [QTk] = {Link ['x-oz://system/wp/QTk.ozf']}
 Pi = 3.14159265
 class PendulumModel
    feat

K

    attr

angle velocity

    meth init(length:L       <= 1.0    %% meters

gravity:G <= 9.81  %% m/s² initialAngle:A <= Pi/2.) %% radians self.K = ~G / L angle := A velocity := 0.0

    end
    meth nextAngle(deltaT:DeltaTMS %% milliseconds

?Angle)  %% radians

       DeltaT = {Int.toFloat DeltaTMS} / 1000.0 %% seconds
       Acceleration = self.K * {Sin @angle}
    in
       velocity := @velocity + Acceleration * DeltaT
       angle := @angle + @velocity * DeltaT
       Angle = @angle
    end
 end
 %% Animates a pendulum on a given canvas.
 class PendulumAnimation from Time.repeat
    feat
       Pend
       Rod
       Bob
       home:pos(x:160 y:50)
       length:140.0

delay

    meth init(Pendulum Canvas delay:Delay <= 25) %% milliseconds

self.Pend = Pendulum self.delay = Delay %% plate and pivot

       {Canvas create(line 0 self.home.y 320 self.home.y width:2 fill:grey50)}
       {Canvas create(oval 155 self.home.y-5 165 self.home.y+5 fill:grey50 outline:black)}

%% the pendulum itself self.Rod = {Canvas create(line 1 1 1 1 width:3 fill:black handle:$)}

       self.Bob = {Canvas create(oval 1 1 2 2 fill:yellow outline:black handle:$)}
       %%
       {self setRepAll(action:Animate delay:Delay)}
    end
    meth Animate

Theta = {self.Pend nextAngle(deltaT:self.delay $)} %% calculate x and y from angle X = self.home.x + {Float.toInt self.length * {Sin Theta}} Y = self.home.y + {Float.toInt self.length * {Cos Theta}}

    in

%% update canvas try {self.Rod setCoords(self.home.x self.home.y X Y)} {self.Bob setCoords(X-15 Y-15 X+15 Y+15)} catch system(tk(alreadyClosed ...) ...) then skip end

    end
 end
  
 Pendulum = {New PendulumModel init}
 Canvas
 GUI = td(title:"Pendulum"
          canvas(width:320 height:210 handle:?Canvas)
          action:proc {$} {Animation stop} {Window close} end
         )
 Window = {QTk.build GUI}
 Animation = {New PendulumAnimation init(Pendulum Canvas)}

in

 {Window show}
 {Animation go}</lang>

Prolog

SWI-Prolog has a graphic interface XPCE. <lang Prolog>:- use_module(library(pce)).

pendulum :- new(D, window('Pendulum')), send(D, size, size(560, 300)), new(Line, line(80, 50, 480, 50)), send(D, display, Line), new(Circle, circle(20)), send(Circle, fill_pattern, colour(@default, 0, 0, 0)), new(Boule, circle(60)), send(Boule, fill_pattern, colour(@default, 0, 0, 0)), send(D, display, Circle, point(270,40)), send(Circle, handle, handle(h/2, w/2, in)), send(Boule, handle, handle(h/2, w/2, out)), send(Circle, connect, Boule, link(in, out, line(0,0,0,0,none))), new(Anim, animation(D, 0.0, Boule, 200.0)), send(D, done_message, and(message(Anim, free), message(Boule, free), message(Circle, free), message(@receiver,destroy))), send(Anim?mytimer, start), send(D, open).



- pce_begin_class(animation(window, angle, boule, len_pendulum), object).

variable(window, object, both, "Display window"). variable(boule, object, both, "bowl of the pendulum"). variable(len_pendulum, object, both, "len of the pendulum"). variable(angle, object, both, "angle with the horizontal"). variable(delta, object, both, "increment of the angle"). variable(mytimer, timer, both, "timer of the animation").

initialise(P, W:object, A:object, B : object, L:object) :->

       "Creation of the object"::
       send(P, window, W),
       send(P, angle, A),
       send(P, boule, B),
       send(P, len_pendulum, L),
       send(P, delta, 0.01),

send(P, mytimer, new(_, timer(0.01,message(P, anim_message)))).

% method called when the object is destroyed % first the timer is stopped % then all the resources are freed unlink(P) :-> send(P?mytimer, stop), send(P, send_super, unlink).


% message processed by the timer anim_message(P) :-> get(P, angle, A), get(P, len_pendulum, L), calc(A, L, X, Y), get(P, window, W), get(P, boule, B), send(W, display, B, point(X,Y)), % computation of the next position get(P, delta, D), next_Angle(A, D, NA, ND), send(P, angle, NA), send(P, delta, ND).

- pce_end_class.

% computation of the position of the bowl. calc(Ang, Len, X, Y) :- X is Len * cos(Ang)+ 250, Y is Len * sin(Ang) + 20.


% computation of the next angle % if we reach 0 or pi, delta change. next_Angle(A, D, NA, ND) :- NA is D + A, (((D > 0, abs(pi-NA) < 0.01); (D < 0, abs(NA) < 0.01))-> ND = - D; ND = D). </lang>

PureBasic

If the code was part of a larger application it could be improved by specifying constants for the locations of image elements. <lang PureBasic>Procedure handleError(x, msg.s)

 If Not x
   MessageRequester("Error", msg)
   End
 EndIf

EndProcedure

  1. ScreenW = 320
  2. ScreenH = 210

handleError(OpenWindow(0, 0, 0, #ScreenW, #ScreenH, "Animated Pendulum", #PB_Window_SystemMenu), "Can't open window.") handleError(InitSprite(), "Can't setup sprite display.") handleError(OpenWindowedScreen(WindowID(0), 0, 0, #ScreenW, #ScreenH, 0, 0, 0), "Can't open screen.")

Enumeration ;sprites

 #bob_spr
 #ceiling_spr
 #pivot_spr

EndEnumeration

TransparentSpriteColor(#PB_Default, RGB(255, 0, 255)) CreateSprite(#bob_spr, 32, 32) StartDrawing(SpriteOutput(#bob_spr))

 Box(0, 0, 32, 32, RGB(255, 0, 255))
 Circle(16, 16, 15, RGB(253, 252, 3))
 DrawingMode(#PB_2DDrawing_Outlined)
 Circle(16, 16, 15, RGB(0, 0, 0))

StopDrawing()

CreateSprite(#pivot_spr, 10, 10) StartDrawing(SpriteOutput(#pivot_spr))

 Box(0, 0, 10, 10, RGB(255, 0, 255))
 Circle(5, 5, 4, RGB(125, 125, 125))
 DrawingMode(#PB_2DDrawing_Outlined)
 Circle(5, 5, 4, RGB(0,0 , 0))

StopDrawing()

CreateSprite(#ceiling_spr,#ScreenW,2) StartDrawing(SpriteOutput(#ceiling_spr))

 Box(0,0,SpriteWidth(#ceiling_spr), SpriteHeight(#ceiling_spr), RGB(126, 126, 126))

StopDrawing()

Structure pendulum

 length.d   ; meters
 constant.d ; -g/l
 gravity.d  ; m/s²
 angle.d    ; radians
 velocity.d ; m/s

EndStructure

Procedure initPendulum(*pendulum.pendulum, length.d = 1.0, gravity.d = 9.81, initialAngle.d = #PI / 2)

 With *pendulum
   \length = length
   \gravity = gravity
   \angle = initialAngle
   \constant = -gravity / length
   \velocity = 0.0
 EndWith

EndProcedure


Procedure updatePendulum(*pendulum.pendulum, deltaTime.d)

 deltaTime = deltaTime / 1000.0 ;ms
 Protected acceleration.d = *pendulum\constant * Sin(*pendulum\angle)
 *pendulum\velocity + acceleration * deltaTime
 *pendulum\angle + *pendulum\velocity * deltaTime

EndProcedure

Procedure drawBackground()

 ClearScreen(RGB(190,190,190))
 ;draw ceiling
 DisplaySprite(#ceiling_spr, 0, 47)
 ;draw pivot
 DisplayTransparentSprite(#pivot_spr, 154,43) ;origin in upper-left

EndProcedure

Procedure drawPendulum(*pendulum.pendulum)

 ;draw rod
 Protected x = *pendulum\length * 140 * Sin(*pendulum\angle) ;scale = 1 m/140 pixels
 Protected y = *pendulum\length * 140 * Cos(*pendulum\angle)
 StartDrawing(ScreenOutput())
   LineXY(154 + 5,43 + 5, 154 + 5 + x, 43 + 5 + y) ;draw from pivot-center to bob-center, adjusting for origins
 StopDrawing()
 
 ;draw bob
 DisplayTransparentSprite(#bob_spr, 154 + 5 - 16 + x, 43 + 5 - 16 + y) ;adj for origin in upper-left

EndProcedure

Define pendulum.pendulum, event initPendulum(pendulum) drawPendulum(pendulum)

AddWindowTimer(0, 1, 50) Repeat

 event = WindowEvent()
 Select event
   Case #pb_event_timer
     drawBackground()
     Select EventTimer()
       Case 1
         updatePendulum(pendulum, 50)
         drawPendulum(pendulum)
     EndSelect
     FlipBuffers() 
   Case #PB_Event_CloseWindow
     Break
 EndSelect

ForEver</lang>

Python

Library: pygame
Translation of: C

<lang python>import pygame, sys from pygame.locals import * from math import sin, cos, radians

pygame.init()

WINDOWSIZE = 250 TIMETICK = 100 BOBSIZE = 15

window = pygame.display.set_mode((WINDOWSIZE, WINDOWSIZE)) pygame.display.set_caption("Pendulum")

screen = pygame.display.get_surface() screen.fill((255,255,255))

PIVOT = (WINDOWSIZE/2, WINDOWSIZE/10) SWINGLENGTH = PIVOT[1]*4

class BobMass(pygame.sprite.Sprite):

   def __init__(self):
       pygame.sprite.Sprite.__init__(self)
       self.theta = 45
       self.dtheta = 0
       self.rect = pygame.Rect(PIVOT[0]-SWINGLENGTH*cos(radians(self.theta)),
                               PIVOT[1]+SWINGLENGTH*sin(radians(self.theta)),
                               1,1)
       self.draw()
   def recomputeAngle(self):
       scaling = 3000.0/(SWINGLENGTH**2)
       firstDDtheta = -sin(radians(self.theta))*scaling
       midDtheta = self.dtheta + firstDDtheta
       midtheta = self.theta + (self.dtheta + midDtheta)/2.0
       midDDtheta = -sin(radians(midtheta))*scaling
       midDtheta = self.dtheta + (firstDDtheta + midDDtheta)/2
       midtheta = self.theta + (self.dtheta + midDtheta)/2
       midDDtheta = -sin(radians(midtheta)) * scaling
       lastDtheta = midDtheta + midDDtheta
       lasttheta = midtheta + (midDtheta + lastDtheta)/2.0
       
       lastDDtheta = -sin(radians(lasttheta)) * scaling
       lastDtheta = midDtheta + (midDDtheta + lastDDtheta)/2.0
       lasttheta = midtheta + (midDtheta + lastDtheta)/2.0
       self.dtheta = lastDtheta
       self.theta = lasttheta
       self.rect = pygame.Rect(PIVOT[0]-
                               SWINGLENGTH*sin(radians(self.theta)), 
                               PIVOT[1]+
                               SWINGLENGTH*cos(radians(self.theta)),1,1)


   def draw(self):
       pygame.draw.circle(screen, (0,0,0), PIVOT, 5, 0)
       pygame.draw.circle(screen, (0,0,0), self.rect.center, BOBSIZE, 0)
       pygame.draw.aaline(screen, (0,0,0), PIVOT, self.rect.center)
       pygame.draw.line(screen, (0,0,0), (0, PIVOT[1]), (WINDOWSIZE, PIVOT[1]))
   def update(self):
       self.recomputeAngle()
       screen.fill((255,255,255))
       self.draw()

bob = BobMass()

TICK = USEREVENT + 2 pygame.time.set_timer(TICK, TIMETICK)

def input(events):

   for event in events:
       if event.type == QUIT:
           sys.exit(0)
       elif event.type == TICK:
           bob.update()

while True:

   input(pygame.event.get())
   pygame.display.flip()</lang>

RLaB

The plane pendulum motion is an interesting and easy problem in which the facilities of RLaB for numerical computation and simulation are easily accessible. The parameters of the problem are , the length of the arm, and the magnitude of the gravity.

We start with the mathematical transliteration of the problem. We solve it in plane (2-D) in terms of describing the angle between the -axis and the arm of the pendulum, where the downwards direction is taken as positive. The Newton equation of motian, which is a second-order non-linear ordinary differential equation (ODE) reads

In our example, we will solve the problem as, so called, initial value problem (IVP). That is, we will specify that at the time t=0 the pendulum was at rest , extended at an angle radians (equivalent to 30 degrees).

RLaB has the facilities to solve ODE IVP which are accessible through odeiv solver. This solver requires that the ODE be written as the first order differential equation,

Here, we introduced a vector , for which the original ODE reads

.

The RLaB script that solves the problem is

<lang RLaB> // // example: solve ODE for pendulum //

// we first define the first derivative function for the solver dudt = function(t, u, p) {

 // t-> time
 // u->[theta, dtheta/dt ]
 // p-> g/L, parameter
 rval = zeros(2,1);
 rval[1] =  u[2];
 rval[2] = -p[1] * sin(u[1]);
 return rval;

};

// now we solve the problem // physical parameters L = 5; // (m), the length of the arm of the pendulum p = mks.g / L; // RLaB has a built-in list 'mks' which contains large number of physical constants and conversion factors T0 = 2*const.pi*sqrt(L/mks.g); // approximate period of the pendulum

// initial conditions theta0 = 30; // degrees, initial angle of deflection of pendulum u0 = [theta0*const.pi/180, 0]; // RLaB has a built-in list 'const' of mathematical constants.

// times at which we want solution t = [0:4:1/64] * T0; // solve for 4 approximate periods with at time points spaced at T0/64

// prepare ODEIV solver optsode = <<>>; optsode.eabs = 1e-6; // relative error for step size optsode.erel = 1e-6; // absolute error for step size optsode.delta_t = 1e-6; // maximum dt that code is allowed optsode.stdout = stderr(); // open the text console and in it print the results of each step of calculation optsode.imethod = 5; // use method No. 5 from the odeiv toolkit, Runge-Kutta 8th order Prince-Dormand method //optsode.phase_space = 0; // the solver returns [t, u1(t), u2(t)] which is default behavior optsode.phase_space = 1; // the solver returns [t, u1(t), u2(t), d(u1)/dt(t), d(u2)/dt]

// solver do my bidding y = odeiv(dudt, p, t, u0, optsode);

// Make an animation. We choose to use 'pgplot' rather then 'gnuplot' interface because the former is // faster and thus less cache-demanding, while the latter can be very cache-demanding (it may slow your // linux system quite down if one sends lots of plots for gnuplot to plot). plwins (1); // we will use one pgplot-window

plwin(1); // plot to pgplot-window No. 1; necessary if using more than one pgplot window plimits (-L,L, -1.25*L, 0.25*L); xlabel ("x-coordinate"); ylabel ("z-coordinate"); plegend ("Arm"); for (i in 1:y.nr) {

 // plot a line between the pivot point at (0,0) and the current position of the pendulum
 arm_line = [0,0; L*sin(y[i;2]), -L*cos(y[i;2])]; // this is because theta is between the arm and the z-coordinate
 plot  (arm_line); 
 sleep (0.1); // sleep 0.1 seconds between plots

}

</lang>

Ruby

Library: Ruby/Tk
Translation of: Tcl

This does not have the window resizing handling that Tcl does -- I did not spend enough time in the docs to figure out how to get the new window size out of the configuration event. Of interest when running this pendulum side-by-side with the Tcl one: the Tcl pendulum swings noticibly faster.

<lang ruby>require 'tk'

$root = TkRoot.new("title" => "Pendulum Animation") $canvas = TkCanvas.new($root) do

 width 320
 height 200
 create TkcLine, 0,25,320,25,   'tags' => 'plate', 'width' => 2, 'fill' => 'grey50'
 create TkcOval, 155,20,165,30, 'tags' => 'pivot', 'outline' => "", 'fill' => 'grey50'
 create TkcLine, 1,1,1,1, 'tags' => 'rod', 'width' => 3, 'fill' => 'black'
 create TkcOval, 1,1,2,2, 'tags' => 'bob', 'outline' => 'black', 'fill' => 'yellow'

end $canvas.raise('pivot') $canvas.pack('fill' => 'both', 'expand' => true)

$Theta = 45.0 $dTheta = 0.0 $length = 150 $homeX = 160 $homeY = 25

def show_pendulum

 angle = $Theta * Math::PI / 180
 x = $homeX + $length * Math.sin(angle)
 y = $homeY + $length * Math.cos(angle)
 $canvas.coords('rod', $homeX, $homeY, x, y)
 $canvas.coords('bob', x-15, y-15, x+15, y+15)

end

def recompute_angle

 scaling = 3000.0 / ($length ** 2)
 # first estimate
 firstDDTheta = -Math.sin($Theta * Math::PI / 180) * scaling
 midDTheta = $dTheta + firstDDTheta
 midTheta = $Theta + ($dTheta + midDTheta)/2
 # second estimate
 midDDTheta = -Math.sin(midTheta * Math::PI / 180) * scaling
 midDTheta = $dTheta + (firstDDTheta + midDDTheta)/2
 midTheta = $Theta + ($dTheta + midDTheta)/2
 # again, first
 midDDTheta = -Math.sin(midTheta * Math::PI / 180) * scaling
 lastDTheta = midDTheta + midDDTheta
 lastTheta = midTheta + (midDTheta + lastDTheta)/2
 # again, second
 lastDDTheta = -Math.sin(lastTheta * Math::PI/180) * scaling
 lastDTheta = midDTheta + (midDDTheta + lastDDTheta)/2
 lastTheta = midTheta + (midDTheta + lastDTheta)/2
 # Now put the values back in our globals
 $dTheta  = lastDTheta
 $Theta = lastTheta

end

def animate

 recompute_angle
 show_pendulum
 $after_id = $root.after(15) {animate}

end

show_pendulum $after_id = $root.after(500) {animate}

$canvas.bind('<Destroy>') {$root.after_cancel($after_id)}

Tk.mainloop</lang>

Scheme

Library: Scheme/PsTk
Translation of: Ruby

This is a direct translation of the Ruby/Tk example into Scheme + PS/Tk.

<lang scheme>#!r6rs

R6RS implementation of Pendulum Animation

(import (rnrs)

       (lib pstk main) ; change this for your pstk installation
       )

(define PI 3.14159) (define *conv-radians* (/ PI 180)) (define *theta* 45.0) (define *d-theta* 0.0) (define *length* 150) (define *home-x* 160) (define *home-y* 25)

estimates new angle of pendulum

(define (recompute-angle)

 (define (avg a b) (/ (+ a b) 2))
 (let* ((scaling (/ 3000.0 (* *length* *length*)))
        ; first estimate
        (first-dd-theta (- (* (sin (* *theta* *conv-radians*)) scaling)))
        (mid-d-theta (+ *d-theta* first-dd-theta))
        (mid-theta (+ *theta* (avg *d-theta* mid-d-theta)))
        ; second estimate
        (mid-dd-theta (- (* (sin (* mid-theta *conv-radians*)) scaling)))
        (mid-d-theta-2 (+ *d-theta* (avg first-dd-theta mid-dd-theta)))
        (mid-theta-2 (+ *theta* (avg *d-theta* mid-d-theta-2)))
        ; again first
        (mid-dd-theta-2 (- (* (sin (* mid-theta-2 *conv-radians*)) scaling)))
        (last-d-theta (+ mid-d-theta-2 mid-dd-theta-2))
        (last-theta (+ mid-theta-2 (avg mid-d-theta-2 last-d-theta)))
        ; again second
        (last-dd-theta (- (* (sin (* last-theta *conv-radians*)) scaling)))
        (last-d-theta-2 (+ mid-d-theta-2 (avg mid-dd-theta-2 last-dd-theta)))
        (last-theta-2 (+ mid-theta-2 (avg mid-d-theta-2 last-d-theta-2))))
   ; put values back in globals
   (set! *d-theta* last-d-theta-2)
   (set! *theta* last-theta-2)))
The main event loop and graphics context

(let ((tk (tk-start)))

 (tk/wm 'title tk "Pendulum Animation")
 (let ((canvas (tk 'create-widget 'canvas)))
   ;;; redraw the pendulum on canvas
   ;;; - uses angle and length to compute new (x,y) position of bob
   (define (show-pendulum canvas)
     (let* ((pendulum-angle (* *conv-radians* *theta*))
            (x (+ *home-x* (* *length* (sin pendulum-angle))))
            (y (+ *home-y* (* *length* (cos pendulum-angle)))))
       (canvas 'coords 'rod *home-x* *home-y* x y)
       (canvas 'coords 'bob (- x 15) (- y 15) (+ x 15) (+ y 15))))
   ;;; move the pendulum and repeat after 20ms
   (define (animate)
     (recompute-angle)
     (show-pendulum canvas)
     (tk/after 20 animate))
   ;; layout the canvas
   (tk/grid canvas 'column: 0 'row: 0)
   (canvas 'create 'line 0 25 320 25 'tags: 'plate 'width: 2 'fill: 'grey50)
   (canvas 'create 'oval 155 20 165 30 'tags: 'pivot 'outline: "" 'fill: 'grey50)
   (canvas 'create 'line 1 1 1 1 'tags: 'rod 'width: 3 'fill: 'black)
   (canvas 'create 'oval 1 1 2 2 'tags: 'bob 'outline: 'black 'fill: 'yellow)
   ;; get everything started
   (show-pendulum canvas)
   (tk/after 500 animate)
   (tk-event-loop tk)))

</lang>

Tcl

Works with: Tcl version 8.5
Library: Tk

<lang tcl>package require Tcl 8.5 package require Tk

  1. Make the graphical entities

pack [canvas .c -width 320 -height 200] -fill both -expand 1 .c create line 0 25 320 25 -width 2 -fill grey50 -tags plate .c create line 1 1 1 1 -tags rod -width 3 -fill black .c create oval 1 1 2 2 -tags bob -fill yellow -outline black .c create oval 155 20 165 30 -fill grey50 -outline {} -tags pivot

  1. Set some vars

set points {} set Theta 45.0 set dTheta 0.0 set pi 3.1415926535897933 set length 150 set homeX 160

  1. How to respond to a changing in size of the window

proc resized {width} {

   global homeX
   .c coords plate 0 25 $width 25
   set homeX [expr {$width / 2}]
   .c coords pivot [expr {$homeX-5}] 20 [expr {$homeX+5}] 30
   showPendulum

}

  1. How to actually arrange the pendulum, mapping the model to the display

proc showPendulum {} {

   global Theta dTheta pi length homeX
   set angle [expr {$Theta * $pi/180}]
   set x [expr {$homeX + $length*sin($angle)}]
   set y [expr {25 + $length*cos($angle)}]
   .c coords rod $homeX 25 $x $y
   .c coords bob [expr {$x-15}] [expr {$y-15}] [expr {$x+15}] [expr {$y+15}]

}

  1. The dynamic part of the display

proc recomputeAngle {} {

   global Theta dTheta pi length
   set scaling [expr {3000.0/$length**2}]
   # first estimate
   set firstDDTheta [expr {-sin($Theta * $pi/180)*$scaling}]
   set midDTheta [expr {$dTheta + $firstDDTheta}]
   set midTheta [expr {$Theta + ($dTheta + $midDTheta)/2}]
   # second estimate
   set midDDTheta [expr {-sin($midTheta * $pi/180)*$scaling}]
   set midDTheta [expr {$dTheta + ($firstDDTheta + $midDDTheta)/2}]
   set midTheta [expr {$Theta + ($dTheta + $midDTheta)/2}]
   # Now we do a double-estimate approach for getting the final value
   # first estimate
   set midDDTheta [expr {-sin($midTheta * $pi/180)*$scaling}]
   set lastDTheta [expr {$midDTheta + $midDDTheta}]
   set lastTheta [expr {$midTheta + ($midDTheta + $lastDTheta)/2}]
   # second estimate
   set lastDDTheta [expr {-sin($lastTheta * $pi/180)*$scaling}]
   set lastDTheta [expr {$midDTheta + ($midDDTheta + $lastDDTheta)/2}]
   set lastTheta [expr {$midTheta + ($midDTheta + $lastDTheta)/2}]
   # Now put the values back in our globals
   set dTheta $lastDTheta
   set Theta $lastTheta

}

  1. Run the animation by updating the physical model then the display

proc animate {} {

   global animation
   recomputeAngle
   showPendulum
   # Reschedule
   set animation [after 15 animate]

} set animation [after 500 animate]; # Extra initial delay is visually pleasing

  1. Callback to handle resizing of the canvas

bind .c <Configure> {resized %w}

  1. Callback to stop the animation cleanly when the GUI goes away

bind .c <Destroy> {after cancel $animation}</lang>