76 lines
1.6 KiB
Go
76 lines
1.6 KiB
Go
package main
|
|
|
|
import (
|
|
"fmt"
|
|
"math"
|
|
)
|
|
|
|
type fun = func(float64) float64
|
|
|
|
func aitken(f fun, p0 float64) float64 {
|
|
p1 := f(p0)
|
|
p2 := f(p1)
|
|
p1m0 := p1 - p0
|
|
return p0 - p1m0*p1m0/(p2-2.0*p1+p0)
|
|
}
|
|
|
|
func steffensenAitken(f fun, pinit, tol float64, maxiter int) float64 {
|
|
p0 := pinit
|
|
p := aitken(f, p0)
|
|
iter := 1
|
|
for math.Abs(p-p0) > tol && iter < maxiter {
|
|
p0 = p
|
|
p = aitken(f, p0)
|
|
iter++
|
|
}
|
|
if math.Abs(p-p0) > tol {
|
|
return math.NaN()
|
|
}
|
|
return p
|
|
}
|
|
|
|
func deCasteljau(c0, c1, c2, t float64) float64 {
|
|
s := 1.0 - t
|
|
c01 := s*c0 + t*c1
|
|
c12 := s*c1 + t*c2
|
|
return s*c01 + t*c12
|
|
}
|
|
|
|
func xConvexLeftParabola(t float64) float64 {
|
|
return deCasteljau(2.0, -8.0, 2.0, t)
|
|
}
|
|
|
|
func yConvexRightParabola(t float64) float64 {
|
|
return deCasteljau(1.0, 2.0, 3.0, t)
|
|
}
|
|
|
|
func implicitEquation(x, y float64) float64 {
|
|
return 5.0*x*x + y - 5.0
|
|
}
|
|
|
|
func f(t float64) float64 {
|
|
x := xConvexLeftParabola(t)
|
|
y := yConvexRightParabola(t)
|
|
return implicitEquation(x, y) + t
|
|
}
|
|
|
|
func main() {
|
|
t0 := 0.0
|
|
for i := 0; i < 11; i++ {
|
|
fmt.Printf("t0 = %0.1f : ", t0)
|
|
t := steffensenAitken(f, t0, 0.00000001, 1000)
|
|
if math.IsNaN(t) {
|
|
fmt.Println("no answer")
|
|
} else {
|
|
x := xConvexLeftParabola(t)
|
|
y := yConvexRightParabola(t)
|
|
if math.Abs(implicitEquation(x, y)) <= 0.000001 {
|
|
fmt.Printf("intersection at (%f, %f)\n", x, y)
|
|
} else {
|
|
fmt.Println("spurious solution")
|
|
}
|
|
}
|
|
t0 += 0.1
|
|
}
|
|
}
|