Version 2 of TclSpringies : A simple mass and spring simulator

Updated 2004-12-07 00:25:34

if 0 {

MBS After visiting the Runge-Kutta page recently, I decided to re-visit something that I was playing with about a year ago...

I got the idea of trying to create a simple mass and spring simulation based on Xspringies [L1 ] after finding a Java based version [L2 ].

I wasn't interested in all the 'bells and whistles' that were included in the Xspringies. I wanted to keep it 'simple'.

The following is what I ended up with. Comments and code fixes are very much welcome.


}

if 0 {

Initialize global values...

I don't remember why I called it 'State', but 'State' is a global array that holds almost everything.

}

 global State

 set State(t)    0.0           ;# current time
 set State(dt)   0.1           ;# time step
 set State(MaxTime) 1000.      ;# Max time value

 # After each time step, the maximum distance that any node travels
 # is calculated.  If this distance, "dmax", is less that "tol", then
 # stop the simulation.

 set State(tol)  0.01

 set State(gravity) -9.8
 set State(curNode) none

if 0 {

Initialize the geometry to be displayed.

'node_def' is a list of nodes containing the following for each node:

* node name * X coordinate * Y coordinate * node mass * free/fixed (whether that node is free to move or fixed) * spring list (the list of springs attached to this node)

'spring_def' is the list of springs containing the following info for each spring:

* spring name * spring constant : k * spring damping constant : kd * rest length for this spring * node list (the list of nodes attached to this spring)

'init_geom' has multiple definitions for 'node_def' and 'spring_def'. Each one is just a 'different' test. Select the one you wish to try. The current setting for 'node_def' and 'spring_def' is sort of a simple suspension bridge type of configuration.

I probably should have added some code so that the user could load 'node_def' and 'spring_def' from a file, but I never got around to doing it...

Oh well.

}

 proc init_geom { } {
     global State


 ########################################################################

     set node_def { \
                        { n1 250.  350.  50. fixed { s1    } } \
                        { n2 300.  350.  50. free  { s1 s2 } } \
                        { n3 250.  300.  50. free  { s2    } } \
                    }

     set spring_def { \
                        { s1 10.0 50.0  50. {n1 n2} } \
                        { s2 10.0 50.0  50. {n2 n3} } \
                      }

 ########################################################################

     set node_def { \
                        { n1  100 350  50. fixed { s1  s5      } } \
                        { n2  150 350  50. free  { s1  s2  s6  } } \
                        { n3  200 350  50. free  { s2  s3  s7  } } \
                        { n4  250 350  50. free  { s3  s4  s8  } } \
                        { n5  300 350  50. fixed { s4  s9      } } \
                        { n6  100 300  50. fixed { s5  s10     } } \
                        { n7  150 300  50. free  { s6  s10 s11 } } \
                        { n8  200 300  50. free  { s7  s11 s12 } } \
                        { n9  250 300  50. free  { s8  s12 s13 } } \
                        { n10 300 300  50. fixed { s9  s13     } } \
                    }

     set spring_def { \
                          { s1  50.0  45.0  50. {n1   n2 } } \
                          { s2  50.0  45.0  50. {n2   n3 } } \
                          { s3  50.0  45.0  50. {n3   n4 } } \
                          { s4  50.0  45.0  50. {n4   n5 } } \
                          { s5  50.0  45.0  50. {n1   n6 } } \
                          { s6  50.0  45.0  50. {n2   n7 } } \
                          { s7  50.0  45.0  50. {n3   n8 } } \
                          { s8  50.0  45.0  50. {n4   n9 } } \
                          { s9  50.0  45.0  50. {n5   n10} } \
                          { s10 50.0  45.0  50. {n6   n7 } } \
                          { s11 50.0  45.0  50. {n7   n8 } } \
                          { s12 50.0  45.0  50. {n8   n9 } } \
                          { s13 50.0  45.0  50. {n9   n10} } \
                      }

 ########################################################################

     set node_def { \
                        { n1  50 350 100. fixed {      s1 } } \
                        { n2 100 350  50. free  { s1   s2 } } \
                        { n3 150 350 100. free  { s2   s3 } } \
                        { n4 200 350  50. free  { s3   s4 } } \
                        { n5 250 350  60. free  { s4   s5 } } \
                        { n6 350 350  70. free  { s5   s6 } } \
                        { n7 450 350 100. fixed { s6      } } \
                    }

     set spring_def { \
                          { s1  30.0  25.0  25. {n1   n2} } \
                          { s2  40.0  35.0  25. {n2   n3} } \
                          { s3  50.0  45.0  25. {n3   n4} } \
                          { s4  50.0  45.0  25. {n4   n5} } \
                          { s5  40.0  35.0  25. {n5   n6} } \
                          { s6  30.0  25.0  25. {n6   n7} } \
                      }

 ########################################################################

     set node_def { \
                        { n1  100 350  50. fixed { s1  s5  s14         } } \
                        { n2  150 350  50. free  { s1  s2  s6  s15 s16 } } \
                        { n3  200 350  50. free  { s2  s3  s7  s17 s18 } } \
                        { n4  250 350  50. free  { s3  s4  s8  s19 s20 } } \
                        { n5  300 350  50. fixed { s4  s9  s21         } } \
                        { n6  100 300  50. fixed { s5  s10 s15         } } \
                        { n7  150 300  50. free  { s6  s10 s11 s14 s17 } } \
                        { n8  200 300  50. free  { s7  s11 s12 s16 s19 } } \
                        { n9  250 300  50. free  { s8  s12 s13 s18 s21 } } \
                        { n10 300 300  50. fixed { s9  s13 s20         } } \
                    }

     set spring_def { \
                          { s1  90.0  45.0  50. {n1   n2 } } \
                          { s2  90.0  45.0  50. {n2   n3 } } \
                          { s3  90.0  45.0  50. {n3   n4 } } \
                          { s4  90.0  45.0  50. {n4   n5 } } \
                          { s5  90.0  45.0  50. {n1   n6 } } \
                          { s6  90.0  45.0  50. {n2   n7 } } \
                          { s7  90.0  45.0  50. {n3   n8 } } \
                          { s8  90.0  45.0  50. {n4   n9 } } \
                          { s9  90.0  45.0  50. {n5   n10} } \
                          { s10 90.0  45.0  50. {n6   n7 } } \
                          { s11 90.0  45.0  50. {n7   n8 } } \
                          { s12 90.0  45.0  50. {n8   n9 } } \
                          { s13 90.0  45.0  50. {n9   n10} } \
                          { s14 90.0  45.0  50. {n1   n7 } } \
                          { s15 90.0  45.0  50. {n2   n6 } } \
                          { s16 90.0  45.0  50. {n2   n8 } } \
                          { s17 90.0  45.0  50. {n3   n7 } } \
                          { s18 90.0  45.0  50. {n3   n9 } } \
                          { s19 90.0  45.0  50. {n4   n8 } } \
                          { s20 90.0  45.0  50. {n4   n10} } \
                          { s21 90.0  45.0  50. {n5   n9 } } \
                      }

 ########################  bridge  ######################################

     set node_def { \
                        { n1   25 300  20. fixed {     s1  s49 } } \
                        { n2   50 300  20. free  { s1  s2  s50 } } \
                        { n3   75 300  20. free  { s2  s3  s51 } } \
                        { n4  100 300  20. free  { s3  s4  s52 } } \
                        { n5  125 300  20. free  { s4  s5  s53 } } \
                        { n6  150 325  20. free  { s5  s6  s54 } } \
                        { n7  175 350  20. free  { s6  s7  s55 } } \
                        { n8  200 375  20. free  { s7  s8  s56 } } \
                        { n9  225 400  20. fixed { s8  s9  s57 } } \
                        { n10 250 380  20. free  { s9  s10 s58 } } \
                        { n11 275 370  20. free  { s10 s11 s59 } } \
                        { n12 300 360  20. free  { s11 s12 s60 } } \
                        { n13 325 350  20. free  { s12 s13 s61 } } \
                        { n14 350 360  20. free  { s13 s14 s62 } } \
                        { n15 375 370  20. free  { s14 s15 s63 } } \
                        { n16 400 380  20. free  { s15 s16 s64 } } \
                        { n17 425 400  20. fixed { s16 s17 s65 } } \
                        { n18 450 375  20. free  { s17 s18 s66 } } \
                        { n19 475 350  20. free  { s18 s19 s67 } } \
                        { n20 500 325  20. free  { s19 s20 s68 } } \
                        { n21 525 300  20. free  { s20 s21 s69 } } \
                        { n22 550 300  20. free  { s21 s22 s70 } } \
                        { n23 575 300  20. free  { s22 s23 s71 } } \
                        { n24 600 300  20. free  { s23 s24 s72 } } \
                        { n25 625 300  20. fixed { s24     s73 } } \

                        { n26  25 200  20. fixed {     s25 s49 } } \
                        { n27  50 200  20. free  { s25 s26 s50 } } \
                        { n28  75 200  20. free  { s26 s27 s51 } } \
                        { n29 100 200  20. free  { s27 s28 s52 } } \
                        { n30 125 200  20. free  { s28 s29 s53 } } \
                        { n31 150 200  20. free  { s29 s30 s54 } } \
                        { n32 175 200  20. free  { s30 s31 s55 } } \
                        { n33 200 200  20. free  { s31 s32 s56 } } \
                        { n34 225 200  20. free  { s32 s33 s57 } } \
                        { n35 250 200  20. free  { s33 s34 s58 } } \
                        { n36 275 200  20. free  { s34 s35 s59 } } \
                        { n37 300 200  20. free  { s35 s36 s60 } } \
                        { n38 325 200  20. free  { s36 s37 s61 } } \
                        { n39 350 200  20. free  { s37 s38 s62 } } \
                        { n40 375 200  20. free  { s38 s39 s63 } } \
                        { n41 400 200  20. free  { s39 s40 s64 } } \
                        { n42 425 200  20. free  { s40 s41 s65 } } \
                        { n43 450 200  20. free  { s41 s42 s66 } } \
                        { n44 475 200  20. free  { s42 s43 s67 } } \
                        { n45 500 200  20. free  { s43 s44 s68 } } \
                        { n46 525 200  20. free  { s44 s45 s69 } } \
                        { n47 550 200  20. free  { s45 s46 s70 } } \
                        { n48 575 200  20. free  { s46 s47 s71 } } \
                        { n49 600 200  20. free  { s47 s48 s72 } } \
                        { n50 625 200  20. fixed { s48     s73 } } \

                    }

     set spring_def { \
                          { s1  90.0  45.0  10. {n1   n2 } } \
                          { s2  90.0  45.0  10. {n2   n3 } } \
                          { s3  90.0  45.0  10. {n3   n4 } } \
                          { s4  90.0  45.0  10. {n4   n5 } } \
                          { s5  90.0  45.0  10. {n5   n6 } } \
                          { s6  90.0  45.0  10. {n6   n7 } } \
                          { s7  90.0  45.0  10. {n7   n8 } } \
                          { s8  90.0  45.0  10. {n8   n9 } } \
                          { s9  90.0  45.0  10. {n9   n10} } \
                          { s10 90.0  45.0  10. {n10  n11} } \
                          { s11 90.0  45.0  10. {n11  n12} } \
                          { s12 90.0  45.0  10. {n12  n13} } \
                          { s13 90.0  45.0  10. {n13  n14} } \
                          { s14 90.0  45.0  10. {n14  n15} } \
                          { s15 90.0  45.0  10. {n15  n16} } \
                          { s16 90.0  45.0  10. {n16  n17} } \
                          { s17 90.0  45.0  10. {n17  n18} } \
                          { s18 90.0  45.0  10. {n18  n19} } \
                          { s19 90.0  45.0  10. {n19  n20} } \
                          { s20 90.0  45.0  10. {n20  n21} } \
                          { s21 90.0  45.0  10. {n21  n22} } \
                          { s22 90.0  45.0  10. {n22  n23} } \
                          { s23 90.0  45.0  10. {n23  n24} } \
                          { s24 90.0  45.0  10. {n24  n25} } \


                          { s25 90.0  45.0  10. {n26  n27} } \
                          { s26 90.0  45.0  10. {n27  n28} } \
                          { s27 90.0  45.0  10. {n28  n29} } \
                          { s28 90.0  45.0  10. {n29  n30} } \
                          { s29 90.0  45.0  10. {n30  n31} } \
                          { s30 90.0  45.0  10. {n31  n32} } \
                          { s31 90.0  45.0  10. {n32  n33} } \
                          { s32 90.0  45.0  10. {n33  n34} } \
                          { s33 90.0  45.0  10. {n34  n35} } \
                          { s34 90.0  45.0  10. {n35  n36} } \
                          { s35 90.0  45.0  10. {n36  n37} } \
                          { s36 90.0  45.0  10. {n37  n38} } \
                          { s37 90.0  45.0  10. {n38  n39} } \
                          { s38 90.0  45.0  10. {n39  n40} } \
                          { s39 90.0  45.0  10. {n40  n41} } \
                          { s40 90.0  45.0  10. {n41  n42} } \
                          { s41 90.0  45.0  10. {n42  n43} } \
                          { s42 90.0  45.0  10. {n43  n44} } \
                          { s43 90.0  45.0  10. {n44  n45} } \
                          { s44 90.0  45.0  10. {n45  n46} } \
                          { s45 90.0  45.0  10. {n46  n47} } \
                          { s46 90.0  45.0  10. {n47  n48} } \
                          { s47 90.0  45.0  10. {n48  n49} } \
                          { s48 90.0  45.0  10. {n49  n50} } \

                          { s49 90.0  45.0 100. {n1   n26} } \
                          { s50 90.0  45.0  86. {n2   n27} } \
                          { s51 90.0  45.0  84. {n3   n28} } \
                          { s52 90.0  45.0  85. {n4   n29} } \
                          { s53 90.0  45.0 102. {n5   n30} } \
                          { s54 90.0  45.0 115. {n6   n31} } \
                          { s55 90.0  45.0 132. {n7   n32} } \
                          { s56 90.0  45.0 165. {n8   n33} } \
                          { s57 90.0  45.0 190. {n9   n34} } \
                          { s58 90.0  45.0 174. {n10  n35} } \
                          { s59 90.0  45.0 160. {n11  n36} } \
                          { s60 90.0  45.0 153. {n12  n37} } \
                          { s61 90.0  45.0 145. {n13  n38} } \
                          { s62 90.0  45.0 153. {n14  n39} } \
                          { s63 90.0  45.0 160. {n15  n40} } \
                          { s64 90.0  45.0 174. {n16  n41} } \
                          { s65 90.0  45.0 190. {n17  n42} } \
                          { s66 90.0  45.0 165. {n18  n43} } \
                          { s67 90.0  45.0 132. {n19  n44} } \
                          { s68 90.0  45.0 115. {n20  n45} } \
                          { s69 90.0  45.0 102. {n21  n46} } \
                          { s70 90.0  45.0  85. {n22  n47} } \
                          { s71 90.0  45.0  84. {n23  n48} } \
                          { s72 90.0  45.0  86. {n24  n49} } \
                          { s73 90.0  45.0 100. {n25  n50} } \
                      }


     set State(node_list)   {}
     set State(spring_list) {}

     foreach s $spring_def {
         foreach {spring k damp rest nlist} $s {
             set State($spring,name)    $spring
             set State($spring,k)       $k
             set State($spring,kd)      $damp
             set State($spring,rest)    $rest
             set State($spring,nlist)   $nlist
             set State($spring,x1)      0.
             set State($spring,y1)      0.
             set State($spring,x2)      0.
             set State($spring,y2)      0.
             set State($spring,dist)    0.
             set State($spring,force)   0.
             foreach {n1 n2} $nlist {
                 set State($spring,n1)  $n1
                 set State($spring,n2)  $n2
             }
             lappend State(spring_list) $spring
         }
     }

     foreach n $node_def {
         foreach {node x y m type slist} $n {
             set State($node,name)  $node
             set State($node,x)     $x
             set State($node,y)     $y
             set State($node,rad)   [expr $m * 0.1]
             set State($node,mass)  $m
             set State($node,type)  $type
             set State($node,vx)    0.0
             set State($node,vy)    0.0
             set State($node,ax)    0.0
             set State($node,ay)    0.0
             set State($node,Fx)    0.0
             set State($node,Fy)    0.0
             set State($node,slist) $slist
             lappend State(node_list) $node
         }
     }

     foreach s $State(spring_list) {
         set n1 $State($s,n1)
         set n2 $State($s,n2)
         set State($s,x1) $State($n1,x)
         set State($s,y1) $State($n1,y)
         set State($s,x2) $State($n2,x)
         set State($s,y2) $State($n2,y)
     }

 }

