use constant PI => 3.141592653589793238; ... sub factorial { my $n = shift; return (sqrt((2*$n + 1/3)*PI)*($n**$n)*exp(-$n)); }