Nugroho's blog.

Saturday, April 16, 2016

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()


.



Friday, April 15, 2016

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


.


Thursday, April 14, 2016

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]


.



                          

Wednesday, April 13, 2016

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

.

 


323f (5) amp (1) android (12) apple (7) arduino (18) art (1) assembler (21) astina (4) ATTiny (23) blackberry (4) camera (3) canon (2) cerita (2) computer (106) crazyness (11) debian (1) delphi (39) diary (286) flash (8) fortran (6) freebsd (6) google apps script (8) guitar (2) HTML5 (10) IFTTT (7) Instagram (7) internet (12) iOS (5) iPad (6) iPhone (5) java (1) javascript (1) keynote (2) LaTeX (6) lazarus (1) linux (29) lion (15) mac (28) macbook air (8) macbook pro (3) macOS (1) Math (3) mathematica (1) maverick (6) mazda (4) microcontroler (35) mountain lion (2) music (37) netbook (1) nugnux (6) os x (36) php (1) Physicist (29) Picture (3) programming (189) Python (109) S2 (13) software (7) Soliloquy (125) Ubuntu (5) unix (4) Video (8) wayang (3) yosemite (3)