Friday, June 26, 2015

June 26: Energy graphs

The last time I tried writing some code for adding energy it didn't really work so today I just focused on the accounting.

There's probably a better way to display the code but I'll just dump it here:

from __future__ import division
from visual import *
from visual.graph import *

## INITIAL CONDITIONS
r0 = 3 ## meters, relaxed length of spring
theta0 = math.radians(106.35) ## free angle
h = 5 ## meters, initial height of ball
vx0 = 6 ## initial horizontal velocity

## DEFINE SPRING
initialaxis = vector(r0*math.cos(theta0),r0*math.sin(theta0),0)
spring = helix(pos=(initialaxis.x,h - initialaxis.y,0), axis=initialaxis, radius=0.6,color=color.yellow)
spring.constant = 90 ## N/m

## DEFINE BALL
ball = sphere(pos=(0,h,0), radius=1,make_trail = True)
ball.mass = 1 # kg
ball.velocity = vector(vx0,0,0)
ball.acceleration = vector(0,0,0)
ball.force = vector(0,0,0)
floor = box(size=(50,.1,2),pos=(0,0,0))

## OTHER DEFINITIONS
accelOfGrav = vector(0,-9.8,0) ## acceleration of gravity

K0 = 0.5 * ball.mass * vx0**2
Ugrav0 = ball.mass * 9.8 * h
energy0 = K0 + Ugrav0 ## joules

dt = 0.005
t = 0

## GRAPHING
graph1 = gdisplay()
energygraph = gcurve(color=color.red)
Kgraph = gcurve(color=color.green)
Ugravgraph = gcurve(color=color.white)
Uspringgraph = gcurve(color = color.yellow)
positiongraph = gcurve(color=color.cyan)

while True:
    rate(100)
    t += dt
    ball.force = ball.mass * accelOfGrav
    
    ## STANCE PHASE
    if ball.pos.y < initialaxis.y:## when spring touches floor
        spring.pos.y = 0.1
        spring.axis = ball.pos - spring.pos##update spring length
        spring.displacement = initialaxis.mag*norm(spring.axis) - spring.axis
        ## ENERGY ADDITIONS/SUBTRACTIONS
        dEnergy = energy0 - energy
        springForce = spring.constant * spring.displacement
        ball.force += springForce ## update force on ball

    ## FLIGHT PHASE
    else:
        spring.pos = ball.pos - initialaxis ##update spring position
        spring.axis = initialaxis
        spring.displacement = vector(0,0,0)

    ball.acceleration = ball.force/ball.mass
    ball.velocity += ball.acceleration * dt ##update ball velocity
    ball.pos += ball.velocity * dt ##update ball position

    ## ENERGY ACCOUNTING
    K = 0.5 * ball.mass * ball.velocity.mag**2
    Ugrav = ball.mass * 9.8 * ball.pos.y
    Uspring = 0.5 * spring.constant * spring.displacement.mag**2
    energy = K + Ugrav + Uspring #calculate current energy

    ## UPDATE GRAPHS
    energygraph.plot(pos=(ball.pos.x,energy))
    positiongraph.plot(pos=(ball.pos.x, ball.pos.y))
    Kgraph.plot(pos=(ball.pos.x,K))
    Ugravgraph.plot(pos=(ball.pos.x,Ugrav))
    Uspringgraph.plot(pos=(ball.pos.x,Ugrav))


I created graphs, did the energy accounting section, and made the graphs update. I also went and commented everything when I got slightly frustrated.

I had problems with spring potential energy. Like regarding placement of the code or whatnot. I suppose if I actually study algorithms or programming or something these kinds of problems won't appear.

Red is total energy (K+U). Yellow is spring potential energy. Green is kinetic energy. Gravitational potential energy is supposed to be plotted but it's not there. Blue is just position of the ball (x, y). Everything is plotted to the x-position of the ball, not time.

No comments:

Post a Comment