61 lines
1.5 KiB
Plaintext
61 lines
1.5 KiB
Plaintext
with javascript_semantics
|
|
function aitken(integer f, atom p0)
|
|
atom p1 = f(p0),
|
|
p2 = f(p1),
|
|
p1m0 = p1 - p0
|
|
return p0 - (p1m0 * p1m0) / (p2 - (2 * p1) + p0)
|
|
end function
|
|
|
|
function steffensenAitken(integer f, maxiter, atom pinit, tol)
|
|
atom p0 = pinit, p = aitken(f, p0)
|
|
integer iter = 1
|
|
while abs(p-p0)>tol and iter<maxiter do
|
|
p0 = p
|
|
p = aitken(f, p0)
|
|
iter += 1
|
|
end while
|
|
return iff(abs(p-p0)>tol?"none":p)
|
|
end function
|
|
|
|
function deCasteljau(atom c0, c1, c2, t)
|
|
atom s = 1 - t,
|
|
c01 = (s * c0) + (t * c1),
|
|
c12 = (s * c1) + (t * c2),
|
|
c012 = (s * c01) + (t * c12)
|
|
return c012
|
|
end function
|
|
|
|
function xConvexLeftParabola(atom t)
|
|
return deCasteljau(2, -8, 2, t)
|
|
end function
|
|
|
|
function yConvexLeftParabola(atom t)
|
|
return deCasteljau(1, 2, 3, t)
|
|
end function
|
|
|
|
function implicitEquation(atom x, y)
|
|
return (5 * x * x) + y - 5
|
|
end function
|
|
|
|
function f(atom t)
|
|
atom x = xConvexLeftParabola(t),
|
|
y = yConvexLeftParabola(t)
|
|
return implicitEquation(x, y) + t
|
|
end function
|
|
|
|
for i=0 to 10 do
|
|
printf(1,"t0 = %.1f : ",i/10)
|
|
object t = steffensenAitken(f, 1000, i/10, 0.00000001)
|
|
if string(t) then
|
|
printf(1,"no answer\n")
|
|
else
|
|
atom x = xConvexLeftParabola(t),
|
|
y = yConvexLeftParabola(t)
|
|
if abs(implicitEquation(x, y)) <= 0.000001 then
|
|
printf(1,"intersection at (%.6f, %.6f)\n",{x,y})
|
|
else
|
|
printf(1,"spurious solution\n")
|
|
end if
|
|
end if
|
|
end for
|