Spring.

 Bahan baru buat ide jahil besok, :)

#code
from visual import *
from random import uniform,random
from visual.controls import *

display(center=(0,0,0),background=(1,1,1), autoscale=False, range=(3,5,3),        
        width=600, height=600,  forward=(-1.4,-1.3,-1)) #arah kamera

distant_light(direction=(1,1,1), color=color.red)
spring = helix(pos=(0,2,1), axis=(5,0,0), radius=0.5, color=color.red)

dl  = .1
l   = 1.
def proses():
    global l,dl
    spring.length = l
    if l>2:
        l   = 2
        dl  *= -1
    elif l<.5:
        l   = .5
        dl  *= -1
    l   += dl

while 1:
    rate (19)
    proses()
    
    
    
    

    

.

Piston Free Expansion Oscillation.

The result of daydreaming, :)

#code
from visual import *
from random import uniform,random
from visual.controls import *

display(center=(0,0,0),background=(1,1,1), autoscale=False, range=(5,7,7),        
        width=600, height=600,  forward=(0.4,-0.3,-1)) #arah kamera

distant_light(direction=(1,1,1), color=color.red)
cylinder(pos=(-4,0,0), axis=(7,0,0), radius=1, color=color.green, opacity=.1)
piston  = cylinder(pos=(-2,0,0), axis=(1,0,0), radius=.7, color=color.red)
piston1  = cylinder(pos=(-2,0,0), axis=(5,0,0), radius=.1, color=color.red)

p0  = 1.
A   = 1.
x0  = 5.        #panjang piston
c   = p0*A*x0
x   = 3
dx  = x0-x       #
m   = 1.
v   = 0.
dt  = 1./8.
def proses():
    global x,v,dx
    p   = c/x
    dp  = p-p0
    print x,dp
    f   = -dp
    a   = f/m
    v   += a*dt
    dx  += v*dt
    x   -= v*dt

    piston.x    = dx
    piston1.x    = dx
while 1:
    rate (19)
    proses()
    
    
    
    

   

.

Ball Oscillation Inside a Sphere

 with Visual Python.

 I use nonlinear solution, :)

#code
from visual import *
from random import uniform,random
from visual.controls import *

display(center=(0,0,0),background=(1,1,1), #autoscale=False,        
        width=600, height=600,  forward=(-0.4,-0.3,-1)) #arah kamera

distant_light(direction=(1,1,1), color=color.red)

dl  = .01
l   = 1.1

box(color=color.white, pos=(0,0,0),length=2*l,height=dl, width=dl, opacity=.3)
box(color=color.white, pos=(0,0,0),length=dl,height=2*l, width=dl, opacity=.3)
box(color=color.white, pos=(0,0,0),length=dl,height=dl, width=2*l, opacity=.3)


bola1 = sphere (pos=(0,0,0), radius=1.1, color=color.green, opacity=.3)
bola = sphere (pos=(1,0,0), radius=.1, color=color.red)

r   = 1.
x   = 1.
m   = 1.
g   = 1.
vx  = 0.
dt  = 1./16.
def proses():
    global x,vx
    f   = -m*g*x/r
    a   = f/m
    vx  += a*dt
    x   += vx*dt
    r2  = r*r
    x2  = x*x
    if x2>r2:
        x2 = r2
    y   = r-sqrt(r2-x2)
    bola.x  = x
    print x
    bola.y  = y-r
    
while 1:
    rate (19)
    proses()
    
    
    
    

    

.

3D Ball Collision in Python with Visual Module.

I only compute the collision between bed and blue one. They're bola[1] and bola[2] respectively

from visual import *
from random import uniform,random
from visual.controls import *

def change(): 
    global jalan
    if b.value:
        jalan  = True
        for ball in(bola):
            ball.vx  = uniform(-7,7)
            ball.vy  = uniform(-7,7)
            ball.vz  = uniform(-7,7)
    else:
        jalan   = False
        
c = controls(title='Tempat Tombol',x=800, y=0, width=300, height=300, range=50) 
b = toggle( pos=(0,0), width=20, height=20, text='Click me', action=lambda: change() )
display(center=(0,0,0),background=(1,1,1), #autoscale=False,        
        width=600, height=600,  forward=(-0.4,-0.3,-1)) #arah kamera
g   = -1.
dt  = .1
e   = 1.
b.value = True
jalan   = True
l       = 17.
dl      = .01
n       = 11

