Laboratoire de l'Informatique et du Parallélisme » MEPlib

#! /usr/bin/perl

# checkSolutions-02.pl - ST - 2006-03-15

#

# Check the quality of and sort a set of solutions with the Maple infnorm.

#

# COMMAND LINE SYNTAX:

# checkSolutions-01.pl  digitsFile solutionsFile function lowerBound upperBound

# NOTE(S):

# Changes from checkSolutions-02.pl: This versions can use several bit sizes

# ---------------------------------- instead of a single one for all 

#                                    coefficients

#

# This version checks the solutions produced by the "gp" PARI utility

# Use only for "monochunk" coefficients. For multiple chunks 

# coefficients use the "checkSimpleSolutions" utilities family.

##############################################################################

# PROGRAM DECLARATIONS

##############################################################################

BEGIN

{

  # Modify @INC here if necessary

  my $scriptDirName = `dirname $0`;

  chomp($scriptDirName);

  unshift(@INC,$scriptDirName);

  unshift(@INC, '../utils');

};

# core modules

use strict; # must be first to force strict in all others

use Carp;

use Getopt::Long;

use Math::BigInt;

use Math::BigRat;

use Math::BigFloat;

# custom modules

# Uncomment if you want to use standard exit error codes.

# the stdExitErrors.pm file must be available around.

use stdExitErrors;

##############################################################################

# PROGRAM CONSTANT DECLARATIONS

##############################################################################

my $bigIntZero = Math::BigInt->new('0');

my $debug = 1;

my $variableName = "x"; 

my $mapleApplication = "/usr/local/maple/bin/maple";

##############################################################################

# PROGRAM VARIABLE DECLARATIONS

##############################################################################

my $digitsFileName    = "";

my $solutionsFileName = "";

my $positionsFileName = "";

my $numDigits         = 0;

my $numDigitsData     = 0;

my @digits            = ();

my $function          = "";

my $lowerBound        = "";

my $upperBound        = "";

my $elem              = "";

my $line              = "";

my $lineNum           = 0;

my @line              = ();

my $isFirstLine       = 1;

my $powerOfTen        = 0;

my $solutionVectorRef = 0;

my @solutionMatrix    = ();

my @positionsVector   = ();

my $powerIndex        = 0;

my $polynomialDegree  = 0;

my $polynomialString  = "";

my @polynomials       = ();

my $mapleProgram      = "";

my $mapleResult       = 0;

my $infNorm           = "";

my @infNorm           = ();

my %sortHash          = ();

my @sortedInfNorms    = ();

#

#

#

##############################################################################

# TEST ENVIRONMENT AND CONFIGURATION VARIABLES

##############################################################################

##############################################################################

# INITIALIZE DEBUG FILE

##############################################################################

##############################################################################

# INITIALIZE LOG FILE

##############################################################################

##############################################################################

# GET COMMAND LINE ARGUMENTS

##############################################################################

##############################################################################

# PROGRAM MAIN

##############################################################################

&usage();

#

# Check the command line parameters

#

# It's a TODO!

#

# Open the files

#

# Open the digits file

#

if (! open(DF, "<$digitsFileName"))

{

  print STDERR "\n\nCan't open $digitsFileName. Aborting the program!\n\n";

  exit($EX_NOINPUT);

}

# Next the solution vector file.

#

if (! open(SVF, "<$solutionsFileName"))

{

  print STDERR "\n\nCan't open $solutionsFileName. Aborting the program!\n\n";

  exit($EX_NOINPUT);

}

#

# Read the digits file

#

while ($line = <DF>)