if 0 {

This routine loops over each spring and calculates the force due to the displacement of each node.

}

 proc solve {func} {
     global State

     if { $State(running) == 0 } {
         return
     }

     set funcq    [uplevel 1 namespace which -command $func]

     set done 0

     while { $done == 0 } {

         set t        $State(t)
         set dmax 0.0

         update_spring_forces
         apply_spring_force_to_nodes

         foreach n $State(node_list) {

             set xvec [list $State($n,vx) $State($n,vy) $State($n,ax) $State($n,ay)]
             set xresult [rungeKuttaStep $n $State(t) $State(dt) $xvec rkfunc]
             foreach {nvx nvy xa ya} $xresult {break}

             set yvec [list $State($n,x) $State($n,y) $nvx $nvy]
             set yresult [rungeKuttaStep $n $State(t) $State(dt) $yvec rkfunc ]
             foreach {nx ny nvx nvy} $yresult {break}

             set m [apply_forces $n $nx $ny $nvx $nvy]

             if {[expr {abs($m)}] > $dmax} {
                 set dmax [expr {abs($m)}]
             }
         }

         set State(dmax) $dmax

         set State(t) [expr {$State(t) + $State(dt)}]

         update_geom

         if {$dmax < $State(tol) || $State(t) > $State(MaxTime)} {
             set done 1
         }
         if { $State(running) == 0 } {
             set done 1
         }
     }
 }