distant_light(direction=(1,1,1), color=color.red)
lantai      = box(color=color.white, pos=(0,0,0),length=l,height=dl, width=l, opacity=.3)
dindingKiri = box(color=color.white, pos=(-l/2,l/2,0),length=dl,height=l, width=l, opacity=.3)
dindingKanan= box(color=color.white, pos=(l/2,l/2,0),length=dl,height=l, width=l, opacity=.3)
dBelakang   = box(color=color.white, pos=(0,l/2,-l/2),length=l,height=l, width=dl, opacity=.3)
atap        = box(color=color.white, pos=(0,l,0),length=l,height=dl, width=l, opacity=.3)
bola = [] 
for i in arange(n):
    ball = sphere (pos=(uniform(1,7),0,uniform(-7,7)), radius=.3, color=color.green)
    ball.v  = vector(uniform(-7,7),uniform(-7,7),uniform(-7,7))
    bola.append(ball)


bola[1].radius = 2.5
bola[2].radius = 2.5

bola[1].color = color.red
bola[2].color = color.blue
def proses():
    for ball in (bola):
        a       = g
        ball.v[1]      += a*dt
        ball.pos+= ball.v*dt
    tumbukan()
    tumbukanBola()
    
def tumbukan():
    for ball in(bola):
        if ball.y<0:
            ball.y  = 0.01
            ball.v[1]  *=-1.*e
        elif ball.y>l:
            ball.y  = l-.01
            ball.v[1]  *= -1
        if ball.x<-l/2:
            ball.x=-l/2+.01
            ball.v[0]  *= -1*e
        if ball.x>l/2:
            ball.x=l/2-.01
            ball.v[0]  *= -1*e
        if ball.z<-l/2:
            ball.z=-l/2+.01
            ball.v[2]  *= -1*e
        if ball.z>l/2:
            ball.z=l/2-.01
            ball.v[2]  *= -1*e

def tumbukanBola():
    #pass
    jarak = mag(bola[2].pos-bola[1].pos)
    if jarak<(bola[1].radius+bola[2].radius):
        arah    = norm(bola[2].pos-bola[1].pos)
        v1      = dot(bola[1].v,arah)
        v2      = dot(bola[2].v,arah)
        dv      = v2-v1
        bola[1].v   += dv*arah
        bola[2].v   -= dv*arah 



while 1:
    rate (51)
    if jalan:
        proses()
    
    
    
    

    

.

Alright, Now Do It in Python, with style, :)

Python version of this flash action script of electron under Lorentz force, :)

from visual import *
from random import uniform

display(center=(0,0,0),background=(1,1,1), autoscale=False,        
        width=600, height=600,  
        #forward=(-0.4,-0.3,-1)
        ) 

distant_light(direction=(1,1,1), color=color.red)


l       = 11
dt      = 1./8.
medan   = box(color=color.white, 
        pos=(l/2,0,0),length=l,height=l, 
        width=l, opacity=.3)

ball = sphere (pos=(-7,0,0), radius=.3, color=(uniform(0,1),uniform(0,1),uniform(0,1)))
 
def awal():
    global q,m,B,v
    q   = 1.
    m   = 1.
    B   = vector(0.,0.,1.)
    v   = vector(3.,0.,0.)


def proses():
    global v,B
    if ball.x > 0:
        B   = vector(0.,0.,1.)
    else:
        B   = vector(0.,0.,0.)
    print B
    F   = q*(cross(v,B))
    a   = F/m
    v   += a*dt
    
    ball.pos += v*dt
    print ball.x
    
awal()
while 1:
    rate (37)
    proses()
    
    

.

Menyapa Senjata Lama.

Sudah lama tidak otak-atik Macromedia Flash.

Yup, masih yang Macromedia, bukan Adobe, :) .

Membuat animasi bola bermuatan (atau elektron) yang bergerak lurus dengan kecepatan konstan tiba-tiba mencapai daerah dengan  medan magnet. Sesuai hukum Lorentz maka bola/elektron akan bergerak melengkung.

Action script hanya ditulis pada action di frame pertama layer background.

