RosettaCodeData/Task/Steffensens-method/FreeBASIC/steffensens-method.basic

70 lines
1.9 KiB
Plaintext

Function aitken(f As Function(As Double) As Double, p0 As Double) As Double
Dim As Double p1 = f(p0)
Dim As Double p2 = f(p1)
Dim As Double p1m0 = p1 - p0
Return p0 - (p1m0 * p1m0) / (p2 - (2 * p1) + p0)
End Function
Function steffensenAitken(f As Function(As Double) As Double, pinit As Double, tol As Double, maxiter As Integer) As Double
Dim As Double p0 = pinit
Dim As Double p = aitken(f, p0)
Dim As Integer iter = 1
While Abs(p-p0) > tol Andalso iter < maxiter
p0 = p
p = aitken (f, p0)
iter += 1
Wend
If Abs (p - p0) > tol Then Return 0 Else Return p
End Function
Function deCasteljau(c0 As Double, c1 As Double, c2 As Double, t As Double) As Double
Dim As Double s = 1 - t
Dim As Double c01 = (s * c0) + (t * c1)
Dim As Double c12 = (s * c1) + (t * c2)
Dim As Double c012 = (s * c01) + (t * c12)
Return c012
End Function
Function xConvexLeftParabola(t As Double) As Double
Return deCasteljau(2, -8, 2, t)
End Function
Function yConvexLeftParabola(t As Double) As Double
Return deCasteljau(1, 2, 3, t)
End Function
Function implicitEquation(x As Double, y As Double) As Double
Return (5 * x * x) + y - 5
End Function
Function f(t As Double) As Double
Dim As Double x = xConvexLeftParabola(t)
Dim As Double y = yConvexLeftParabola(t)
Return implicitEquation(x, y) + t
End Function
Dim As Double t0 = 0
For i As Integer = 0 To 10
Print Using "t0 = #.# : "; t0;
Dim As Double t = steffensenAitken(@f, t0, 0.00000001, 1000)
If t = 0 Then
Print "no answer"
Else
Dim As Double x = xConvexLeftParabola(t)
Dim As Double y = yConvexLeftParabola(t)
If Abs(implicitEquation(x, y)) <= 0.000001 Then
Print Using "intersection at (##.######, ##.######)"; x; y
Else
Print "spurious solution"
End If
End If
t0 += 0.1
Next
Sleep