# Mass Spring System:

In this example, we are going to study a simple mass-spring system. This example serves to demonstrate the visualization handling abilities of mechanics-visualization module which is still in development process

## Dependencies:

For this system, following python modules are required:

1) sympy=0.7.2

2) numpy=1.6.2

``````

In [49]:

from sympy.physics.mechanics import *
from sympy import symbols,cos,Eq,sqrt,sin
import viz                       #This module is included with the notebook
from IPython.display import *

``````

That is all we are going to need to carry out this visualization.

Now, we will create symbols for coordinates, speed, force etc of the mass(bob).

``````

In [50]:

q = dynamicsymbols('q')                  # Generalized coordinate
u = dynamicsymbols('u')                  # Generalized speed
f = dynamicsymbols('f')                  # Net Force applied to the box
m = symbols('m')                         # Mass of body
g, k, t, L = symbols('g k t L')               # Gravity,spring constant, time and length of unstretched spring

``````

Now we create an inertial reference frame, assign an origin, and set origin velocity to zero.

``````

In [51]:

I = ReferenceFrame('I')            # Inertial reference frame
O = Point('O')                     # Origin point
O.set_vel(I, 0)                    # Origin's velocity is zero

``````

Now we define the bob as a Particle ..

``````

In [52]:

P = Point('P')                      # Center Of Mass of the box
P.set_pos(O, q * I.y)               # Set the position of P
P.set_vel(I, u * I.y )              # Set the velocity of P
bob = Particle('bob',P,m)

``````

We assign forces acting on the Particle. Coordinate axes are defined to be as : x-axis: running from left to right, y-axis: running from top to bottom, z-axis: running from behind the plane to above the plane.

Following the above notion, Gravity is taken to be in (+I.y) direction, and force due to spring(F=Kx), in -I.y direction ..

``````

In [53]:

f = m*g* I.y - k*(q)*I.y    #Net force on bob.
force = [(P,f)]

``````

Setting up equation of motion for the body:

``````

In [54]:

diff_eq = [q.diff(t) - u]

``````

Now, we call the Kane's Method to generate the equation of motion of the bob

``````

In [55]:

kane = KanesMethod(I, q_ind=[q], u_ind=[u], kd_eqs=diff_eq) # Initialize the object
fr, frstar = kane.kanes_equations(force, [bob])

``````

The Differential equation of the motion of bob is given by the expression: fr + frstar = 0, i.e. fr = frstar

``````

In [56]:

fr, frstar

``````
``````

Out[56]:

([-g*m + k*q(t)], [m*Derivative(u(t), t)])

``````

Now, We need to solve this differential equation, for this we would be utilizing code generator sub-module in sympy mechanics, which is under development, and is not merged in the package.

First packing the parameters in a list.

``````

In [57]:

parameters = [m, k, g, L]     #Mass, stiffness, gravity
parameter_vals = [1, 1, 9.8, 1]

``````

Now we import the code generator module.

``````

In [58]:

import code_gen
right_hand_side = code_gen.numeric_right_hand_side(kane, parameters)

``````

The numeric_right_hand_side method takes the Kane's object and parameters in symbolic form, and returns the differential equation in a numeric form, which can be integrated using odeint() method in SciPy.

Now we integrate right_hand_side using odeint

``````

In [65]:

from scipy.integrate import odeint
x0 = [1,0]        #Initial conditions, q=1,u=0.
t = [i*0.1 for i in range(0,10)]    #Taking 10 time intervals of 0.1 sec
numeric_vals = odeint(right_hand_side, x0, t, args=(parameter_vals,))
numeric_vals

``````
``````

Out[65]:

array([[ 1.        ,  0.        ],
[ 1.04396334,  0.87853405],
[ 1.17541411,  1.7482901 ],
[ 1.39303891,  2.6005778 ],
[ 1.69466326,  3.42688139],
[ 2.07727345,  4.21894471],
[ 2.53704657,  4.96885373],
[ 3.06938873,  5.66911561],
[ 3.66898093,  6.31273356],
[ 4.32983223,  6.89327675]])

``````

The variable contains the positions and velocities, we need only positions for this 1D simulation, so we separate them, and pack them alongwith time.

``````

In [63]:

Q = []
for entries in numeric_vals:
Q.append(entries[0])

simulation_coords = [Q,t]     #Packing positions and time in a single array

``````

## Simulation:

