Spel7: Planetsystem

newton-kimdir                                      Om Newton och Newtons ekvationer.

Vi skriver datorprogram för lösning av Newton’s rörelseekvationer, med användning av 2d vektorn vec2() och tabeller för position, hastighet och acceleration

  • pos = pos +hast*dt,
  • hast = vel + acc*dt
  • acc = kraft/massa

där kraften mellan två kroppar av massa1 och massa2 med position pos1 och pos2 med avstånd r, är given av

  • kraft = (pos1 -pos2)*massa1*massa2/r^3.

Vi simulerar olika planetsystem med olika startvärden för position och hastighet  och då speciellt vårt eget med sol-jord-måne + övriga planeter. Jfr med detta. Jämför med Mechanics1 på App Store.

 

Template 1: (Följ planeters rörelse under ömsesidig gravitation och här i 3d)

function setup()
I=20

pos = {}
vel = {}
acc= {}
force = {}
m={}
for i=1,I do
pos[i] = vec2(0,0)
vel[i] = vec2(0,0)
acc[i] = vec2(0,0)
force[i] = vec2(0,0)
m[i]= 1/i
end

print(“Planet”)
print(“force~1/r^2”)
print(“acc=force/mass”)
print(“vel=vel+acc*dt”)
print(“pos=pos+vel*dt”)
dt =0.001
for i=1,I do
pos[i]=vec2(math.random(200,400),math.random(300,600))
vel[i]=vec2(math.random(-100,100),math.random(-100,100))
m[i] = math.random(1,100)
end

end

function draw()

font(“AmericanTypewriter-Condensed”)
fontSize(30)
fill(0, 22, 255, 255)
text(“Motion-Force: Newton’s 2nd Law”,250, 800)
text(“force = mass x acceleration, f = m x a”,250, 750)
text(“acc=force/mass”,250, 700)
text(“vel=vel+acc*dt”,250, 650)
text(“pos=pos+vel*dt”,250, 600)
stroke(43, 203, 22, 108)
strokeWidth(1)
for i=1,I do
force[i]=vec2(0,0)
for j=1,I do
r= (pos[i]-pos[j]):len()
force[i]=force[i]-(1000000/(r+10)^3)*(pos[i]-pos[j])
end
acc[i] = force[i]
vel[i]=vel[i]+acc[i]*dt
pos[i]=pos[i]+vel[i]*dt
end
fill(42, 0, 255, 255)
for i=1,I do
ellipse(pos[i].x,pos[i].y,10,10)
end
end

Spel8: Kvadratroten ur 2

sqrt2

Försök beräkna kvadratroten ur 2, noterad som \sqrt{2}, dvs försök bestämma ett (positivt) tal x  sådant att

  • x*x=2.

Prova:

  • x=1 ger x*x=1 < 2, dvs x=1 är för litet
  • x=2 ger x*x = 4 > 2, dvs x=2 är för stort
  • x=1.5 ger x*x=2.25, dvs x=2.25 är för stort
  • x=1.4 ger x*x=1.96, dvs x=1.4 är för litet,
  • x=1.45 ger x*x=2.1025, dvs 1.45 är för stort,
  • fortsätt att prova med hjälp av miniräknare och bestäm fler decimaler!

Skriv sedan ett datorprogram som gör denna provning och bestäm så många decimaler som möjligt.

Prova sedan med iterationen x = 1/x+ x/2 med start x=1.  Observera att x*x=2 kan skrivas som 2*x*x = 2+x*x, dvs efter division med 2*x, som  x=1/x +x/2.

Utvidga till lösning av ekvationen x*x= a för godtyckligt a>0, tex via iterationen

  • x=a/(2*x) +x/2.

Utvidga till lösning av x^n=a för n=3,4,5,…, tex via iterationen

  • x=(n-1)*x/n+a/(n*x^(n-1)).

Jämför med Pythagoras på App Store.

Spel9: Derivata och Integral

midpointrule

Programmera tidsstegningen

  • x=x+v*dt eller dx = v*dt

för positionen x med olika val av hastigheten v och tidssteget dt.  Börja med v=1 och dt = 1 och plotta x på skärmen. Välj sedan v = t, v = t^2, v = t^3 osv. Jämför med Möte 4.

Notera att v*dt kan tolkas som ytan av en rektangel med bas dt och höjd v* och att därför x kan tolkas som en summa av rektangelytor som tillsammans bildar ytan under hastighetskurvan.