if 0 {

Calculate the derivatives

}

 proc rkfunc {n t xvec} {

     foreach {x y vx vy} $xvec {break}
     return [list $vx $vy 0.0 0.0]
 }

if 0 {

This routine loops over each spring and calculates the force due to the displacement of each node.

}

 proc update_spring_forces { } {
     global State

     zero_node_forces

     foreach s $State(spring_list) {

         set dx     [expr { $State($s,x2) - $State($s,x1) } ]
         set dy     [expr { $State($s,y2) - $State($s,y1) } ]
         set dist   [expr { sqrt($dx * $dx + $dy * $dy) } ]

         set State($s,dist)  $dist
         set State($s,force) [expr { -$State($s,k) * ($dist - $State($s,rest))}]
     }
 }

if 0 {

Set the forces applied to each node to zero.

}

 proc zero_node_forces { } {
     global State

     foreach n $State(node_list) {
         set State($n,Fx) 0.0
         set State($n,Fy) 0.0
     }
 }

if 0 {

Calculate the forces,accelerations,velocities,etc. to be applied to each node.

}

 proc apply_spring_force_to_nodes { } {
     global State

     foreach n $State(node_list) {

         set State($n,Fy) [expr {$State(gravity) * $State($n,mass)}]

         foreach s $State($n,slist) {

             if {$n == $State($s,n1)} {
                 set n1 $State($s,n1)
                 set n2 $State($s,n2)
             } else {
                 set n1 $State($s,n2)
                 set n2 $State($s,n1)
             }

             set dx    [expr { $State($n1,x) - $State($n2,x) } ]
             set dy    [expr { $State($n1,y) - $State($n2,y) } ]
             set dvx   [expr { $State($n1,vx) - $State($n2,vx) } ]
             set dvy   [expr { $State($n1,vy) - $State($n2,vy) } ]
             set VdotL [expr { $dvx * $dx + $dvy * $dy } ]

             set dist  $State($s,dist)
             set damp  [expr { $State($s,kd) * $VdotL / $dist }]
             set force [expr { $State($s,force) - $damp } ]

             set Fx [expr { $force * $dx / $dist } ]
             set Fy [expr { $force * $dy / $dist} ]

             set State($n,Fx) [expr {$State($n,Fx) + $Fx}]
             set State($n,Fy) [expr {$State($n,Fy) + $Fy}]
         }
         set State($n,ax) [expr {$State($n,Fx) / $State($n,mass)}]
         set State($n,ay) [expr {$State($n,Fy) / $State($n,mass)}]
     }
 }