px0 = bola._x;
py0 = bola._y;
awal();
_root.onEnterFrame = function() {
 if (jalan == true) {
  proses();
 }
 //trace(jalan)                         
};
function awal() {
 jalan = false;
 bola._x = px0;
 bola._y = py0;
 q = 1;
 B = 1;
 tB.text = B;
 //massa
 m = 1;
 tm.text = m;
 //kecepatan
 vx = 50;
 vy = -10;
 tvx.text = vx;
 tvy.text = vy;
 //percepatan
 ax = 0;
 ay = 0;
 //posisi
 px = 0;
 py = 0;
 dt = 1/8;
}
function bacaInput() {
 vx = Number(tvx.text);
 vy = Number(tvy.text);
 m = Number(tm.text);
 B = Number(tB.text);
}
function updateNilai() {
 tB.text = B;
 tm.text = m;
 tvx.text = vx;
 tvy.text = vy;
}
function proses() {
 v = Math.sqrt(vx*vx+vy*vy);
 if (bola._x>200) {
  F = B*q*v;
 } else {
  F = 0;
 }
 a = F/m;
 //arah vektor normal v
 nvx = vx/v;
 nvy = vy/v;
 /*arah vektor percepatan
 karena vektor yang tegak lurus a=(ax,ay) adalah at=(-ay,ax)
 */
 nax = -nvy;
 nay = nvx;
 ax = a*nax;
 ay = a*nay;
 //hitung kecepatan baru
 vx += ax*dt;
 vy += ay*dt;
 px += vx*dt;
 py -= vy*dt;
 trace(F);
 //update posisi bola
 bola._x = px0+px;
 bola._y = py0+py;
 updateNilai();
}
//tombol-tombol
//tombol tbJalan
tbJalan.onRelease = function() {
 bacaInput();
 jalan = true;
};
tbStop.onRelease = function() {
 jalan = false;
};
tbReset.onRelease = function() {
 awal();
};

.

Iterasi di Python

Jika kita punya sebuah list bernama bola yang didefinisikan sebagai:

bola = []

kemudian kita membuat obyek bernama ball

ball = sphere (pos=(uniform(1,7),0,uniform(-7,7)), radius=.3, color (uniform(0,1),uniform(0,1),uniform(0,1)))
ball.vx  = uniform(-7,7)
ball.vy  = uniform(-7,7)
ball.vz  = uniform(-7,7)
 
Obyek ini kita masukkan ke dalam bola dengan perintah

bola.append(ball)

Kita dapat melakukannya berkali-kali sehingga pada list bola terdapat beberapa obyek bernama ball

Jika kita ingin mengakses obyek tersebut, kita dapat menggunakan perintah semacam

bola[0].vx = 1

jika kita ingin mengakses vx di semua obyek bola, kita dapat menggunakan iterasi

for i in arange (n):
  bola[i].vx = 1


Namun di python ada cara lain yang juga mudah

for ball in (bola):
  ball.vx = 1

Berikut contoh kode yang menggunakan iterasi seperti itu

from visual import *
from random import uniform,random
from visual.controls import *

def change(): 
    global jalan,vx,vy,vz
    if b.value:
        jalan  = True
        for ball in(bola):
            ball.vx  = uniform(-7,7)
            ball.vy  = uniform(-7,7)
            ball.vz  = uniform(-7,7)
    else:
        jalan   = False
        
c = controls(title='Tempat Tombol',x=800, y=0, width=300, height=300, range=50) 
b = toggle( pos=(0,0), width=20, height=20, text='Click me', action=lambda: change() )
display(center=(0,0,0),background=(1,1,1), #autoscale=False,        
        width=600, height=600,  forward=(-0.4,-0.3,-1)) #arah kamera
g   = -1.
dt  = .1
e   = 1.
b.value = True
jalan   = True
l       = 17.
dl      = .01
n       = 11

distant_light(direction=(1,1,1), color=color.red)
lantai      = box(color=color.white, pos=(0,0,0),length=l,height=dl, width=l, opacity=.3)
dindingKiri = box(color=color.white, pos=(-l/2,l/2,0),length=dl,height=l, width=l, opacity=.3)
dindingKanan= box(color=color.white, pos=(l/2,l/2,0),length=dl,height=l, width=l, opacity=.3)
dBelakang   = box(color=color.white, pos=(0,l/2,-l/2),length=l,height=l, width=dl, opacity=.3)
atap        = box(color=color.white, pos=(0,l,0),length=l,height=dl, width=l, opacity=.3)
bola = [] 
for i in arange(n):
    ball = sphere (pos=(uniform(1,7),0,uniform(-7,7)), radius=.3, color=(uniform(0,1),uniform(0,1),uniform(0,1)))
    ball.vx  = uniform(-7,7)
    ball.vy  = uniform(-7,7)
    ball.vz  = uniform(-7,7)
    bola.append(ball)

def proses():
    global vx,vy,vz
    for ball in (bola):
        a       = g
        ball.vy      += a*dt
        ball.pos+= vector(ball.vx*dt,ball.vy*dt,ball.vz*dt)
    tumbukan()
    
def tumbukan():
    global vx,vy,vz
    for ball in(bola):
        if ball.y<0:
            ball.y  = 0.01
            ball.vy  *=-1.*e
        elif ball.y>l:
            ball.y  = l-.01
            ball.vy  *= -1
        if ball.x<-l/2:
            ball.x=-l/2+.01
            ball.vx  *= -1*e
        if ball.x>l/2:
            ball.x=l/2-.01
            ball.vx  *= -1*e
        if ball.z<-l/2:
            ball.z=-l/2+.01
            ball.vz  *= -1*e
        if ball.z>l/2:
            ball.z=l/2-.01
            ball.vz  *= -1*e

