Gas Model Python

Pasting this code into Visual Python will create a 3D simulation of our Stirling Engine, with visible gas molecules moving in random paths in a region of the pressure vessel below the displacer.

from visual import *
home = cylinder(pos=(0,-2.5,0), axis=(0,5,0), radius=3, color=color.cyan, opacity=0.2)
displacer = cylinder(pos=(0,1.5,0), axis=(0,1,0), radius=2.9, color=color.blue, opacity=1)
ball1 = sphere(pos=(0,0,0), radius=0.1, color=color.yellow)
ball1.velocity = vector(20,20,20)
ball2 = sphere(pos=(0,1.2,2), radius=0.1, color=color.red)
ball2.velocity = vector(20,-20,20)
ball3 = sphere(pos=(-2,-2.3,0), radius=0.1, color=color.green)
ball3.velocity = vector(20,20,20)
ball4 = sphere(pos=(0.9,1.4,-2), radius=0.1, color=color.magenta)
ball4.velocity = vector(20,-20,20)
ball5 = sphere(pos=(-1.0,0.5,2), radius=0.1, color=color.orange)
ball5.velocity = vector(20,-20,20)
ball6 = sphere(pos=(0,1,-0.9), radius=0.1, color=color.blue)
ball6.velocity = vector(-20,20,20)
ball7 = sphere(pos=(2.1,-1,0.5), radius=0.1, color=color.white)
ball7.velocity = vector(20,20,-20)
ball10 = sphere(pos=(2.3,-0.3,1.7), radius=0.1, color=color.cyan)
ball10.velocity = vector(20,-20,20)
ballA = sphere(pos=(2.3,1.7,1.7), radius=0.1, color=color.blue, opacity=0.5)
ballA.velocity = vector(-20,20,-20)
ballB = sphere(pos=(-0.3,-0.4,0.9), radius=0.1, color=color.orange, opacity=0.5)
ballB.velocity = vector(20,20,20)
deltat = 0.005
t = 0
while True:
rate(100)
ball1.pos = ball1.pos + ball1.velocity*deltat
ball2.pos = ball2.pos + ball2.velocity*deltat
ball3.pos = ball3.pos + ball3.velocity*deltat
ball4.pos = ball4.pos + ball4.velocity*deltat
ball5.pos = ball5.pos + ball5.velocity*deltat
ball6.pos = ball6.pos + ball6.velocity*deltat
ball7.pos = ball7.pos + ball7.velocity*deltat
ball10.pos = ball10.pos + ball10.velocity*deltat
ballA.pos = ballA.pos + ballA.velocity*deltat
ballB.pos = ballB.pos + ballB.velocity*deltat
if sqrt(math.pow(ball1.pos.x, 2) + math.pow(ball1.pos.z, 2))>home.radius:
ball1.velocity.x=-ball1.velocity.x
ball1.velocity.z=-ball1.velocity.z
if ball1.pos.y < -1.5:
ball1.velocity.y = -ball1.velocity.y
if ball1.pos.y > 1.5:
ball1.velocity.y = -ball1.velocity.y
if sqrt(math.pow(ball2.pos.x, 2) + math.pow(ball2.pos.z, 2))>home.radius:
ball2.velocity.x=-ball2.velocity.x
ball2.velocity.z=-ball2.velocity.z
if ball2.pos.y < -1.5:
ball2.velocity.y = -ball2.velocity.y
if ball2.pos.y > 1.5:
ball2.velocity.y = -ball2.velocity.y
if sqrt(math.pow(ball3.pos.x, 2) + math.pow(ball3.pos.z, 2))>home.radius:
ball3.velocity.x=-ball3.velocity.x
ball3.velocity.z=-ball3.velocity.z
if ball3.pos.y < -1.5:
ball3.velocity.y = -ball3.velocity.y
if ball3.pos.y > 1.5:
ball3.velocity.y = -ball3.velocity.y
if sqrt(math.pow(ball4.pos.x, 2) + math.pow(ball4.pos.z, 2))>home.radius:
ball4.velocity.x=-ball4.velocity.x
ball4.velocity.z=-ball4.velocity.z
if ball4.pos.y < -1.5:
ball4.velocity.y = -ball4.velocity.y
if ball4.pos.y > 1.5:
ball4.velocity.y = -ball4.velocity.y
if sqrt(math.pow(ball5.pos.x, 2) + math.pow(ball5.pos.z, 2))>home.radius:
ball5.velocity.x=-ball5.velocity.x
ball5.velocity.z=-ball5.velocity.z
if ball5.pos.y < -1.5:
ball5.velocity.y = -ball5.velocity.y
if ball5.pos.y > 1.5:
ball5.velocity.y = -ball5.velocity.y
if sqrt(math.pow(ball6.pos.x, 2) + math.pow(ball6.pos.z, 2))>home.radius:
ball6.velocity.x=-ball6.velocity.x
ball6.velocity.z=-ball6.velocity.z
if ball6.pos.y < -1.5:
ball6.velocity.y = -ball6.velocity.y
if ball6.pos.y > 1.5:
ball6.velocity.y = -ball6.velocity.y
if sqrt(math.pow(ball7.pos.x, 2) + math.pow(ball7.pos.z, 2))>home.radius:
ball7.velocity.x=-ball7.velocity.x
ball7.velocity.z=-ball7.velocity.z
if ball7.pos.y < -1.5:
ball7.velocity.y = -ball7.velocity.y
if ball7.pos.y > 1.5:
ball7.velocity.y = -ball7.velocity.y
if sqrt(math.pow(ball10.pos.x, 2) + math.pow(ball10.pos.z, 2))>home.radius:
ball10.velocity.x=-ball10.velocity.x
ball10.velocity.z=-ball10.velocity.z
if ball10.pos.y < -1.5:
ball10.velocity.y = -ball10.velocity.y
if ball10.pos.y > 1.5:
ball10.velocity.y = -ball10.velocity.y
if sqrt(math.pow(ballA.pos.x, 2) + math.pow(ballA.pos.z, 2))>home.radius:
ballA.velocity.x=-ballA.velocity.x
ballA.velocity.z=-ballA.velocity.z
if ballA.pos.y < -1.5:
ballA.velocity.y = -ballA.velocity.y
if ballA.pos.y > 1.5:
ballA.velocity.y = -ballA.velocity.y
if sqrt(math.pow(ballB.pos.x, 2) + math.pow(ballB.pos.z, 2))>home.radius:
ballB.velocity.x=-ballB.velocity.x
ballB.velocity.z=-ballB.velocity.z
if ballB.pos.y < -1.5:
ballB.velocity.y = -ballB.velocity.y
if ballB.pos.y > 1.5:
ballB.velocity.y = -ballB.velocity.y