if 0 {

I borrowed this routine from tcllib. Specifically I borrowed the math::calculus::rungeKuttaStep routine so that I could modify it slightly. In the process, I believe that I found a "typo/bug" in the tcllib routine (search MBS below for more info). The bug has been reported and will be fixed in the next release of tcllib (Note: It was in version 1.7).

}

 ########################################################################

 # rungeKuttaStep --
 #    Integrate a system of ordinary differential equations of the type
 #    x' = f(x,t), where x is a vector of quantities. Integration is
 #    done over a single step according to Runge-Kutta 4th order.
 #
 # Arguments:
 #    t           Start value of independent variable (time for instance)
 #    tstep       Step size of interval
 #    xvec        Vector of dependent values at the start
 #    func        Function taking the arguments t and xvec to return
 #                the derivative of each dependent variable.
 # Return value:
 #    List of values at the end of the step
 #
 proc rungeKuttaStep { n t tstep xvec func } {

    set funcq    [uplevel 1 namespace which -command $func]

    #
    # Four steps:
    # - k1 = tstep*func(t,x0)
    # - k2 = tstep*func(t+0.5*tstep,x0+0.5*k1)
    # - k3 = tstep*func(t+0.5*tstep,x0+0.5*k2)
    # - k4 = tstep*func(t+    tstep,x0+    k3)
    # - x1 = x0 + (k1+2*k2+2*k3+k4)/6
    #
    set tstep2   [expr {$tstep/2.0}]
    set tstep6   [expr {$tstep/6.0}]

    set xk1      [$funcq $n $t $xvec]
    set xvec2    {}
    foreach x1 $xvec xv $xk1 {
       lappend xvec2 [expr {$x1+$tstep2*$xv}]
    }

    set xk2      [$funcq $n [expr {$t+$tstep2}] $xvec2]
    set xvec3    {}
    foreach x1 $xvec xv $xk2 {
       lappend xvec3 [expr {$x1+$tstep2*$xv}]
    }

    set xk3      [$funcq $n [expr {$t+$tstep2}] $xvec3]
    set xvec4    {}
    foreach x1 $xvec xv $xk3 {
       lappend xvec4 [expr {$x1+$tstep *$xv}]
    }
    #***************************************************************************
    # MBS :
    # Previously the above line had:
    #   lappend xvec4 [expr {$x1+$tstep2*$xv}]
    # The "bug" is that $tstep2 in the above line should actually be $tstep
    # From the description above
    #           # - k4 = tstep*func(t+    tstep,x0+    k3)
    # tstep is used rather than tstep*0.5 (ie tstep2)
    #***************************************************************************


    set xk4      [$funcq $n [expr {$t+$tstep}] $xvec4]
    set result   {}
    foreach x0 $xvec k1 $xk1 k2 $xk2 k3 $xk3 k4 $xk4 {
       set dx [expr {$k1+2.0*$k2+2.0*$k3+$k4}]
       lappend result [expr {$x0+$dx*$tstep6}]
    }
    return $result
 }