while 1:
    rate (100)
    if jalan:
        proses()
    
    

.

Lorenz Attractor in Python with Visual Module.

 So much faster than matplotlib 3d projection

from visual import *
from random import uniform
display(center=(0,0,0),         #pusat display
        background=(1,1,1), 
        #autoscale=False,        #agar display tidak otomatis mengikuti obyek 
        width=600, 
        height=600,
        forward=(-0.4,-0.3,-1)) #arah kamera
x   = 1.
y   = 1.
z   = 1.

dt  = 1./64.
s   = 10.
b   = 8./3.
r   = 28
w   = 0.
dw  = 0.01
    
n   = 0
while 1:
    rate(1)
    while n<3000:
        n   += 1

        xdot    = s * (y-x)
        ydot    = x*r -x*z -y
        zdot    = x*y -b*z

        x       = x+xdot*dt
        y       = y+ydot*dt
        z       = z+zdot*dt

        w       +=dw
        if w>1 or w<0:
            dw   = -dw

        sphere(pos=(x,y,z),radius=.7,color=(0.,w,0))

Forced to get Force using Center of Mass

 I know. Using solely center of mass as main contributor in force calculation is bad idea in nearly homogenous n-body system.

 Well, I'm curious about that. So, just bite it.

 And, yeah, bad idea.

 The system is fine, generally.

 Each body move as usual, at a glimpse.

 But if we go detail about that. There's strange behavior here and there.


 On previous code with brute force calculation, we could saw, binary system consist of two body orbiting each other. But, with this center-of-mass-code, we see no one. A body that have very close neighbor don't make a twin system. It continues orbiting center of mass.

 As we have system with -body with same mass, it should be natural if two close-body will have greater gravity force than with any other body. It should have "priority"; orbiting each other.

 Ok, this "forced"-center-mass-code is weird, but it'll useful in other system. :)

from visual import *
from random import uniform,random

l       = 17.
dl      = .01
display(center=(0,0,0),background=(1,1,1),autoscale=False, width=600, height=600,forward=(-0.7,-0.7,-1))
distant_light(direction=(1,1,1), color=color.red)
#local_light(pos=(0,0,0), color=(0,0,1))
#sphere(pos=(0,0,0), color=(0,0,0), material=materials.emissive, opacity=0.9)

box(color=color.white, pos=(0,l/2,0),length=l,height=dl, width=l, opacity=0.3)
box(color=color.white, pos=(0,-l/2,0),length=l,height=dl, width=l, opacity=0.3)
box(color=color.white, pos=(l/2,0,0),length=dl,height=l, width=l, opacity=0.3)
box(color=color.white, pos=(-l/2,0,0),length=dl,height=l, width=l, opacity=0.3)
box(color=color.white, pos=(0,0,l/2),length=l,height=l, width=dl, opacity=0.3)

bola    = [] 
n       = 11
dv      = .1        #kecepatan setelah menabrak dinding
G       = 31.
for i in range(n):
    ball            = sphere (pos=(uniform(-7,7),uniform(-7,7),uniform(-7,7)), radius=.2, color=(random(),random(),random()))
    ball.v          = vector(uniform(-1,1),uniform(-1,1),uniform(-1,1))
    
    bola.append(ball)


dt = 1./32.
bola[0].pos=(1,1,1)


def hitungGaya(i,rcm):
    r           = bola[i].pos-rcm
    jarak       = mag(r)
    arah        = -norm(r)
    gaya        = G*1./(pow(jarak,2)+.1)*arah*n
    return gaya
def hitungPusatMassa():
    xcm = 0.
    ycm = 0.
    zcm = 0.
    for i in arange(n):
        xcm += bola[i].x
        ycm += bola[i].y
        zcm += bola[i].z
    xcm /= n
    ycm /= n
    zcm /= n
    return vector(xcm,ycm,zcm)
        
