COLLISIONS EN 3D AVEC VISUAL PYTHON

#! /usr/bin/env python 
# -*- coding: Latin-1 -*- 

"COLLISIONS_VPYTHON"
"Ce script est démarqué de l'exemple fourni dans 'GAS' de l'extension VISUAL de PYTHON"

from visual import *
from time import clock
from visual.graph import *
from random import random
from visual.controls import*

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

win=500
# Valeurs typiques
L=1.0 
gray=(0.7,0.7,0.7) # couleur des bordures
Raxes=0.01 # rayon des bordures
scene=display(title="COLLISIONS_VPYTHON",width=win,height=win,background=(0,1,1),x=0,y=0,range=L/4,center=(L/2.,L/2.,L/2.))
scene.lights=[vector(1,1,1)]
scene.ambient=0
xaxis=curve(pos=[(0,0,0),(L,0,0)],color=gray,radius=Raxes)
yaxis=curve(pos=[(0,0,0),(0,L,0)],color=gray,radius=Raxes)
xaxis2=curve(pos=[(0,L,0),(L,L,0)],color=gray,radius=Raxes)
yaxis2=curve(pos=[(L,0,0),(L,L,0)],color=gray,radius=Raxes)
dt=0.002

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

def boules():
    ""
    global pos,v,radius,boule,N,R,Nsteps,t,dt
    N =10 # valeur à modifier pour augmenter ou diminuer le nombre de boules
    R=0.05 
    boule=[]
    colors=[color.red,color.green,color.blue,color.yellow,color.cyan,color.magenta,color.white]
    poslist,vlist,rlist=[],[],[]
    for i in range(N):
        Lmin=R
        Lmax=L-Lmin
        x,y,z,r=R+(Lmax-Lmin)*random(),R+(Lmax-Lmin)*random(),0,R
        boule=boule+[sphere(pos=(x,y,z),radius=r,color=colors[i % 7])]
        theta,phi=pi*random(),2*pi*random()
        vx,vy,vz=sin(theta)*cos(phi),sin(theta)*sin(phi),0
        poslist.append((x,y,z))
        vlist.append((vx,vy,vz))
        rlist.append(r)
    pos=array(poslist)
    v=array(vlist)
    radius=array(rlist)
    t=0.0
    Nsteps=0
    pos=pos+v*(dt/2.) 
    time=clock()

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

def move_boules():
    ""
    global pos,v,radius,boule,Nsteps,t,dt
    while 1:
        rate(100)
        pos =pos+ v*dt
        r = pos-pos[:,NewAxis] 
        rmag = sqrt(add.reduce(r*r,-1))
        hit = less_equal(rmag,radius+radius[:,NewAxis])-identity(N)
        hitlist = sort(nonzero(hit.flat)).tolist()
       # si collisions:
        for ij in hitlist:
            i,j= divmod(ij,N)
            hitlist.remove(j*N+i)
            vtot=v[i]+v[j]
            vi,vj=v[i],v[j]
            ri,rj=boule[i].radius,boule[j].radius
            a=mag(vj-vi)**2
            if a==0: continue 
            b=2*dot(pos[i]-pos[j],vj-vi)
            c=mag(pos[i]-pos[j])**2-(ri+rj)**2
            d=b**2-4.*a*c
            if d<0: continue 
            deltat=(-b+sqrt(d))/(2.*a) # t-deltat lorsque il y a contact
            pos[i],pos[j]=pos[i]-v[i]*deltat,pos[j]-v[j]*deltat # retour à la configuration de contact
            vcmi,vcmj=v[i]-vtot/2,v[j]-vtot/2 
            rrel=norm(pos[j]-pos[i])
            vcmi,vcmj=vcmi-2*dot(vcmi,rrel)*rrel,vcmj-2*dot(vcmj,rrel)*rrel
            v[i],v[j]=vcmi+vtot/2,vcmj+vtot/2
            pos[i],pos[j]=pos[i]+v[i]*deltat,pos[j]+v[j]*deltat
        # rebonds sur les bordures
        outside=less_equal(pos,R+Raxes)
        v1=v*outside
        v=v-v1+abs(v1)
        outside=greater_equal(pos,L-R-Raxes)
        v1=v*outside
        v=v-v1-abs(v1)
        # mise à jour des positions
        for i in range(N):
            boule[i].pos=pos[i]
        Nsteps=Nsteps+1
        t=t+dt

def lancer():
    ""
    boules()
    move_boules()
lancer()