Math::Taylor - Taylor Polynomials and remainders

  use Math::Taylor;

  # Create new approximation
  my $approximation = Math::Taylor->new(
    function       => "sin(y) * cos(x)",
    point          => 2,
    variable       => 'y',
    remainder_type => 'cauchy',
  );

  # Calculate Taylor Polynomial of degree 2
  my $poly = $approximation->taylor_poly(2);
  print "$poly\n";

  # Upper bounds of the remainder are also availlable:
  my $remainder = $approximation->remainder(2, 'cauchy');

package Math::Taylor;

use 5.006;
use strict;
use warnings;

our $VERSION = '1.00';

use Math::Symbolic qw/parse_from_string/;
use Math::Symbolic::MiscCalculus;
use Carp qw/confess cluck/;

our $Default_Point      = 0;
our $Default_Variable   = Math::Symbolic::Variable->new('x');
our $Default_Remainder_Type = 'lagrange';

sub new {
    my $proto = shift;
    my $class = ref($proto) || $proto;

    confess "Math::Taylor called with uneven number of arguments."
      if @_ % 2;

    my %args = @_;

    # function to approximate,
    # variable of the function
    # point to approximate about
    my $self = {
        function       => undef,
        variable       => $Default_Variable,
        point          => $Default_Point,
        remainder_type => $Default_Remainder_Type,
    };

    # Clone prototype if applicable.
    if ( ref($proto) ) {
        $self->{function} = $proto->{function}->new();
        $self->{variable} = $proto->{variable}->new()
          if defined $proto->{variable};
        $self->{point} = $proto->{point} if defined $proto->{point};
        $self->{remainder_type} = $proto->{remainder_type}
          if defined $proto->{remainder_type};
    }

    bless $self => $class;

    $self->function( $args{function} ) if exists $args{function};
    $self->variable( $args{variable} ) if exists $args{variable};
    $self->point( $args{point} )       if exists $args{point};

    confess "Cannot create a Math::Taylor object without at least a function."
      if not defined $self->function();

    return $self;
}

sub function {
    my $self = shift;
    if ( not @_ ) {
        return $self->{function};
    }
    my $function = shift;

    if ( not defined $function ) {
        confess "Won't set the function of a Math::Taylor object to 'undef'.\n";
    }
    elsif ( ref($function) =~ /^Math::Symbolic/ ) {
        $self->{function} = $function;
    }
    elsif ( not ref($function) ) {
        my $parsed;
        eval { $parsed = parse_from_string($function); };
        if (   $@
            or not defined $parsed
            or not ref($parsed) =~ /^Math::Symbolic/ )
        {
            confess <<"HERE"
Could not parse function of Math::Taylor object as Math::Symbolic tree.
Argument was: '$function'
Error (if any) was: '$@'
HERE
        }
        $self->{function} = $parsed;
    }
    return $self->{function};
}

sub point {
    my $self = shift;
    if ( not @_ ) {
        return $self->{point};
    }
    my $point = shift;
    if ( defined $point ) {
        $self->{point} = $point;
    }
    else {
        confess
"Cannot set the 'point' attribute of a Math::Taylor object to 'undef'.";
    }
    return $self->{point};
}

sub variable {
    my $self = shift;
    if ( not @_ ) {
        return $self->{variable};
    }
    my $variable = shift;
    if (
        defined $variable
        and ( ref($variable) eq 'Math::Symbolic::Variable'
            or not ref($variable) )
      )
    {
        if ( not ref($variable) ) {
            my $parsed;
            eval { $parsed = Math::Symbolic::Variable->new($variable); };
            if (   $@
                or not defined $parsed
                or not ref($parsed) =~ /^Math::Symbolic::Variable/ )
            {
                confess <<"HERE"
Could not parse variable of Math::Taylor object as Math::Symbolic::Variable.
Argument was: '$variable'
Error (if any) was: '$@'
HERE
            }
            $variable = $parsed;
        }

        # Is the variable contained in the function at all?
        $self->_is_variable_in_function( $self->{function}, $variable );
    }
    else {
        confess
          "Tried to create a variable for Math::Taylor object from dubious"
          . "source. Source: '"
          . ( !defined($variable) ? 'undef' : $variable ) . "'";
    }
    $self->{variable} = $variable;
    return $variable;
}