def proses():
    global bola
    rcm = hitungPusatMassa()
    for i in arange(n):
        gaya        = hitungGaya(i,rcm)
        r           = bola[i].pos
        #tambah gaya dari pusat
        '''
        jarak       = mag(bola[i].pos)
        arah        = -norm(bola[i].pos)
        gaya        += 10*G*1./(pow(jarak,2)+.05)*arah
        '''
        a           = gaya
        v           = bola[i].v
        v           += a*dt
        r           += v*dt
        #cek pantul ke tembok
        if bola[i].x > l/2.:
            bola[i].x   = l/2.
            v.x         *= -dv
        elif bola[i].x < -l/2.:
            bola[i].x   = -l/2.
            v.x         *= -dv
        elif bola[i].y > l/2.:
            bola[i].y   = l/2.
            v.y         *= -dv
        if bola[i].y < -l/2.:
            bola[i].y   = -l/2.
            v.y         *= -dv
        if bola[i].z > l/2.:
            bola[i].z   = l/2.
            v.z         *= -dv
        if bola[i].z < -l/2.:
            bola[i].z   = -l/2.
            v.z         *= -dv

        bola[i].pos = r
        
while 1:
    rate (10)
    proses()

.

Menjinakkan Gravitasi.

 Cari cara agar mereka jadi enak untuk dilihat.

 Di kode sebelumnya, bola sudah dicegah agar tidak melarikan diri dengan menambahkan kotak, mereka akan memantul jika menabrak kotak.

 Beres, sekarang tambahkan interaksi antar bola, interaksi sederhana, gravitasi.

 Ok, pake metode brutal, untuk tiap-tiap bola, hitung satu-satu interaksi gravitasi dengan seluruh bola, kemudian jumlahkan.





(ide, cari pusat masa seluruh bola, perlakukan interaksi gaya seperti kode sebelumnya)

 (saat nulis tadi, sambil mikir, ternyata itu ide buruk, tapi boleh juga, lihat saja nanti, bagus atau tidak)

 (tapi tidak sekarang, tidak masuk kode yang sekarang)

 (kode di posting ini tetap pake metode brutal)

 Hasilnya adalah, ..., ok, mereka tetap memantul saat menabrak tembok, tetapi di dalam kotak mereka bergerak dengan sangat cepat, bola liar.

 (yep, metode brutal menghasilkan hasil brutal, tetapi tidak selalu...)

 Ide licik, tambahkan koefisien restitusi di dinding, tumbukan dengan dinding tidak lenting sempurna, kecepatannya tinggal 10% setelah memantul

 Dengan demikian mereka lebih jinak di dalam kotak, :)

from visual import *
from random import uniform,random

l       = 17.
dl      = .01
display(center=(0,0,0),background=(1,1,1),autoscale=False, width=600, height=600,forward=(-0.7,-0.7,-1))
distant_light(direction=(1,1,1), color=color.red)
#local_light(pos=(0,0,0), color=(0,0,1))
#sphere(pos=(0,0,0), color=(0,0,0), material=materials.emissive, opacity=0.9)

box(color=color.white, pos=(0,l/2,0),length=l,height=dl, width=l, opacity=0.3)
box(color=color.white, pos=(0,-l/2,0),length=l,height=dl, width=l, opacity=0.3)
box(color=color.white, pos=(l/2,0,0),length=dl,height=l, width=l, opacity=0.3)
box(color=color.white, pos=(-l/2,0,0),length=dl,height=l, width=l, opacity=0.3)
box(color=color.white, pos=(0,0,l/2),length=l,height=l, width=dl, opacity=0.3)

bola    = [] 
n       = 11
dv      = .1        #kecepatan setelah menabrak dinding
G       = 31.
for i in range(n):
    ball            = sphere (pos=(uniform(-7,7),uniform(-7,7),uniform(-7,7)), radius=.3, color=(random(),random(),random()))
    ball.v          = vector(uniform(-1,1),uniform(-1,1),uniform(-1,1))
    
    bola.append(ball)


dt = 1./32.
bola[0].pos=(1,1,1)


def hitungGaya(i):
    gaya    = vector(0,0,0)
    for j in arange(n):
        if i!=j:  
            r           = bola[i].pos-bola[j].pos
            jarak       = mag(r)
            arah        = -norm(r)
            gaya        += G*1./(pow(jarak,2)+.1)*arah
    return gaya
def proses():
    global bola
    for i in arange(n):
        gaya        = hitungGaya(i)
        #tambah gaya dari pusat
        r           = bola[i].pos
        jarak       = mag(bola[i].pos)
        arah        = -norm(bola[i].pos)
        gaya        += 10*G*1./(pow(jarak,2)+.05)*arah

        a           = gaya
        v           = bola[i].v
        v           += a*dt
        r           += v*dt
        #cek pantul ke tembok
        if bola[i].x > l/2.:
            bola[i].x   = l/2.
            v.x         *= -dv
        elif bola[i].x < -l/2.:
            bola[i].x   = -l/2.
            v.x         *= -dv
        elif bola[i].y > l/2.:
            bola[i].y   = l/2.
            v.y         *= -dv
        if bola[i].y < -l/2.:
            bola[i].y   = -l/2.
            v.y         *= -dv
        if bola[i].z > l/2.:
            bola[i].z   = l/2.
            v.z         *= -dv
        if bola[i].z < -l/2.:
            bola[i].z   = -l/2.
            v.z         *= -dv

        bola[i].pos = r
        