Notera att tidsstegningen dx = v*dt kan skrivas dx/dt = v där dx/dt kallas derivatan av x med avseende på t, dvs

  • derivatan av positionen med avseende på tiden är lika med hastigheten.

Omvänt säger man att

  • Positionen x är integralen av hastigheten v om dx/dt = v.

Tidsstegningen dv =a*dt med v hastighet och a acceleration kan alltså uttryckas på följande sätt:

  • derivatan av hastigheten med avseende på tiden är lika med accelerationen
  • hastigheten är integralen av accelerationen.

Derivatan dx/dt anger hur snabbt positionen ändras (dx) med ändring av tiden (dt), vilket ju är hastigheten v=dx/dt, och dv/dt anger hur snabbt hastigheten ändras (dv) med tiden (dt), vilket är accelerationen a=dv/dt.

Notera att om dx = v*dt, eller dx/dt = v, så gäller att

  • positionen x är integralen av derivatan dx/dt av x
  • hastigheten v är derivatan av integralen x av v.

Jämför med:

  • barnet V är avkomma (barn) till föräldrarna X
  • föräldrarna X är upphov (föräldrar) till barnet V
  • X = upphov till avkomma till  X
  • V = avkomma till upphov  till V.

Jämför med Calculus2 på App Store.

Template 1: (Följ tidsstegningen för dx=vdt med v=t,  för dt=1 och för dt=0.2)

function setup()

x = 0
v = 0
dt = 1
n=0
print(“Touch för att tidsstega dx=vdt med v=t”)

end

function draw()

fill(255,0,0,0,255)
rect(n*dt*20,200,20*dt,v*5)
ellipse(n*dt*20,200+x*5,10)
fill(0,255,0,255)
ellipse(n*dt*20, 200+v*5,10)
end

function touched(touch)

n=n+1
v=n*dt
x=x+v*dt

end

Spel10: Rävar och Kaniner

predator-prey-model

Låt R vara antalet rävar och K antalet kaniner i ett visst område och antag att R och K förändras enligt tidsstegningen

  • R = R +R*K*dt-10*R*dt
  • K = K – R*K*dt +10*K*dt,

där R=10 och K=20 som utgångsvärden och tidsteget dt=0.002. Skriv ett program som utför tidsstegningen och studera hur R och K ändras med tiden. Vilka lagar för förändringen av R och K uttrycks genom tidsstegningen? Ändra parametrar i modellen och studera utfallet. Utvidga till fler rovdjur och bytesdjur. Utforma och programmera ett datorspel som bygger på räv-kanin modellen. Jämför med Biology1 på App Store.

Template 1 (Titta på jakten)

function setup()
print(“Rävar och Kaniner”)
end

K=20
R=10
i=1
dt=0.002

function draw()
i=i+1
R=R+K*R*dt-10*R*dt
K=K-R*K*dt+10*K*dt
fill(255, 0, 74, 255)
ellipse(i/2+20,15*R+300,10)
fill(0, 255, 9, 255)
ellipse(i/2+20,15*K+300,10)
end

Spel11: Kampen mellan Färgerna

Aspidites_ramsayi-RD-example

Kampen mellan två motstridiga sidor, säga röda och blå, kan simuleras med tidsstegningen:

u[i]=u[i]+A*(u[i+1]-2*u[i]+u[i-1])+u[i]*(1-u[i]^2)  för i=1,…,M,

där u[1], u[2],…u[M] är en lista av värden med utgångsvärden mellan -1 och +1 (tex med u[i]=1-2*math.random()), där u[i]=-1 betyder röd och u[i]=1 betyder blå och u[1]=1 och u[M]=-1.

Programmera tidsstegningen och studera kampen mellan röda och blå. Välj tex M=100 och A=0.1 till att börja med och se vad som händer då Du ändrar A (tex genom att lägga upp den som en parameter i Codea och sedan styra den interaktivt under tidsstegningen). Koppla till mönstret på ormens skinn i bilden. Försök ge en fysikalisk mening till tidsstegningen. Jämför med Biology1 på App Store.

Template 1: (Titta på kampen)

function setup()

uold = {}
unew = {}
um = {}
M =100

for i=1,M do
uold[i]=1-2*math.random()
unew[i]=uold[i]
um[i]=0
end

uold[1]=1
unew[1]=1
uold[M]=-1
unew[M]=-1

print(“Wave Motion 1d”)
parameter.number(“A”,1,10,1,function(n) print(“A set to “..n) end)
parameter.number(“N”,1,200,200,function(n) print(“A set to “..n) end)

end

function draw()

background(237, 237, 241, 255)

k=1/N

