1# pendulum.tcl -- 2# 3# This demonstration illustrates how Tcl/Tk can be used to construct 4# simulations of physical systems. 5# 6# RCS: @(#) $Id$ 7 8if {![info exists widgetDemo]} { 9 error "This script should be run from the \"widget\" demo." 10} 11 12package require Tk 13 14set w .pendulum 15catch {destroy $w} 16toplevel $w 17wm title $w "Pendulum Animation Demonstration" 18wm iconname $w "pendulum" 19positionWindow $w 20 21label $w.msg -font $font -wraplength 4i -justify left -text "This demonstration shows how Tcl/Tk can be used to carry out animations that are linked to simulations of physical systems. In the left canvas is a graphical representation of the physical system itself, a simple pendulum, and in the right canvas is a graph of the phase space of the system, which is a plot of the angle (relative to the vertical) against the angular velocity. The pendulum bob may be repositioned by clicking and dragging anywhere on the left canvas." 22pack $w.msg 23 24## See Code / Dismiss buttons 25set btns [addSeeDismiss $w.buttons $w] 26pack $btns -side bottom -fill x 27 28# Create some structural widgets 29pack [panedwindow $w.p] -fill both -expand 1 30$w.p add [labelframe $w.p.l1 -text "Pendulum Simulation"] 31$w.p add [labelframe $w.p.l2 -text "Phase Space"] 32 33# Create the canvas containing the graphical representation of the 34# simulated system. 35canvas $w.c -width 320 -height 200 -background white -bd 2 -relief sunken 36$w.c create text 5 5 -anchor nw -text "Click to Adjust Bob Start Position" 37# Coordinates of these items don't matter; they will be set properly below 38$w.c create line 0 25 320 25 -tags plate -fill grey50 -width 2 39$w.c create oval 155 20 165 30 -tags pivot -fill grey50 -outline {} 40$w.c create line 1 1 1 1 -tags rod -fill black -width 3 41$w.c create oval 1 1 2 2 -tags bob -fill yellow -outline black 42pack $w.c -in $w.p.l1 -fill both -expand true 43 44# Create the canvas containing the phase space graph; this consists of 45# a line that gets gradually paler as it ages, which is an extremely 46# effective visual trick. 47canvas $w.k -width 320 -height 200 -background white -bd 2 -relief sunken 48$w.k create line 160 200 160 0 -fill grey75 -arrow last -tags y_axis 49$w.k create line 0 100 320 100 -fill grey75 -arrow last -tags x_axis 50for {set i 90} {$i>=0} {incr i -10} { 51 # Coordinates of these items don't matter; they will be set properly below 52 $w.k create line 0 0 1 1 -smooth true -tags graph$i -fill grey$i 53} 54# FIXME: UNICODE labels 55$w.k create text 0 0 -anchor ne -text "q" -font {Symbol 8} -tags label_theta 56$w.k create text 0 0 -anchor ne -text "dq" -font {Symbol 8} -tags label_dtheta 57pack $w.k -in $w.p.l2 -fill both -expand true 58 59# Initialize some variables 60set points {} 61set Theta 45.0 62set dTheta 0.0 63set pi 3.1415926535897933 64set length 150 65set home 160 66 67# This procedure makes the pendulum appear at the correct place on the 68# canvas. If the additional arguments "at $x $y" are passed (the 'at' 69# is really just syntactic sugar) instead of computing the position of 70# the pendulum from the length of the pendulum rod and its angle, the 71# length and angle are computed in reverse from the given location 72# (which is taken to be the centre of the pendulum bob.) 73proc showPendulum {canvas {at {}} {x {}} {y {}}} { 74 global Theta dTheta pi length home 75 if {$at eq "at" && ($x!=$home || $y!=25)} { 76 set dTheta 0.0 77 set x2 [expr {$x - $home}] 78 set y2 [expr {$y - 25}] 79 set length [expr {hypot($x2, $y2)}] 80 set Theta [expr {atan2($x2, $y2) * 180/$pi}] 81 } else { 82 set angle [expr {$Theta * $pi/180}] 83 set x [expr {$home + $length*sin($angle)}] 84 set y [expr {25 + $length*cos($angle)}] 85 } 86 $canvas coords rod $home 25 $x $y 87 $canvas coords bob \ 88 [expr {$x-15}] [expr {$y-15}] [expr {$x+15}] [expr {$y+15}] 89} 90showPendulum $w.c 91 92# Update the phase-space graph according to the current angle and the 93# rate at which the angle is changing (the first derivative with 94# respect to time.) 95proc showPhase {canvas} { 96 global Theta dTheta points psw psh 97 lappend points [expr {$Theta+$psw}] [expr {-20*$dTheta+$psh}] 98 if {[llength $points] > 100} { 99 set points [lrange $points end-99 end] 100 } 101 for {set i 0} {$i<100} {incr i 10} { 102 set list [lrange $points end-[expr {$i-1}] end-[expr {$i-12}]] 103 if {[llength $list] >= 4} { 104 $canvas coords graph$i $list 105 } 106 } 107} 108 109# Set up some bindings on the canvases. Note that when the user 110# clicks we stop the animation until they release the mouse 111# button. Also note that both canvases are sensitive to <Configure> 112# events, which allows them to find out when they have been resized by 113# the user. 114bind $w.c <Destroy> { 115 after cancel $animationCallbacks(pendulum) 116 unset animationCallbacks(pendulum) 117} 118bind $w.c <1> { 119 after cancel $animationCallbacks(pendulum) 120 showPendulum %W at %x %y 121} 122bind $w.c <B1-Motion> { 123 showPendulum %W at %x %y 124} 125bind $w.c <ButtonRelease-1> { 126 showPendulum %W at %x %y 127 set animationCallbacks(pendulum) [after 15 repeat [winfo toplevel %W]] 128} 129bind $w.c <Configure> { 130 %W coords plate 0 25 %w 25 131 set home [expr %w/2] 132 %W coords pivot [expr $home-5] 20 [expr $home+5] 30 133} 134bind $w.k <Configure> { 135 set psh [expr %h/2] 136 set psw [expr %w/2] 137 %W coords x_axis 2 $psh [expr %w-2] $psh 138 %W coords y_axis $psw [expr %h-2] $psw 2 139 %W coords label_dtheta [expr $psw-4] 6 140 %W coords label_theta [expr %w-6] [expr $psh+4] 141} 142 143# This procedure is the "business" part of the simulation that does 144# simple numerical integration of the formula for a simple rotational 145# pendulum. 146proc recomputeAngle {} { 147 global Theta dTheta pi length 148 set scaling [expr {3000.0/$length/$length}] 149 150 # To estimate the integration accurately, we really need to 151 # compute the end-point of our time-step. But to do *that*, we 152 # need to estimate the integration accurately! So we try this 153 # technique, which is inaccurate, but better than doing it in a 154 # single step. What we really want is bound up in the 155 # differential equation: 156 # .. - sin theta 157 # theta + theta = ----------- 158 # length 159 # But my math skills are not good enough to solve this! 160 161 # first estimate 162 set firstDDTheta [expr {-sin($Theta * $pi/180)*$scaling}] 163 set midDTheta [expr {$dTheta + $firstDDTheta}] 164 set midTheta [expr {$Theta + ($dTheta + $midDTheta)/2}] 165 # second estimate 166 set midDDTheta [expr {-sin($midTheta * $pi/180)*$scaling}] 167 set midDTheta [expr {$dTheta + ($firstDDTheta + $midDDTheta)/2}] 168 set midTheta [expr {$Theta + ($dTheta + $midDTheta)/2}] 169 # Now we do a double-estimate approach for getting the final value 170 # first estimate 171 set midDDTheta [expr {-sin($midTheta * $pi/180)*$scaling}] 172 set lastDTheta [expr {$midDTheta + $midDDTheta}] 173 set lastTheta [expr {$midTheta + ($midDTheta + $lastDTheta)/2}] 174 # second estimate 175 set lastDDTheta [expr {-sin($lastTheta * $pi/180)*$scaling}] 176 set lastDTheta [expr {$midDTheta + ($midDDTheta + $lastDDTheta)/2}] 177 set lastTheta [expr {$midTheta + ($midDTheta + $lastDTheta)/2}] 178 # Now put the values back in our globals 179 set dTheta $lastDTheta 180 set Theta $lastTheta 181} 182 183# This method ties together the simulation engine and the graphical 184# display code that visualizes it. 185proc repeat w { 186 global animationCallbacks 187 188 # Simulate 189 recomputeAngle 190 191 # Update the display 192 showPendulum $w.c 193 showPhase $w.k 194 195 # Reschedule ourselves 196 set animationCallbacks(pendulum) [after 15 [list repeat $w]] 197} 198# Start the simulation after a short pause 199set animationCallbacks(pendulum) [after 500 [list repeat $w]] 200