# Internal method to test whether a given variable is actually part of
# a function's signature.
sub _is_variable_in_function {
    my $self      = shift;
    my $function  = shift;
    my $var       = shift;
    my %signature = map { ( $_, undef ) } $function->explicit_signature();
    if ( not exists $signature{ $var->to_string() } ) {
        confess <<"HERE";
Variable not contained in function to approximate.
Variable: '$var'
Function: '$function'
HERE
    }
}

sub remainder_type {
    my $self = shift;
    if ( not @_ ) {
        return $self->{remainder_type};
    }
    my $err = shift;
    if ( defined $err and $err eq 'lagrange' or $err eq 'cauchy' ) {
        $self->{remainder_type} = $err;
    }
    else {
        confess
"Cannot set the 'remainder_type' attribute of a Math::Taylor object to anything\n"
          . "but 'cauchy' or 'lagrange'.";
    }
    return $self->{remainder_type};
}

sub taylor_poly {
	my $self = shift;
	my $degree = shift;
	$degree = 1 if not defined $degree;
	confess("The degree of a Taylor approximation has to be >= 0.")
		unless $degree >= 0;

	# Get all necessary data.
	my $function = $self->function()->new();
	my $variable = $self->variable();
	my $pos = $self->point();

	# Make sure we don't have any vars in the function that clash with
	# the nomenclature of the output of TaylorPolynomial:
	# If TaylorPolynomial uses "x" as variable, it will include "x_0" as
	# a new variable in the output.
	my $posname = $variable->to_string() . '_0';	
	my %sig = map {($_, undef)} $function->explicit_signature();
	my @replace;
	if (exists $sig{$posname}) {
		my $newname = $posname;
		while (exists $sig{$newname}) {
			$newname .= '_';
		}
		$function->implement($posname => $newname);
		@replace = ($newname, $posname);
	}

	my $poly = Math::Symbolic::MiscCalculus::TaylorPolynomial(
		$function,
		$degree,
		$variable,
		$pos
	);
	
	if (not defined $poly or not ref($poly) =~ /^Math::Symbolic/) {
		confess(
		"Could not calculate Taylor approximation of degree $degree using\n"
		."function '$function',\n position '$pos',\n and variable '$variable'."
		);
	}

	# Insert value for x_0
	$poly->implement($posname => $pos);
	
	# Undo all changes to variable names.
	if (@replace) {
		$poly->implement(@replace);
	}

	return $poly;
}

sub remainder {
	my $self = shift;
	my $degree = shift;
	$degree = 1 if not defined $degree;
	confess("The degree of a Taylor approximation has to be >= 0.")
		unless $degree >= 0;
	
	my $tvar = shift;
	$tvar = Math::Symbolic::Variable->new('theta') if not defined $tvar;
	$tvar = Math::Symbolic::Variable->new($tvar);
	
	# Get all necessary data.
	my $function = $self->function()->new();
	my $variable = $self->variable();
	my $pos = $self->point();
	my $type = $self->{remainder_type};

	# Make sure we don't have any vars in the function that clash with
	# the nomenclature of the output of TaylorPolynomial:
	# If TaylorPolynomial uses "x" as variable, it will include "x_0" as
	# a new variable in the output.
	my $posname = $variable->to_string() . '_0';	
	my %sig = map {($_, undef)} $function->explicit_signature();
	my @replace;
	if (exists $sig{$posname}) {
		my $newname = $posname;
		while (exists $sig{$newname}) {
			$newname .= '_';
		}
		$function->implement($posname => $newname);
		@replace = ($newname, $posname);
	}

	# get remainder;
	my $rem;
	if ($type eq 'lagrange') {
		$rem = Math::Symbolic::MiscCalculus::TaylorErrorLagrange(
			$function, $degree, $variable, $pos, $tvar
		);
	} else {
		$rem = Math::Symbolic::MiscCalculus::TaylorErrorCauchy(
			$function, $degree, $variable, $pos, $tvar
		);
	}
	
	# Insert value for x_0
	$rem->implement($posname => $pos);
	
	# Undo all changes to variable names.
	if (@replace) {
		$rem->implement(@replace);
	}

	return $rem;
}





1;
__END__