[2020 Oct 13] - fast, accurate, framerate-independent
Code: Select all
/' -- fast terminal velocity - 2020 Oct 13 - by dafhi
goal: fps-independence and simple physics
currently
1. acceleration is not fps-independent
2. necessary pre-calcs aren't straight-forward
here's the process:
dim sng dens = .6 '' > 0 and less than 1
dim sng v = velocity_const / velo_sum(dens, phys_fps)
dim sng a = accel_const / accel_sum(dens, phys_fps)
dim sng p '' pos
dim sng k = dens ^ (1 / phys_fps)
for i int = 0 to phys_fps*5 '' extend the curve
'' horizontal and vertical scaling
pset (i * 130 / phys_fps, h-1 - p * 99), col
v += a
v *= k
p += v
Next
'/
'' replaces Int()
#Define flo(x) (((x)*2.0-0.5)shr 1) '' http://www.freebasic.net/forum/viewtopic.php?p=118633
'' hacks
#undef int
#define int as integer
#define sng as single
function gs(k sng = .5, n sng = 2) sng
/'
-- geometric sum
sum = k + k^2 .. k^n
sum - k*sum = k + k^2 .. k^n
- k^2 - k^3 .. k^(n+1)
= k + k^2 .. k^n
- k^2 .. k^n - k^(n+1)
'/
return k*(1-k^n) / (1-k)
end function
function sgs(k sng, fps sng) sng
/'
-- sum of geometric sum
inputs particle_density, physics fps
output
position at time t
- theory
can summarize with a loop:
[ vel+=accel, vel*=k, pos+=vel ]
vel frame 1: k 'velocity = (velocity + accel) * k
vel frame 2: k+kk 'notice k in front each time (acceleration)
vel frame 3: k kk kkk 'geometric sum
pos frame 1: k 'position = velocity sum = sum of geometric sum
pos frame 2: k k kk
pos frame 3: k k kk k kk kkk
the pattern becomes visible to the experienced human
geometric sum = k*(1-k^n) / (1-k)
sum of geometric sum = k*(1-k^1)/(1-k) + .. k*(1-k^n)/(1-k)
representation 'r'
k/(1-k)
series
1-k^1 + .. 1-k^n
= n - [k^1 + .. k^n]
= n - geometric sum
solution
= r*( n - k*(1-k^n)/(1-k) )
= r*( n - r*(1-k^n) )
'/
static sng r
r=k/(1-k)
return r*(fps - r*(1-k^fps))
end function
function velo_sum(k sng = .5, fps sng = 2) sng
return gs( k^(1/fps), fps ) / k
end function
function accel_sum(k sng, fps sng) sng
return sgs( k^(1/fps), fps ) / k
end function
const w = 800, wm=w-1
const h = 600, hm=h-1
sub fps_plot( phys_fps sng, col as ulong = -1 )
dim sng dens = .6
dim sng v = 2 / velo_sum(dens, phys_fps)
dim sng a = 0.1 / accel_sum(dens, phys_fps)
dim sng p '' pos
dim sng k = dens ^ (1 / phys_fps)
for i int = 0 to phys_fps*5 '' extend the curve
'' horizontal and vertical scaling
pset (i * 130 / phys_fps, h-1 - p * 99), col
v += a
v *= k
p += v
Next
end sub
sub main
screenres w,h,32
const golden = (1+5^.5)/2
dim sng phys_fps = .25
for i int = 1 to 10
fps_plot phys_fps, rgb(rnd*255,rnd*255,rnd*255)
phys_fps *= golden
next
end sub
main
sleep
[2017 Sep 17] - accurate, framerate-independent
Code: Select all
#ifndef pi
const TwoPi = 8*atn(1)
const Pi = 4*atn(1)
#EndIf
type float as single
type v3
as float x,y,z
declare function rndf as v3
declare sub rand
End Type
sub v3.rand
y=2*(rnd-.5): var r=sqr(1-y*y)
z=rnd*twopi: x=r*cos(z): z=r*sin(z)
End Sub
function v3.rndf as v3
var y=2*(rnd-.5): var r=sqr(1-y*y)
var z=rnd*twopi: return type(r*cos(z), y, r*sin(z))
End function
operator -(r as v3) as v3: return type(-r.x, -r.y, -r.z): end operator
operator -(l as v3,r as v3) as v3: return type(l.x-r.x,l.y-r.y,l.z-r.z): end operator
operator +(l as v3,r as v3) as v3: return type(l.x+r.x, l.y+r.y, l.z+r.z): end operator
operator /(l as v3,r as float) as v3: dim as float s = 1/r: return type(l.x*s,l.y*s,l.z*s): end operator
operator *(l as v3,r as float) as v3: return type(l.x*r,l.y*r,l.z*r): end operator
operator *(l as float, r as v3) as v3: return type(l*r.x,l*r.y,l*r.z): end operator
operator *(l as v3,r as v3) as v3: return type(l.x*r.x,l.y*r.y,l.z*r.z): end operator
type tTermVel
as single v, a
as single k=.5, fps=67
declare sub ini(k as single = .5, fps as single=67)
declare sub dpos(i as single)
declare sub new_kinetics
as single ifps=1/67
declare constructor(k as single = 5, fps as single=67)
private:
as single ik, iok, i, j, kj
End Type
constructor.tTermVel(k as single = 5, fps as single=67): ini k, fps
end constructor
sub tTermVel.ini(_k as single, _fps as single): fps = _fps
ifps = 1/fps: k=_k
new_kinetics: j=0: kj=1
End Sub
sub tTermVel.new_kinetics
ik = 1/(1-k)
iok = 1/(1-1/k)
End Sub
sub tTermVel.dpos(i as single)
var ki = k^i
v = k*(kj-ki)
a = ((i-j)*k + v*iok)*ik
v *= ik
j=i: kj=ki
End Sub
const w = 800, wm=w-1, midx = wm/2
const h = 600, hm=h-1, midy = hm/2
screenres w,h,32
sub insta_plot(fps as single, col as ulong=-1, frames as ubyte=30)
dim as v3 p = type(100, 200, 0)
dim as v3 v = type(20,20) * 3
dim as v3 a = type(0,-1,0) * 3
var k = .85
dim as tTermVel t=type(k, fps)
var ifps = 1/fps, age = rnd*ifps
for n as long = 1 to fps*frames
t.dpos age
p += v*t.v + a*t.a
pset ( p.x, hm - p.y ), col
age += ifps
Next
end sub
var l = 1
var u = 7
for fps as single = u to l step -1
var col0 = 70
var col = col0 + (255-col0)*(u+1-fps)/(u+1-l)
insta_plot fps, rgb(col,col,col)
Next
sleep 2000
Code: Select all
#ifndef pi
const TwoPi = 8*atn(1)
const Pi = 4*atn(1)
#EndIf
type float as single
type v3
as float x,y,z
declare function rndf as v3
declare sub rand
End Type
sub v3.rand
y=2*(rnd-.5): var r=sqr(1-y*y)
z=rnd*twopi: x=r*cos(z): z=r*sin(z)
End Sub
function v3.rndf as v3
var y=2*(rnd-.5): var r=sqr(1-y*y)
var z=rnd*twopi: return type(r*cos(z), y, r*sin(z))
End function
operator -(r as v3) as v3: return type(-r.x, -r.y, -r.z): end operator
operator -(l as v3,r as v3) as v3: return type(l.x-r.x,l.y-r.y,l.z-r.z): end operator
operator +(l as v3,r as v3) as v3: return type(l.x+r.x, l.y+r.y, l.z+r.z): end operator
operator /(l as v3,r as float) as v3: dim as float s = 1/r: return type(l.x*s,l.y*s,l.z*s): end operator
operator *(l as v3,r as float) as v3: return type(l.x*r,l.y*r,l.z*r): end operator
operator *(l as float, r as v3) as v3: return type(l*r.x,l*r.y,l*r.z): end operator
operator *(l as v3,r as v3) as v3: return type(l.x*r.x,l.y*r.y,l.z*r.z): end operator
function gsum(k as single=.5, n as single=2) as double
/' --- geometric sum ---
1. s = k+k^2 .. k^n
2. s-ks = k+k^2 .. k^n
- k^2 .. k^n - k^(n+1)
____________________
= k - k^(n+1) '/
return (k-k^(n+1)) / (1-k)
end function
function vgsum(k as single, i as single) as single
/' --- velo grav sum ---
Loop
[ vel+=accel, vel*=k, pos+=vel ]
vel frame 1: k
pos frame 1: k
vel frame 2: kk+k
pos frame 2: kk+k+k
- solving for pos -
p = 2k1 + 1k2
pk = 2k2 + 1k3
p(1-k) = 2k1 + 1k2 - 2k2 - 1k3
f = 2k1
g = - 1k2 - 1k3
g/k = - 1k - 1k2
g(1-1/k) = - k3 + k
g = (k-k3) / (1 - 1/k)
f+g = 2k + (k-k3) / (1-1/k)
p =(2k + (k-k3) / (1-1/k)) / (1-k)
'/
return ( i*k + (k-k^(i+1)) / (1-1/k) ) / (1-k)
end function
const w = 800, wm=w-1, midx = wm/2
const h = 600, hm=h-1, midy = hm/2
screenres w,h,32
sub insta_plot(fps as single, col as ulong=-1, frames as ubyte=30)
dim as v3 p = type(100, 200, 0)
dim as v3 v = type(20,20) * 3
dim as v3 a = type(0,-1,0) * 3
var k = .85
var ifps = 1/fps
var age = rnd*ifps
var ik = 1/(1-k)
var iok = 1/(1-1/k)
var j=csng(0)
for n as long = 1 to fps*frames
p += v*(gsum(k,age)-gsum(k,j)) + _ 'pos delta is based upon absolute pos: current - previous
a*(vgsum(k,age)-vgsum(k,j))
'logic for absolute pos is shown in function tvel(),
'at the end of this source code
pset ( p.x, hm - p.y ), col
j = age
age += ifps
Next
end sub
var l = 1
var u = 7
for fps as single = u to l step -1
var col0 = 70
var col = col0 + (255-col0)*(u+1-fps)/(u+1-l)
insta_plot fps, rgb(col,col,col)
Next
sleep
function tvel(a as single, k as single, t as single, v0 as single=0) as single
/' --- terminal velocity - solved using "recurrence relations"
'pos(time+1) = pos(time) + (vel(time) + accel)*k
fn+1 = fn+(vn+a)k <--- shorthand
fn+2 = fn+1 + (vn+1 + a)k
= fn+(vn+a)k + (vn+1 + a)k
vn+1 = (vn+a)k
fn+2 = fn+(vn+a)k + ( (vn+a)k + a)k
= fn + (vn+a)k + (vn+a)kk + ak
fn+3 = fn+1 + (vn+1+a)k + (vn+1+a)kk + ak
= fn+(vn+a)k + ((vn+a)k+a)k + ((vn+a)k+a)kk + ak
= fn + vnk + ak + vnkk + akk + ak + vnkkk + akkk + akk + ak
= fn + vn*(k+kk+kkk) + a*(3k + 2kk + 1kkk)
= fn + vn*geom sum + a*poly
'/
'return v0*gsum(k,j) + a*vgsum(k,j)
dim as single kpow = k-k^(t+1)
return ( v0*kpow + a*( t*k + kpow/(1-1/k) ) )/(1-k)
End Function