while 1:
    rate (10)
    proses()

.

Bermain dengan Gaya Sentral di Python.

 Merangkak dari kode sebelumnya, dengan tambahan kotak yang mengungkung bola-bola yang mengitari titik asal dengan berbagai kecepatan. Bola-bola tersebut tak bisa menembus kotak, jika menabrak salah satu dinding kotak, maka dia akan memantul

from visual import *
from random import uniform,random

l       = 17.
dl      = .01
display(center=(0,0,0),background=(1,1,1),autoscale=False, width=600, height=600,forward=(-0.7,-0.7,-1))
distant_light(direction=(1,1,1), color=color.red)

box(color=color.white, pos=(0,l/2,0),length=l,height=dl, width=l, opacity=0.3)
box(color=color.white, pos=(0,-l/2,0),length=l,height=dl, width=l, opacity=0.3)
box(color=color.white, pos=(l/2,0,0),length=dl,height=l, width=l, opacity=0.3)
box(color=color.white, pos=(-l/2,0,0),length=dl,height=l, width=l, opacity=0.3)
box(color=color.white, pos=(0,0,l/2),length=l,height=l, width=dl, opacity=0.3)

bola    = [] 
n       =   11

G       = 31.
for i in range(n):
    ball            = sphere (pos=(uniform(1,7),0,uniform(-7,7)), radius=.3, color=(random(),random(),random()))
    ball.v          = vector(0,uniform(1,3),0)
    
    bola.append(ball)


dt = 1./32.
bola[0].pos=(1,1,1)
while 1:
    rate (100)
    for i in arange(n):
        r           = bola[i].pos
        v           = bola[i].v
        jarak       = mag(bola[i].pos)
        arah        = -norm(bola[i].pos)
        gaya        = G*1./(pow(jarak,2)+.1)*arah
        a           = gaya
        v           += a*dt
        r           += v*dt
        #cek pantul ke tembok
        if bola[i].x > l/2.:
            bola[i].x   = l/2.
            v.x         *= -1
        elif bola[i].x < -l/2.:
            bola[i].x   = -l/2.
            v.x         *= -1
        elif bola[i].y > l/2.:
            bola[i].y   = l/2.
            v.y         *= -1
        if bola[i].y < -l/2.:
            bola[i].y   = -l/2.
            v.y         *= -1
        if bola[i].z > l/2.:
            bola[i].z   = l/2.
            v.z         *= -1
        if bola[i].z < -l/2.:
            bola[i].z   = -l/2.
            v.z         *= -1

        bola[i].pos = r
        

.

Playing with Central Force in Pyton with Visual Module

 I don't use n body calculation, every sphere only attracted by force from origin, dependent upon its position.

from visual import *
from random import uniform,random

l       = 17.
dl      = .01
display(center=(0,0,0),background=(1,1,1),autoscale=False, width=600, height=600,forward=(-0.4,-0.3,-1))
distant_light(direction=(1,1,1), color=color.red)
box(color=color.white, pos=(0,0,0),length=l,height=dl, width=l, opacity=0.3)
box(color=color.white, pos=(0,0,0),length=dl,height=l, width=l, opacity=0.3)
box(color=color.white, pos=(0,0,0),length=l,height=l, width=dl, opacity=0.3)

bola    = [] 
n       =   11

G       = 31.
for i in range(n):
    ball            = sphere (pos=(uniform(1,7),0,uniform(-7,7)), radius=.3, color=(random(),random(),random()))
    ball.v          = vector(0,uniform(1,3),0)
    
    bola.append(ball)


dt = 1./32.
bola[0].pos=(1,1,1)
while 1:
    rate (100)
    for i in arange(n):
        r           = bola[i].pos
        v           = bola[i].v
        jarak       = mag(bola[i].pos)
        arah        = -norm(bola[i].pos)
        gaya        = G*1./(pow(jarak,2)+.1)*arah
        a           = gaya
        v           += a*dt
        r           += v*dt

        bola[i].pos = r
        

.

Wave Simulation Using Visual Python.

 I used previous code (with matplotlib animation).

 Removed the matplotlib part, swap it with vpython, :)

"""
Gelombang
"""
import numpy as np                          #untuk operasi array
from visual import *

#variabel
n   = 39
x   = np.arange(0., 1.,1./n )   #array dari 0 s.d 1 berjarak 1/n
y   = np.zeros(n)               #array sejumlah n isinya 0
y1  = np.zeros(n)               # y1 : array [0..n] of real
y2  = np.zeros(n)
y1  = np.exp(-1*np.power(10*x-3,2))
y2  = np.exp(-1*np.power(10*x-3,2))

