Steffensen's method: Difference between revisions
m
→{{header|Phix}}: use pygments
(Add MATLAB implementation) |
m (→{{header|Phix}}: use pygments) |
||
| (3 intermediate revisions by 3 users not shown) | |||
Line 717:
<pre>
Same as C example
</pre>
=={{header|Haxe}}==
{{Trans|ATS|via ALGOL 68}}
<syntaxhighlight lang="haxe">
// Steffensens's method - translated from the ATS sample via Algol 68
class OptionalReal
{
public var v:Float;
public var hasValue:Bool;
public function new( value:Float, notNull:Bool )
{
v = value;
hasValue = notNull;
};
}
class Main {
// Aitken's extrapolation
static function aitken( f:Float->Float // function Float to Float
, p0:Float // initial fixed point estimate
):Float
{
var p1 = f( p0 );
var p2 = f( p1 );
var p1m0 = p1 - p0;
return p0 - ( p1m0 * p1m0 ) / ( p2 - ( 2 * p1 ) + p0 );
};
// finds fixed point p such that f(p) = p
static function steffensenAitken( f:Float->Float // function Float to Float
, pinit:Float // initial estimate
, tol:Float // tolerance
, maxiter:Float // maximum number of iterations
):OptionalReal // return a Float, IF tolerance is met
{
var p0 = pinit;
var p = aitken( f, p0 ); // first iteration
var i = 1;
while( i < ( maxiter - 1 ) && Math.abs( p - p0 ) > tol ) {
p0 = p;
p = aitken( f, p0 );
i += 1;
}
if( Math.abs( p - p0 ) > tol ) {
return new OptionalReal( 0, false );
} else {
return new OptionalReal( p, true );
}
};
static function deCasteljau( c0:Float // control point coordinates (one axis)
, c1:Float
, c2:Float
, t :Float // the indep}ent parameter
) :Float // value of x(t) or y(t)
{
var s = 1 - t;
var c01 = ( s * c0 ) + ( t * c1 );
var c12 = ( s * c1 ) + ( t * c2 );
return ( s * c01 ) + ( t * c12 );
};
static function xConvexLeftParabola( t:Float ):Float { return deCasteljau( 2, -8, 2, t ); }
static function yConvexLeftParabola( t:Float ):Float { return deCasteljau( 1, 2, 3, t ); }
static function implicitEquation( x:Float, y:Float ):Float { return ( 5 * x * x ) + y - 5; }
static function f( t:Float ):Float // find fixed points of this function
{
var x = xConvexLeftParabola( t );
var y = yConvexLeftParabola( t );
return implicitEquation( x, y ) + t;
};
static function main()
{
var t0:Float = 0;
for( i in 1...12 ) {
Sys.print( 't0 = $t0 : ' );
var result = steffensenAitken( f, t0, 0.00000001, 1000 );
if( ! result.hasValue ) {
Sys.println( 'no answer' );
} else {
var x = xConvexLeftParabola( result.v );
var y = yConvexLeftParabola( result.v );
if( Math.abs( implicitEquation( x, y ) ) <= 0.000001 ) {
Sys.println( 'intersection at ($x, $y)' );
} else {
Sys.println( 'spurious solution' );
}
}
t0 += 0.1;
}
}
}
</syntaxhighlight>
{{out}}
<pre>
t0 = 0 : no answer
t0 = 0.1 : intersection at (0.881024967590665, 1.11897503240933)
t0 = 0.2 : no answer
t0 = 0.3 : no answer
t0 = 0.4 : no answer
t0 = 0.5 : no answer
t0 = 0.6 : no answer
t0 = 0.7 : no answer
t0 = 0.8 : no answer
t0 = 0.9 : intersection at (-0.68102496759051, 2.68102496759069)
t0 = 1 : no answer
</pre>
Line 1,125 ⟶ 1,238:
=={{header|Phix}}==
{{trans|C}}
<!--
<syntaxhighlight lang="phix">
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
</syntaxhighlight>
Output same as C
Line 1,425 ⟶ 1,539:
t0 = .9 : intersection at ( -.681024968, 2.68102497)
t0 = 1. : no answer
</pre>
=={{header|Scala}}==
{{trans|Java}}
<syntaxhighlight lang="Scala">
object Steffensen {
def aitken(p0: Double): Double = {
val p1 = f(p0)
val p2 = f(p1)
val p1m0 = p1 - p0
p0 - p1m0 * p1m0 / (p2 - 2.0 * p1 + p0)
}
def steffensenAitken(pinit: Double, tol: Double, maxiter: Int): Option[Double] = {
var p0 = pinit
var p = aitken(p0)
var iter = 1
while (math.abs(p - p0) > tol && iter < maxiter) {
p0 = p
p = aitken(p0)
iter += 1
}
if (math.abs(p - p0) > tol) None else Some(p)
}
def deCasteljau(c0: Double, c1: Double, c2: Double, t: Double): Double = {
val s = 1.0 - t
val c01 = s * c0 + t * c1
val c12 = s * c1 + t * c2
s * c01 + t * c12
}
def xConvexLeftParabola(t: Double): Double = deCasteljau(2.0, -8.0, 2.0, t)
def yConvexRightParabola(t: Double): Double = deCasteljau(1.0, 2.0, 3.0, t)
def implicitEquation(x: Double, y: Double): Double = 5.0 * x * x + y - 5.0
def f(t: Double): Double = {
val x = xConvexLeftParabola(t)
val y = yConvexRightParabola(t)
implicitEquation(x, y) + t
}
def main(args: Array[String]): Unit = {
var t0 = 0.0
for (i <- 0 until 11) {
print(f"t0 = $t0%3.1f : ")
steffensenAitken(t0, 0.00000001, 1000) match {
case None => println("no answer")
case Some(t) =>
val x = xConvexLeftParabola(t)
val y = yConvexRightParabola(t)
if (math.abs(implicitEquation(x, y)) <= 0.000001) {
println(f"intersection at ($x%f, $y%f)")
} else {
println("spurious solution")
}
}
t0 += 0.1
}
}
}
</syntaxhighlight>
{{out}}
<pre>
t0 = 0.0 : no answer
t0 = 0.1 : intersection at (0.881025, 1.118975)
t0 = 0.2 : no answer
t0 = 0.3 : no answer
t0 = 0.4 : no answer
t0 = 0.5 : no answer
t0 = 0.6 : no answer
t0 = 0.7 : no answer
t0 = 0.8 : no answer
t0 = 0.9 : intersection at (-0.681025, 2.681025)
t0 = 1.0 : no answer
</pre>
Line 1,431 ⟶ 1,624:
{{libheader|Wren-fmt}}
Wren's only built-in numeric type is a 'double' (in C) so no surprise that we're getting the same answers as ATS when using the same tolerances and maximum number of iterations.
<syntaxhighlight lang="
var aitken = Fn.new { |f, p0|
| |||