if 0 {

Apply the newly calculated positions and velocities of each node.

}

 proc apply_forces {n nx ny nvx nvy} {
     global State

     set State($n,dx) 0.0
     set State($n,dy) 0.0

     if { $State($n,type) == "fixed" } {
         return 0.0
     }

     set dx [expr { $nx - $State($n,x) } ]
     set dy [expr { $ny - $State($n,y) } ]

     set State($n,dx) $dx
     set State($n,dy) $dy

     set State($n,x) $nx
     set State($n,y) $ny

     set State($n,vx) $nvx
     set State($n,vy) $nvy

     foreach s $State($n,slist) {
         if {$n == $State($s,n1)} {
             set State($s,x1)     $State($n,x)
             set State($s,y1)     $State($n,y)
         } else {
             set State($s,x2)     $State($n,x)
             set State($s,y2)     $State($n,y)
         }
     }
     return [expr {sqrt($dx*$dx + $dy*$dy)} ]
 }

if 0 {

Initialize the interface :

Hopefully this is pretty straight forward...

}

 proc init_display { } {
     global State

     # Window things
     wm title . "Springy balls"
     wm geometry . "+0+0"
     wm protocol . WM_DELETE_WINDOW "file_quit"

     # Create canvas
     set State(canvas) [canvas .c -width 700 -height 700 -background black]
     pack $State(canvas) -fill both -expand true

     bind $State(canvas) <ButtonPress-1>   {c_select %W %x %y}
     bind $State(canvas) <B1-Motion>       {c_drag   %W %x %y}
     bind $State(canvas) <ButtonRelease-1> {c_drop   %W %x %y}

     frame .fr -bd 2 -relief raised
     pack .fr -side bottom -fill x
     button .fr.quit  -text "Quit"  -command {exit}
     button .fr.start -text "Start" -command {set State(running) 1; solve rkfunc}
     button .fr.stop  -text "Stop"  -command {set State(running) 0}
     button .fr.reset -text "Reset" -command {reset}
     label .fr.mtlbl -text "Max Time "
     entry .fr.maxt -textvariable State(MaxTime) -width 5

     label .fr.ltol -text "Tolerance "
     entry .fr.etol -textvariable State(tol) -width 8

     pack .fr.quit  -side left
     pack .fr.start -side left
     pack .fr.stop  -side left
     pack .fr.reset -side left
     pack .fr.mtlbl -side left
     pack .fr.maxt  -side left
     pack .fr.ltol  -side left
     pack .fr.etol  -side left


     frame .f1 -bd 2 -relief raised
     label .f1.tlab -text "Time : "
     label .f1.tval -textvariable State(t) -background white -width 20

     pack .f1.tlab -side left
 #    pack .f1.tval -side left -fill x -expand yes
     pack .f1.tval -side left -fill x 

     frame .f2 -bd 2 -relief raised
     label .f2.dlab -text "dMax : "
     label .f2.dval -textvariable State(dmax) -background white -width 20

     pack .f2.dlab -side left
 #    pack .f2.dval -side left -fill x -expand yes
     pack .f2.dval -side left -fill x

     pack .f1 .f2 -side left -fill x -expand yes

     update
 }