r2  = 1./512   

display(center=(.5,0,0),background=(1,1,1))
#kotak
dindingKiri   = box (pos=(0.,0,0), length=.01, height=1, width=1, color=color.green)
dindingKanan   = box (pos=(1.,0,0), length=.01, height=1, width=1, color=color.blue)
bola    = []
for i in range(n):
    ball            = sphere (pos=(x[i],y2[i],0), radius=.01, color=color.red)
    bola.append(ball)

#nilai awal fungsi gaussian

#print y2


def proses():
    #hitung nilai baru
    for i in np.arange(1,n-1):
        y[i]    =   2*(1-r2)*y1[i]-y2[i]+r2*(y1[i+1]+y1[i-1])
    #geser
    y2[:]   =   y1[:]
    y1[:]   =   y[:]
    #print y2

    return y2

while 1:
    rate(100)
    proses()
    for i in np.arange(n):
        bola[i].y   = y2[i]

.

                          

Hello (Again) Visual Python, :D .

Menyapa kembali mainan lama, :)

from visual import *

floor   = box (pos=(0,0,0), length=4, height=0.5, width=4, color=color.blue)
bola    = [] 
n   =   2
for i in range(n):
    ball            = sphere (pos=(0,4,0), radius=1, color=color.red)
    ball.v          = vector(0,-1,0)

    bola.append(ball)

dt = 0.01
bola[1].pos=(1,1,1)
while 1:
    rate (100)
    bola[0].pos = bola[0].pos + bola[0].v*dt
    if bola[0].y < bola[0].radius:
        bola[0].v.y = abs(bola[0].v.y)
    else:
        bola[0].v.y = bola[0].v.y - 9.8*dt

.

 

Lorenz Attractor 3D Scatter Plot using Python with Matplotlib

 After failed with regular plot() syntax, I have luck with the scatter() one, :)

"""
Cluster
"""
import numpy as np                          #untuk operasi array
import matplotlib.pyplot as plt             #untuk gambar grafik
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.animation as animation    #untuk menggerakkan grafik
from itertools import cycle
randColor = cycle('bgrcmk').next
randMarker= cycle('.,ov^<>12348sp*hH+xDd|_').next

#fig, ax = plt.subplots()
fig = plt.figure()
ax  = fig.add_subplot(111,projection='3d')

x   = 1.
y   = 1.
z   = 1.

ax.scatter(x, y, z )

def animate(i):
    global x,y,z,n
    
    dt  = 1./64.
    s   = 10.
    b   = 8./3.
    r   = 28
    
    
    
    xdot    = s * (y-x)
    ydot    = x*r -x*z -y
    zdot    = x*y -b*z

    x       = x+xdot*dt
    y       = y+ydot*dt
    z       = z+zdot*dt


    ax.scatter(x,y,z, c=randColor(), marker=randMarker())

ani = animation.FuncAnimation(fig, animate, frames=2000, interval=10, blit=False)
#ani.save('LorenzAtrractor3D.mp4',bitrate=1024)
plt.show()

.

 

Lorenz Attractor in Python with Matplotlib.

 Using three mapping.

import numpy as np                          #untuk operasi array
import matplotlib.pyplot as plt             #untuk gambar grafik
import matplotlib.animation as animation    #untuk menggerakkan grafik

fig, (ax, ay, az) = plt.subplots(3,sharex=True)

x   = 1.
y   = 1.
z   = 1.
#plt.ylim(-43,43)
#plt.xlim(-43,43)

#membuat garis/kurva dengan sumbu-x adalah x, sumbu-y adalah y
line, = ax.plot(x, y, 'o' )
line, = ay.plot(x, z, 'o' )
line, = az.plot(y, z, 'o' )

def animate(i):
    global x,y,z
    
    
    dt  = 1./64.
    s   = 10.
    b   = 8./3.
    r   = 28
    
    xdot    = s * (y-x)
    ydot    = x*r -x*z -y
    zdot    = x*y -b*z

    x       = x+xdot*dt
    y       = y+ydot*dt
    z       = z+zdot*dt

    plt.figure(1)

    line, = ax.plot(x, y, 'o' )
    line, = ay.plot(x, z, 'o' )
    line, = az.plot(y, z, 'o' )

    return line,

ani = animation.FuncAnimation(fig, animate, frames=2000, interval=100, blit=False)
#ani.save('Lorenz.mp4',bitrate=1024)
plt.show()

.

Game of Life.

I use Conway's model in Python with numpy and matplotlib module.

