REAL FUNCTION sekante (polynom, x0, abbruch)
 IMPLICIT NONE
 REAL    :: polynom(:), abbruch, x0
 REAL    :: x1, a, b, y0, y1, delta
 INTEGER :: it, maxit

 INTERFACE
  REAL FUNCTION funktionswert(f,x)
   IMPLICIT NONE
   REAL :: f(0:), x
  END FUNCTION
 END INTERFACE

 abbruch = MAX(abbruch,EPSILON(abbruch))
 maxit = 1000
 it = 0
 delta = abbruch/10  ! beliebig gewaehlter Wert, koennte uU nicht gut sein
 DO
  it = it + 1
  x1 = x0 + delta
  y0 = funktionswert(polynom,x0)
  y1 = funktionswert(polynom,x1)
  a = (y1 - y0) / delta
  b = y0 - a*x0
  IF (ABS(a) < EPSILON(x0)) THEN   ! sonst Division durch 0
   delta = delta / 2.0
   CYCLE
  ENDIF
  x1 = - b / a
  IF (ABS(x1-x0) < abbruch) THEN
   sekante = x0
   RETURN         ! Nullstelle gefunden
  ENDIF
  IF (it == maxit) THEN
   sekante = x0   ! letzter Wert - aber:
   RETURN         ! keine Nullstelle gefunden
  ENDIF
  x0 = x1
  y0 = funktionswert(polynom,x0)
 END DO
END FUNCTION