for m=1,4 do
for i = 1,M do
um[i]=(uold[i]+unew[i])/2
end
for i =2,M-1 do
D2=A*(um[i+1]-2*um[i]+um[i-1])+(um[i]-um[i]^3)
unew[i] = uold[i]+ k*D2
end
end

for i=1,M do
uold[i]=unew[i]
fill(255*uold[i], 0,255*(1-uold[i]),255)
rect(i*5,500,10,200)
ellipse(i*5,300+100*uold[i],10)
end

end

Spel12: Svängande Strängen

313751.image0

En svängande gitarrsträng, eller ett rep som man skakar i ändarna, kan simuleras med tidsstegningen med tidssteget dt:

  • v[i]=v[i]+dt*(u[i+1]-2*u[i]+u[i-1])/dx^2
  • u[i]=u[i]+dt*v[i]

där u[1], u[2],…u[M] är en lista av vertikala ändringar av position av strängen/repet och v[1], v[2],…,v[M], är en lista av motsvarande hastigheter med vissa utgångsvärden och randvärden för i=1 och i=M givna. Låt tex M=100,  motsvarande en indelning av en sträng av längden 1 i delar av längden dx=1/M.

Programmera tidsstegningen och studera resultaten för olika val av M och dt, utgångsvärden och randvärden. Försök ge tidsstegningen en fysikalisk mening.

Jämför med Mechanics1 på App Store.

Template 1: (Titta på den svängande strängen)

function setup()

print(“Wave Motion 1d”)

uold = {}
unew = {}
vold = {}
vnew = {}
um = {}
M = 50

for i=1,M do
uold[i]=0
unew[i]=0
vold[i]=0
vnew[i]=0
um[i]=0
end

for i=20,30 do
vold[i]=20
end
parameter.number(“F”,-10,10,0,function(n) print(“F set to “..n) end)
parameter.number(“A”,1,10,5,function(n) print(“A set to “..n) end)
parameter.number(“U”,-10,10,0,function(n) print(“U set to “..n) end)

end

— This function gets called once every frame
function draw()

background(40, 40, 50)
k=1/(M*A)
h=1/M
h2=h^2

for m=1,4 do
uold[100]=0.1*U
for i = 1,M do
um[i]=(uold[i]+unew[i])/2
end
for i =2,M-1 do
D2=A*(um[i+1]-2*um[i]+um[i-1])/h2+10*F
vnew[i] = vold[i]+ k*D2
unew[i] = uold[i]+ (vold[i]+vnew[i])*k/2
end
end

for i=1,M do
uold[i]=unew[i]
vold[i]=vnew[i]
end

fill(100*uold[50], 0, 255, 255)
for i = 1,M-1 do
ellipse(8*i,500+200*uold[i],10,10)
end

end

Spel13: Sinus och Cosinus och talet Pi

Betrakta en punkformig kropp som rör sig med hastighet 1 i (motsols) cirkelrörelse med radie kring punkten (0,0) i ett vanligt rätvinkligt (x,y) koordinatsystem. Cirkelrörelsen beskrivs av tidsstegningen (försök förklara varför)

  • x = x – y*dt  (eller dx = -y*dt eller dx/dt = -y)
  • y = y + x*dt  (eller dy = x*dt eller dy/dt = x).

Programmera tidsstegningen med olika startvärdena för x och y och olika tidssteg dt. Försök bestämma sambandet till sinus och cosinus.

Hint: Med startvärdena x=1 och y=0 för t = 0 gäller att x=cos(t) och y=sin(t)  som funktioner av tiden t.

Bestäm sedan hur lång tid det tar för ett varv och därur ett värde på talet Pi.

Jämför med Calculus1 på App Store.

Template 1: (Följ cirkelrörelsen som beskrivs av funktionerna sin(t) och cos(t))

–Sinus och Cosinus

function setup()
x=1
y=0
dt=0.01
end

function draw()
background(255, 255, 255, 255)
strokeWidth(3)
stroke(0, 0, 0, 255)
line(200,0,200,1000)
line(0,500,800,500)

x=x-y*dt
y=y+x*dt

strokeWidth(0)
fill(0, 0, 0, 255)
ellipse(200*x+200,200*y+500,20)
fill(0, 255, 0, 255)
ellipse(200*x+200,500,20)
fill(255, 0, 0, 255)
ellipse(200,200*y+500,20)
stroke(0, 0, 255, 255)
strokeWidth(5)
line(200,500,200*x+200,200*y+500)
line(200,500,200*x+200,500)
line(200,500,200,500+200*y)
end