For simulation and visualization, we are going to utilize, the Javascript rendering capabilities of IPython notebooks to render WebGL's. We are using Three.js Javascript library for 3D visualizations.

``````

In [64]:

viz.visualize(simulation_coords)

``````
``````

//Modifying AnimFrame for smart animations for different browsers
window.requestAnimFrame = (function(){
return  window.requestAnimationFrame       ||
window.webkitRequestAnimationFrame ||
window.mozRequestAnimationFrame    ||
window.oRequestAnimationFrame      ||
window.msRequestAnimationFrame     ||
function(/* function */ callback, /* DOMElement */ element){
window.setTimeout(callback, 1000 / 60);
};
})();

// fetching coorinates from Python output
var coordinates = [[1.0, 1.043963344984979, 1.1754141133673666, 1.393038914946831, 1.6946632604568606, 2.0772734459187014, 2.5370465725650542, 3.0693887293725495, 3.6689809279475774, 4.3298322287782103, 5.0453396329483411, 5.8083540269903864, 6.6112516230633744, 7.4460101341520106, 8.3042889290908075, 9.1775123710715363, 10.05695550508608, 10.933831221694406, 11.799378083071952, 12.644947870758436, 13.462091881024476, 14.242645466951734, 14.978809604297208, 15.663228782059802, 16.289064520686569, 16.850063669767206, 17.340620914519064, 17.755834770309377, 18.091556555358487, 18.344431847834439, 18.511934002381899, 18.592389393806581, 18.584994138551849, 18.4898221240379, 18.307824279078766, 18.040819067522001, 17.69147431722892, 17.263280562295861, 16.760516169565314, 16.188204598152172, 15.552064198144958, 14.858451075166629, 14.114295582854174, 13.327033076139031, 12.504529621363243, 11.655003400373227, 10.786942600151775, 9.9090206348995338, 9.0300093260149481, 8.1586915020008028, 7.303773076883961, 6.4737961345647053, 5.6770535301135059, 4.9215060430446975, 4.2147028534381468, 3.56370610475483, 2.9750203422783539, 2.4545275208507893, 2.0074282314625376, 1.6381897415622193, 1.3505013604415244, 1.1472375770664169, 1.0304293357170493, 1.0012437446905, 1.059972416445093, 1.2060285540442111, 1.4379528151063226, 1.7534278887878063, 2.1493016519961436, 2.621618662853384, 3.1656596880458734, 3.775988848905417, 4.4465079390525402, 5.1705173525768693, 5.9407830266204043, 6.7496087214287481, 7.5889129180462751, 8.4503095585569792, 9.3251918617874043, 10.204818292734096, 11.080399920848826, 11.943188218052438, 12.784562497122385, 13.596116018183057, 14.369740008750897, 15.097704675556272, 15.772736431569612, 16.388090590239063, 16.937618724835971, 17.415830140490293, 17.817946708118495, 18.139950597083427, 18.378624470443459, 18.531583563738522, 18.597299557760884, 18.575115858989967, 18.465254090430427, 18.268811973271564, 17.987752296222297, 17.624883268536969], [0.0, 0.1, 0.2, 0.30000000000000004, 0.4, 0.5, 0.6000000000000001, 0.7000000000000001, 0.8, 0.9, 1.0, 1.1, 1.2000000000000002, 1.3, 1.4000000000000001, 1.5, 1.6, 1.7000000000000002, 1.8, 1.9000000000000001, 2.0, 2.1, 2.2, 2.3000000000000003, 2.4000000000000004, 2.5, 2.6, 2.7, 2.8000000000000003, 2.9000000000000004, 3.0, 3.1, 3.2, 3.3000000000000003, 3.4000000000000004, 3.5, 3.6, 3.7, 3.8000000000000003, 3.9000000000000004, 4.0, 4.1000000000000005, 4.2, 4.3, 4.4, 4.5, 4.6000000000000005, 4.7, 4.800000000000001, 4.9, 5.0, 5.1000000000000005, 5.2, 5.300000000000001, 5.4, 5.5, 5.6000000000000005, 5.7, 5.800000000000001, 5.9, 6.0, 6.1000000000000005, 6.2, 6.300000000000001, 6.4, 6.5, 6.6000000000000005, 6.7, 6.800000000000001, 6.9, 7.0, 7.1000000000000005, 7.2, 7.300000000000001, 7.4, 7.5, 7.6000000000000005, 7.7, 7.800000000000001, 7.9, 8.0, 8.1, 8.200000000000001, 8.3, 8.4, 8.5, 8.6, 8.700000000000001, 8.8, 8.9, 9.0, 9.1, 9.200000000000001, 9.3, 9.4, 9.5, 9.600000000000001, 9.700000000000001, 9.8, 9.9]];