if 0 {

Display the current geometry :

}

 proc draw_geom { } {
     global State

     set canvasHeight [winfo height $State(canvas)]
     set canvasWidth  [winfo width $State(canvas)]

     foreach s $State(spring_list) {
         set x1 $State($s,x1)
         set y1 [expr $canvasHeight - $State($s,y1)]
         set x2 $State($s,x2)
         set y2 [expr $canvasHeight - $State($s,y2)]
         set State($s,id) [$State(canvas) create line $x1 $y1 $x2 $y2 \
                               -fill white -tag $s]
     }

     foreach n $State(node_list) {
         set r  $State($n,rad)
         set x1 [expr                 $State($n,x) - $r]
         set y1 [expr $canvasHeight - $State($n,y) - $r]
         set x2 [expr                 $State($n,x) + $r]
         set y2 [expr $canvasHeight - $State($n,y) + $r]
         if { $State($n,type) == "fixed" } {
             set State($n,id) [$State(canvas) create oval $x1 $y1 $x2 $y2 \
                                   -outline white -fill red -tag $n ]
         } else {
             set State($n,id) [$State(canvas) create oval $x1 $y1 $x2 $y2 \
                                   -outline white -fill blue -tag $n ]
         }
     }
     update
 }

if 0 {

Update the geometry after each time step

}

 proc update_geom { } {
     global State

     set canvasHeight [winfo height $State(canvas)]
     set canvasWidth  [winfo width $State(canvas)]

     foreach s $State(spring_list) {
         set x1 $State($s,x1)
         set y1 [expr $canvasHeight - $State($s,y1)]
         set x2 $State($s,x2)
         set y2 [expr $canvasHeight - $State($s,y2)]
         $State(canvas) coords $State($s,id) $x1 $y1 $x2 $y2 
     }

     foreach n $State(node_list) {
         set dy [expr -1.0 * $State($n,dy)]
         $State(canvas) move $State($n,id) $State($n,dx) $dy
     }
     update
 }