"""
Cluster
"""
import numpy as np                          #untuk operasi array
import matplotlib.pyplot as plt             #untuk gambar grafik
import matplotlib.animation as animation    #untuk menggerakkan grafik

fig, ax = plt.subplots()

n   = 19

plt.ylim(0,n)
plt.xlim(0,n)
a   = np.zeros((n,n))
a0  = np.zeros((n,n))
#buat nilai awal
for i in np.arange(n):
    for j in np.arange(n):
        r   = np.random.randint(100)
        if r<50:
            a0[i,j] = 1
            line,   = ax.plot(i,j,'o')
        else:
            a0[i,j] = 0
            line,   = ax.plot(i,j,'wo')

a[:,:]=a0[:,:]


def animate(i):
    global line
    for i in np.arange(1,n-1):
        for j in np.arange(1,n-1):
            #hitung tetangga
            t = a0[i-1,j-1]+a0[i-1,j]+a0[i-1,j+1]+\
                    a0[i,j-1]+a0[i,j+1]+\
                    a0[i+1,j-1]+a0[i+1,j]+a0[i+1,j+1]
            #hidup atau mati?
            if a0[i,j]==1:
                if t<2 or t>3:
                    a[i,j]  = 0

            else:
                if t==3:
                    a[i,j]  =1
            #gambar
            if a[i,j]==1:
                line,   = ax.plot(i,j,'o')
            else:
                line,   = ax.plot(i,j,'wo')
    a0[:,:] = a[:,:]
    return line,

ani = animation.FuncAnimation(fig, animate, frames=2000, interval=100, blit=False)
#ani.save('cluster.mp4',bitrate=1024)
plt.show()

.

Interference

Try to make a simple code

#mentah

from pylab import *

n   = 193
t   = 0
dy  = 1./64.

w   = 1./8.         #wavelength

L   = 1.            #distance to screen

d   = .001/L          #half slit distance, weird it is, but anyway...
y   = zeros(n)
S   = zeros(n)
A   = 1.
f   = 1./w          # v = 1. ,:)


for i in arange(n):
    y[i]= dy*i
    l1  = sqrt(pow((y[i]-d),2)+pow(L,2))
    l2  = sqrt(pow((y[i]+d),2)+pow(L,2))
    print '***'
    print l1
    print l2
    y1  = A*sin(f*l1)
    y2  = A*sin(f*l2)
    #S[i]= y1+y2
    S[i]= pow((y1+y2),2)   


plot(y,S)
xlabel('y')
ylabel('V')
title('Interferensi')
grid(True)

show()

.

Python Turtle Assymetrical Branch

Ok, it's my last fractal this month, :D

import turtle
import numpy

#buat pola di sini
#kura-kura menghadap ke atas
turtle.shape("turtle")
turtle.speed(10)
turtle.left(90)


lv  = 11
l   = 100
dl  = 3./4.
sl  = 17
sr  = 39
bl  = 1./2.
br  =1./4.

turtle.penup()
turtle.backward(l)
turtle.pendown()
turtle.forward(l)


def maju(l,level):
    level   += 1
    turtle.backward(l*bl)
    l       = l*dl
    turtle.left(sl)
    turtle.forward(l)
    if level<=lv:
        maju(l,level)
    #mundur, tengok kanan
    turtle.backward(l)
    turtle.right(sl)
    turtle.forward(l*bl/dl)
    
    turtle.backward(l*br/dl)
    turtle.right(sr)
    turtle.forward(l)
    if level<lv:
        maju(l,level)
    
    
    turtle.backward(l)
    turtle.left(sr) 
    turtle.forward(l*br/dl)
    
    level   -= 1
    

maju(l,2)

#agar gambar tak langsung hilang
turtle.exitonclick()



 




 

.

Pohon Asimetris, :)

Modifikasi dari kode sebelumnya

import turtle
import numpy

#buat pola di sini
#kura-kura menghadap ke atas
turtle.shape("turtle")
turtle.left(90)


lv  = 11
l   = 100
sl  = 47
sr  = 17

turtle.penup()
turtle.backward(l)
turtle.pendown()
turtle.forward(l)


def maju(l,level):
    l           = 3./4.*l
    #turtle.backward(l)
    turtle.left(sl)
    turtle.forward(l)
    level       += 1
    if level<lv:
        maju(l,level)
    
    turtle.backward(l)
    turtle.right(sl)
    turtle.right(sr)
    turtle.forward(l)
    if level    <lv:
        maju(l,level)
    turtle.backward(l)
    turtle.left(sr)
    level       -= 1

maju(l,2)


#agar gambar tak langsung hilang
turtle.exitonclick()



 

.