var time = coordinates[1][0];
var i=0;
//Setting other varaibles
var controls,scene,camera,renderer,reset;
var \$reset_button;
var \$canvas
var sphere_geometry;

var animation_request_id;

//initial Position, dimensions of bob
var radius = 0.1, segments = 100, rings = 100;

init_canvas();
init();
visualize();

function init_canvas(){

//This function initiates a basic canvas, with trackball controls,
//all drawing work occurs in init() function.

// first of all , a renderer ...
renderer = new THREE.WebGLRenderer();

//show the IPython handle container
container.show();

// create a canvas div
\$canvas = \$("<div/>").attr("id", "#canvas");
//giving background color

\$canvas.attr("style","background-color:rgb(104,104,104)");

//For this particular, giving a top left, position time index

\$canvas.append("<p>click left and move mouse to rotate camera, hit reset button to reset camera</p>");
\$canvas.append('<div id=\"pos_index\" style=\"margin-left:10px;margin-top:20px;\">T= ' + time + '<br /> Pos..='+(-link_length)+'</div>');

// Adding our canvas to IPython UI

// Now lets add a scene ..

// set the scene size
var WIDTH = 400,
HEIGHT = 300;

scene = new THREE.Scene();

//Add a camera to the scene..

// set some camera attributes
var VIEW_ANGLE = 45,
ASPECT = WIDTH / HEIGHT,
NEAR = 0.1,
FAR = 1000;

camera = new THREE.PerspectiveCamera(  VIEW_ANGLE,
ASPECT,
NEAR,
FAR  );

// the camera starts at 0,0,0 so pull it back
camera.position.z = 10;

controls = new THREE.TrackballControls(camera, renderer.domElement);

reset = function(){ controls.reset();}

\$reset_button = \$('<button/>').attr('style','margin-left:40px;').click(reset);
\$reset_button.append('Reset Camera');
\$canvas.append(\$reset_button)

\$canvas.append(renderer.domElement);
element.append(\$canvas);
// start the renderer
renderer.setSize(WIDTH, HEIGHT);

var axesMaterial = new THREE.MeshLambertMaterial(
{
color: 0xFFFFFF

});
var x_axis = new THREE.Mesh(
new THREE.CubeGeometry(WIDTH, 0.03, 0.03),
axesMaterial);

var y_axis = new THREE.Mesh(
new THREE.CubeGeometry(0.03, HEIGHT, 0.03),
axesMaterial);

var z_axis = new THREE.Mesh(
new THREE.CubeGeometry(0.03, 0.03, 10),
axesMaterial);

// create a point light
var pointLight = new THREE.PointLight( 0xFFFFFF );

// set its position
pointLight.position.x = 100;
pointLight.position.y = 100;
pointLight.position.z = 100;

renderer.render(scene, camera);

}

function init(){

// create the sphere's material
var sphereMaterial = new THREE.MeshLambertMaterial(
{
color: 0xCC0000
});

// create a new mesh with sphere geometry

sphere_geometry = new THREE.SphereGeometry(radius, segments, rings)
sphere = new THREE.Mesh(sphere_geometry
,
sphereMaterial);
sphere.position.needsUpdate = true;
sphere.geometry.dynamic = true;

// add the sphere to the scene

// draw!
renderer.render(scene, camera);

}

function visualize() {

controls.update();

renderer.render(scene, camera);

requestAnimationFrame(visualize);

}

var animate  =function() {

time = coordinates[1][i];

console.log(sphere.position.y);
\$("#pos_index").html('T= ' + time + '<br /> Pos..='+(-link_length));

renderer.render(scene, camera);
i+=1;
if(i == coordinates[0].length)
i=0;

animation_request_id = requestAnimationFrame(animate);
}

var stop_animate  =function() {

cancelAnimationFrame(animation_request_id);

}

var \$animate_button = \$('<button/>').attr('style','margin-left:40px;').click(animate);
\$animate_button.append('Animate');
\$canvas.append(\$animate_button)

var \$stop_animate_button = \$('<button/>').attr('style','margin-left:40px;').click(stop_animate);
\$stop_animate_button.append('Stop Animation');
\$canvas.append(\$stop_animate_button)

``````