Roots of a function

Create a program that finds an outputs the roots of a given function, range and (if applicable) step width. The program should identify whether the root is exact or approximate.

For this example, use f(x)=x^3-3x^2+2x.

Ada

with Ada.Text_Io; use Ada.Text_Io;

procedure Roots_Of_Function is
   package Real_Io is new Ada.Text_Io.Float_Io(Long_Float);
   use Real_Io;
   
   function F(X : Long_Float) return Long_Float is
   begin
      return (X**3 - 3.0*X*X + 2.0*X);
   end F;
   
   Step  : constant Long_Float := 1.0E-6;
   Start : constant Long_Float := -1.0;
   Stop  : constant Long_Float := 3.0;
   Value : Long_Float := F(Start);
   Sign  : Boolean := Value > 0.0;
   X     : Long_Float := Start + Step;
 
begin
   if Value = 0.0 then
      Put("Root found at ");
      Put(Item => Start, Fore => 1, Aft => 6, Exp => 0);
      New_Line;
   end if;
   while X <= Stop loop
      Value := F(X);
      if (Value > 0.0) /= Sign then
         Put("Root found near ");
         Put(Item => X, Fore => 1, Aft => 6, Exp => 0);
         New_Line;
      elsif Value = 0.0 then
         Put("Root found at ");
         Put(Item => X, Fore => 1, Aft => 6, Exp => 0);
         New_Line;
      end if;
      Sign := Value > 0.0;
      X := X + Step;
   end loop;
end Roots_Of_Function;

C++

#include <iostream>

double f(double x)
{
	return (x*x*x - 3*x*x + 2*x);
}

int main()
{
	double step = 0.001; // Smaller step values produce more accurate and precise results
	double start = -1;
	double stop = 3;
	double value = f(start);
	double sign = (value > 0);
	
	// Check for root at start
	if ( 0 == value )
		std::cout << "Root found at " << start << std::endl;

	for(	double x = start + step;
			x <= stop;
			x += step )
	{
		value = f(x);
		
		if ( ( value > 0 ) != sign )
			// We passed a root
			std::cout << "Root found near " << x << std::endl;
		else if ( 0 == value )
			// We hit a root
			std::cout << "Root found at " << x << std::endl;
		
		// Update our sign
		sign = ( value > 0 );
	}
}

D

module findroot ;
import std.stdio ;
import std.math ;

void report(T)(T[] r, T function(T) f, T tolerance = cast(T) 1e-4L) {
  if (r.length) {
    writefln("Root found (tolerance = %1.4g) :", tolerance) ;
    foreach(x ; r) {
      T y = f(x) ;    
      if (nearZero(y)) 
        writefln("... EXACTLY at %+1.20f, f(x) = %+1.4g", x, y) ;
      else if (nearZero(y, tolerance))
        writefln(".... MAY-BE at %+1.20f, f(x) = %+1.4g", x, y) ;
      else
        writefln("Verify needed, f(%1.4g) = %1.4g > tolerance in magnitude", x, y) ;
    }
  } else 
    writefln("No root found.") ;
}

bool nearZero(T)(T a, T b = T.epsilon * 4) { return abs(a) <= b ; }

T[] findroot(T)(T function(T) f, T start, T end, T step = cast(T) 0.001L, 
                                 T tolerance = cast(T) 1e-4L) {
  T[T] result ;

  if (nearZero(step))
    writefln("WARNING: step size may be too small.") ;

  T searchRoot(T a, T b) { // search root by simple bisection
    T root ;
    int limit = 49 ;
    T gap = b - a ;
    while (!nearZero(gap) && limit--) {
      if (nearZero(f(a))) return a ;
      if (nearZero(f(b))) return b ;
      root = (b + a)/2.0L ;
      if (nearZero(f(root))) return root ;
      if (f(a) * f(root) < 0)
        b = root ;
      else
        a = root ;
      gap = b - a ;           
    }
    return root ;
  }
  
  T dir = cast(T) (end > start ? 1.0 : -1.0) ;
  step = (end > start) ? abs(step) : - abs(step) ;
  for(T x = start ; x*dir <= end*dir ; x = x + step)
    if (f(x)*f(x + step) <= 0) {
      T r = searchRoot(x, x+ step) ;        
      result[r] = f(r) ;
    } 
  return result.keys.sort ; // reduce duplacated root, if any
}

real f(real x){
  return x*x*x - 3*x*x + 2*x ;
}

void main(){
    findroot(&f, -1.0L, 3.0L, 0.001L).report(&f) ;
}

Output ( NB:smallest increment for real type in D is real.epsilon = 1.0842e-19 ):

Root found (tolerance = 0.0001) :
.... MAY-BE at -0.00000000000000000080, f(x) = -1.603e-18
... EXACTLY at +1.00000000000000000020, f(x) = -2.168e-19
.... MAY-BE at +1.99999999999999999950, f(x) = -8.674e-19

Perl

sub f
{
        my $x = shift;

        return ($x * $x * $x - 3*$x*$x + 2*$x);
}

my $step = 0.001; # Smaller step values produce more accurate and precise results
my $start = -1;
my $stop = 3;
my $value = &f($start);
my $sign = $value > 0;

# Check for root at start

print "Root found at $start\n" if ( 0 == $value );

for(    my $x = $start + $step;
        $x <= $stop;
        $x += $step )
{
        $value = &f($x);

        if ( 0 == $value )
        {
                # We hit a root
                print "Root found at $x\n";
        }
        elsif ( ( $value > 0 ) != $sign )
        {
                # We passed a root
                print "Root found near $x\n";
        }

        # Update our sign
        $sign = ( $value > 0 );
}