Another solution based on the fact that a floating point number is actually a sum of powers of 2.

For 0 <= y <= 0, it holds that y = b₀*2⁰ + b₁*2^-1 + b₂*2^-2 + b₃*2^-3 + ...

And so, x^y = x^{b₀*20 + b₁*2-1 + b₂*2-2 + b₃*2-3 + ...} = x^b₀*20 * x^b₁*2-1 * x^b₂*2-2 * x^b₃*2-3 * ... = Π_{i, bi≠0} x^2-i

Which in Perl becomes...

#!/usr/bin/perl

use strict;
use warnings;

my ($X, $Y, $n) = @ARGV;
$n ||= 20;

my $Z = $X ** $Y;

my ($neg, $y) = ($Y < 0 ? (1, -$Y) : (0, $Y));

$y > 1 and die "Y is out of range";
my $bit = 1;
my $x = $X;
my $z = 1;
while ($y) {
    if ($y >= $bit) {
        $y -= $bit;
        $z *= $x;
    }
    $bit /= 2;
    $x = sqrt($x);
}

$z = 1/$z if $neg;
my $e = abs($z - $X ** $Y);
print "z=$z, e=$e\n";
[download]

sqrt is also an expensive operation. If you have to calculate x^y for different values of y while x stays constant, you may be able to speed up the process creating a table with the values of x^2-i.

Alternatively, you can write a function to calculate e^x (using a table with the values of e²ⁱ), and then calculate x^y as e^log(x)*y.

update: there are also several implementations of exp(x) freely available, for instance, the one in OpenBSD is here.