{

  # Remove the '\n'.

  chomp($line);

  # Coalescing all consecutive spaces to one.

  $line =~ s/(\s)+/ /g;

  # Removing leading and trailling spaces.

  $line =~ s/(^\s*)||(\s*$)//g;

  # Forget about comment lines

  if ($line =~ /^\s*#/) { next; }

  @line = split(/ /, $line);

  if ($isFirstLine)

  {

    if (@line != 1)

    {

      print STDERR "\n\nInvalid number (", (@line * 1), ") of digits data.\n";

      print STDERR "Aborting the program!\n\n";

      exit($EX_DATAERR);

    }

    $numDigitsData  = $line;

    $isFirstLine    = 0;

  }

  else

  {

    if (@line != $numDigitsData)

    {

      print STDERR "\n\nThe number of announced digits data ($numDigitsData) and ";

      print STDERR "the number of actual digits data", (@line * 1), " do not match.\n";

      print STDERR "Aborting the program!\n\n";

      exit($EX_DATAERR);

    }

    @digits = @line;

    last;

  }

} # End while ($line = <DF>)

close(DF);

#

# Scan the solutions file.

#

while ($line = <SVF>)

{

  # Remove the '\n'.

  chomp($line);

  # Coalescing all consecutive spaces to one.

  $line =~ s/(\s)+/ /g;

  # Removing leading and trailling spaces.

  $line =~ s/(^\s*)||(\s*$)//g;

  # Only take into account the lines starting with "A[".

  if (! ($line =~ /^A\[/))

  {

    next;

  }

  # Remove syntatic sugar.

  $line =~ s/.*=\[//;

  $line =~ s/\]://;

  $solutionVectorRef = [];

  @$solutionVectorRef = split(/,/, $line);

  for (my $i = 0 ; $i < @$solutionVectorRef ; $i++)

  {

    $$solutionVectorRef[$i] = Math::BigInt->new($$solutionVectorRef[$i]);

    if($$solutionVectorRef[$i]->is_nan())

    {

      print STDERR "\n\n\"$line\" from the solution vector file holds an invalid";

      print STDERR " number at index $i.";

      print STDERR " Aborting the program!\n\n";

      exit($EX_DATAERR);

    }

  } # End for

  push(@solutionMatrix, $solutionVectorRef);

} # End while ($line = <SVF>)

#

# Degree of the polynomial

#

$polynomialDegree = @{$solutionMatrix[0]} - 1;

if ($polynomialDegree < 0)

{

  print STDERR "\n\n\"Polynomial degree ($polynomialDegree) is invalid:";

  print STDERR " it should be superior or equal to 0.";

  print STDERR " Aborting the program!\n\n";

  exit($EX_DATAERR);

}

#

# Compute the chunk multipliers

#

#

# Shape a polynomial for Maple

#

for (my $i = 0 ; $i < @solutionMatrix ; $i++)

{

  $polynomialString = "P:=$variableName->";

  for (my $j = 0 ; $j <= $polynomialDegree ; $j++)

  {

    if ($j != 0)

    {

      $polynomialString .= " + $variableName * ("

    }

    if (${$solutionMatrix[$i]}[$j] < 0)

    {

      $polynomialString .= "((${$solutionMatrix[$i]}[$j])";

    }

    else

    {

      $polynomialString .= "(${$solutionMatrix[$i]}[$j]";

    }

      #$polynomialString .= "* 2^$decPointPos)";

      $polynomialString .= "/2^$digits[$j]";

      #$polynomialString .= "+";

    if ($j != 0)

    {

      $polynomialString .= ") ";

    }

  } # End or $j

  $polynomialString .= (')' x ($polynomialDegree + 1));

  $polynomialString .= ";";

  push(@polynomials, $polynomialString);

} # End for $i

print STDERR $polynomialString, "\n";

#!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

#exit(0);

#!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

#

if ($debug)

{

  print STDERR "\nNum digits: $numDigitsData\n";

  for (my $i = 0 ; $i < @digits ; $i++)

  {

    print STDERR "$digits[$i] ";

  }

  print SDTERR "\n";

  print STDERR "\n Solutions:\n";

  for (my $i = 0 ; $i < @solutionMatrix ; $i++)

  {

    for (my $j = 0 ; $j < @{$solutionMatrix[$i]} ; $j++)

    {

      print STDERR "${$solutionMatrix[$i]}[$j] ";

    } # End for $j

    print STDERR "\n"

  } # End for $i

  print STDERR "\n";

  print STDERR "\n";

  print STDERR "Degree of the polynomial: $polynomialDegree";

  print STDERR "\n";

  print STDERR "Polynomials:\n";

  foreach $elem (@polynomials)

  {

    print STDERR "$elem\n"

  }

  print STDERR "\n\n";

} # End debug

#

# Run Maple and analyse the output.

#

#

# The interesting part is made of the three last lines.

# The last line is made of a simple prompt at the begining of the line: ">"

# 

for (my $i = 0 ; $i < @polynomials ; $i++)

{

  $mapleProgram = "with(numapprox);\n";

  $mapleProgram .= "Digits:=80;\n";

  $mapleProgram .= "interface(screenwidth=1000);";

  $mapleProgram .= $polynomials[$i];

  $mapleProgram .= "infnorm(($function-P(x)), x=$lowerBound..$upperBound, 'xmax');";

  $mapleResult = &runMaple($mapleProgram);

  for (my $j = 0 ; $j <  @$mapleResult; $j++)

  {

    if ($$mapleResult[$j] =~ /^> quit/)

    {

      #

      # Removing leading and trailling spaces for the exponent line. 

      $$mapleResult[$j - 3] =~ s/(^\s*)||(\s*$)//g;

      # Check if the exponent line is empty and, if not, that it is made of 

      # digits

      if ($$mapleResult[$j - 3] ne "" && $$mapleResult[$j - 3] !~ /[:alpha:]/)

      {

        $infNorm = $$mapleResult[$j - 3];

        $powerOfTen = '(' . $$mapleResult[$j - 3] . '))';

      }

      else

      {

        $infNorm = "";

        $powerOfTen = "";

      }

      # Removing leading and trailling spaces.

      $$mapleResult[$j - 2] =~ s/(^\s*)||(\s*$)//g;

      print STDERR '.';

      if (length($powerOfTen) > 0)

      {

        @infNorm = split(/ /, $$mapleResult[$j - 2]);

        $infNorm = sprintf("%.30f", 

	                eval("$infNorm[0] * ($infNorm[1]**$powerOfTen"));

      }

      else

      {

        $infNorm = sprintf("%.30f",($$mapleResult[$j - 2] * 1));

      }

      while(defined($sortHash{$infNorm}))

      {

        $infNorm .= ' ';

      }

      $sortHash{$infNorm} = $solutionMatrix[$i];

    } # End if ($$mapleResult[$j] =~ /^> quit/)

  } # End for $j

} # End for $i

print STDERR "\n";

@sortedInfNorms = sort {$a <=> $b} (keys(%sortHash));

foreach $elem (@sortedInfNorms)

{

  for (my $i = 0 ; $i < @{$sortHash{$elem}} ; $i ++)

  {

    print  "${$sortHash{$elem}}[$i] ";

  }

  print "-\> ";

  print sprintf("%.10e",eval($elem)),"\n";

}

# Uncomment this line and remove the next one if using stdExitErrors.

exit($EX_OK);

##############################################################################

# PROGRAM SUBS

##############################################################################

#

##############################################################################

# sub runMaple

##############################################################################

# Performed task: runs Maple with a program and returns the answer

#

# input   : $mapleProgram

# output  : $mapleResultRef

# globals : none

# uses    : none

# notes   : none

#

sub runMaple {

  my ($mapleProgram) = @_;

  if (! defined($mapleProgram))

  {

    my @call_info = caller(0);

    print "\n\n", $call_info[3], " missing parameter. Aborting program!\n\n";

    # Uncomment this line and remove the next one if using stdExitErrors.

    #exit($EX_SOFTWARE);

    exit(1);

  }

  my $line        = "";

  my $mapleResult = [];

  #print STDERR $mapleProgram, "\n";

  open(MR, "echo \"$mapleProgram\"  | $mapleApplication |");

  while($line = <MR>)

  {

    chomp($line);

    push(@$mapleResult,$line); 

  }

  return($mapleResult);

} # End runMaple

##############################################################################

# sub mySub

##############################################################################

# Performed task:

#

# input   :

# output  :

# globals :

# uses    :

# notes   :

#

sub mySub {

  my ($someParameter) = @_;

  if (! defined($someParameter))

  {

    my @call_info = caller(0);

    print "\n\n", $call_info[3], " missing parameter. Aborting program!\n\n";

    # Uncomment this line and remove the next one if using stdExitErrors.

    #exit($EX_SOFTWARE);

    exit(1);

  }

} # End mySub

##############################################################################

# sub $retv = _privateSubName($argv)

##############################################################################

sub _privateSubName

{

	my $argv = shift();

	my $retv = undef;

	return($retv);

};

#

##############################################################################

# sub usage

##############################################################################

# Performed task: checks the command line parameters

#

# input   : none

# output  : none

# globals : $ARGV,others...

# uses    : none

# notes   :

#

sub usage {

  if (! defined($ARGV[4]))

  {

    my $scriptName = `basename $0`;

    chomp $scriptName;

    print STDERR "\n\nUsage: $scriptName numDigits solutionsFile \n";

    print STDERR "       " , (" " x length("$scriptName "));

    print STDERR " function lowerBound upperBound\n\n";

    print STDERR "  - digitsFile: the name of the file holding the number of digits\n";

    print STDERR "    for each coefficient;\n";

    print STDERR "  - solutionFile: the name of the solution vector file;\n";

    print STDERR "  - function (e.g. \"cos(x)\", don't forget quoting);\n";

    print STDERR "  - lowerBound (e.g. 0, don't forget quoting for complex expressions;\n";

    print STDERR "  - upperBound (e.g. \"Pi/4\" , don't forget quoting for complex expressions;\n";

    # Uncomment this line and remove the next one if using stdExitErrors.

    #exit($EX_USAGE);

    exit(1);

  } # End if

# Some initializations  

$digitsFileName           = $ARGV[0];

$digitsFileName           =~ s/"//g;

$solutionsFileName        = $ARGV[1];

$solutionsFileName        =~ s/"//g;

$function                 = $ARGV[2];

$function                 =~ s/"//g;

$lowerBound               = $ARGV[3];

$lowerBound               =~ s/"//g;

$upperBound               = $ARGV[4];

$upperBound               =~ s/"//g;

} # End usage

__END__