if 0 {

Some canvas related routines...

Allows the user to select / drag / drop a node using Button-One

}

 proc c_closest_node {w x y} {
     global State

     set x [$w canvasx $x]
     set y [$w canvasy $y]

     set all_list [$w find all]

     set mdist 9999999.
     set pt_num {}

     foreach item $all_list {

         set ptr  [$w itemcget $item -tags]
         set lndx [lsearch $ptr "n*"]

         if {$lndx != -1} {
             set cord [$w coords $item]
             set nx [expr ([lindex $cord 0] + [lindex $cord 2]) / 2.0]
             set ny [expr ([lindex $cord 1] + [lindex $cord 3]) / 2.0]
             set dx [expr $x - $nx]
             set dy [expr $y - $ny]
             set dst [expr sqrt($dx*$dx + $dy*$dy)]
             if {$dst < $mdist} {
                 set it    $ptr
                 set mdist $dst
             }
         }
     }
     set node [lindex $it 0]
     return $node
 }

 proc c_select {w x y} {
     global State

     set canvasHeight [winfo height $State(canvas)]
     set canvasWidth  [winfo width $State(canvas)]

     set node [c_closest_node $w $x $y]

     set State(curNode) $node
     if {$State(curNode) != "none"} {
         c_move_node $State(curNode) $x $y
     } else {
         set State(curNode) "none"
     }
 }

 proc c_drag {w x y} {
     global State

     set canvasHeight [winfo height $State(canvas)]
     set canvasWidth  [winfo width $State(canvas)]

     if {$State(curNode) != "none"} {
         c_move_node $State(curNode) $x $y
     }
 }

 proc c_drop {w x y} {
     global State

     set canvasHeight [winfo height $State(canvas)]
     set canvasWidth  [winfo width $State(canvas)]

     if {$State(curNode) != "none"} {
         c_move_node $State(curNode) $x $y
     }
     set State(curNode) "none"

 }

 proc c_move_node {n x y} {
     global State

     set canvasHeight [winfo height $State(canvas)]
     set canvasWidth  [winfo width $State(canvas)]

     set r           $State($n,rad)
     set State($n,x) $x
     set State($n,y) [expr {$canvasHeight - $y}]
     set x1 [expr                 $State($n,x) - $r]
     set y1 [expr $canvasHeight - $State($n,y) - $r]
     set x2 [expr                 $State($n,x) + $r]
     set y2 [expr $canvasHeight - $State($n,y) + $r]
     $State(canvas) coords $State($n,id) $x1 $y1 $x2 $y2 

     foreach s $State($n,slist) {
         if {$n == $State($s,n1)} {
             set scoords [$State(canvas) coords $State($s,id)]
             foreach {x1 y1 x2 y2} $scoords {break}
             set State($s,x1)     $x
             set State($s,y1)     [expr {$canvasHeight - $y}]
             $State(canvas) coords $State($s,id) $x $y $x2 $y2
         } else {
             set scoords [$State(canvas) coords $State($s,id)]
             foreach {x1 y1 x2 y2} $scoords {break}
             set State($s,x2)     $x
             set State($s,y2)     [expr {$canvasHeight - $y}]
             $State(canvas) coords $State($s,id) $x1 $y1 $x $y
         }
     }

 }

if 0 {

Reset to the initial geometry.

}

 proc reset { } {
     global State

     set State(t) 0.0
     init_geom

     $State(canvas) delete all
     draw_geom
     update_geom
 }


 ########################################################################

 init_display

 init_geom

 draw_geom

 set State(running) 0

 solve rkfunc

if 0 {


}