Statistiques
| Branche: | Révision :

meplib / binaries / bin / checkSimpleSolutions-01.pl @ master

Historique | Voir | Annoter | Télécharger (13,25 ko)

1
#! /usr/bin/perl
2

    
3
# checkSimpleSolutions-01.pl - ST - 2006-01-24
4
#
5
# Check the quality of a solution with the Maple infnorm.
6
#
7
# COMMAND LINE SYNTAX:
8
# checkSimpleSolutions-01.pl  digitsNum decPointPos solutionsFile positionsFile
9

    
10
# NOTE(S):
11
# This version checks the solutions produced by the "gp" PARI utility
12
# The number of chunks is given by the positionsFile
13

    
14
##############################################################################
15
# PROGRAM DECLARATIONS
16
##############################################################################
17

    
18
BEGIN
19
{
20
  # Modify @INC here if necessary
21
  push(@INC, '../utils');
22
};
23

    
24
# core modules
25
use strict; # must be first to force strict in all others
26
use Carp;
27
use Getopt::Long;
28
use Math::BigInt;
29
use Math::BigRat;
30
use Math::BigFloat;
31

    
32
# custom modules
33
# Uncomment if you want to use standard exit error codes.
34
# the stdExitErrors.pm file must be available around.
35
use stdExitErrors;
36

    
37
##############################################################################
38
# PROGRAM CONSTANT DECLARATIONS
39
##############################################################################
40
my $bigIntZero = Math::BigInt->new('0');
41
my $debug = 1;
42
my $variableName = "x"; 
43
my $mapleApplication = "/usr/local/maple/bin/maple";
44
##############################################################################
45
# PROGRAM VARIABLE DECLARATIONS
46
##############################################################################
47
my $solutionsFileName = "";
48
my $positionsFileName = "";
49
my $numDigits         = 0;
50
my $decPointPos       = 0;
51
my $numChunks         = 0;
52
my $function          = "";
53
my $lowerBound        = "";
54
my $upperBound        = "";
55
my $elem              = "";
56
my $line              = "";
57
my $lineNum           = 0;
58
my @line              = ();
59
my $powerOfTen        = 0;
60
my $solutionVectorRef = 0;
61
my @solutionMatrix    = ();
62
my @positionsVector   = ();
63
my $powerIndex        = 0;
64
my $polynomialDegree  = 0;
65
my $mainDivisor       = 0;
66
my @chunkMultipliers  = ();
67
my $polynomialString  = "";
68
my @polynomials       = ();
69
my $mapleProgram      = ();
70
my $mapleResult       = 0;
71
my $infNorm           = "";
72
my @infNorm           = ();
73
my %sortHash          = ();
74
my @sortedInfNorms    = ();
75
#
76
#
77
#
78
##############################################################################
79
# TEST ENVIRONMENT AND CONFIGURATION VARIABLES
80
##############################################################################
81

    
82
##############################################################################
83
# INITIALIZE DEBUG FILE
84
##############################################################################
85

    
86
##############################################################################
87
# INITIALIZE LOG FILE
88
##############################################################################
89

    
90
##############################################################################
91
# GET COMMAND LINE ARGUMENTS
92
##############################################################################
93

    
94
##############################################################################
95
# PROGRAM MAIN
96
##############################################################################
97
&usage();
98
#
99
# Check the command line parameters
100
#
101
# It's a TODO!
102
#
103
# Open the files
104
#
105
# First the solution vector file.
106
#
107
if (! open(SVF, "<$solutionsFileName"))
108
{
109
  print STDERR "\n\nCan't open $solutionsFileName. Aborting the program!\n\n";
110
  exit($EX_NOINPUT);
111
}
112
# Second the positions file.
113
#
114
if (! open(POF, "<$positionsFileName"))
115
{
116
  print STDERR "\n\nCan't open $positionsFileName. Aborting the program!\n\n";
117
  exit($EX_NOINPUT);
118
}
119
#
120
# Scan the solutions file.
121
#
122
while ($line = <SVF>)
123
{
124
  # Remove the '\n'.
125
  chomp($line);
126
  # Coalescing all consecutive spaces to one.
127
  $line =~ s/(\s)+/ /g;
128
  # Removing leading and trailling spaces.
129
  $line =~ s/(^\s*)||(\s*$)//g;
130
  # Only take into account the lines starting with "A[".
131
  if (! ($line =~ /^A\[/))
132
  {
133
    next;
134
  }
135
  # Remove syntatic sugar.
136
  $line =~ s/.*=\[//;
137
  $line =~ s/\]://;
138
  $solutionVectorRef = [];
139
  @$solutionVectorRef = split(/,/, $line);
140
  for (my $i = 0 ; $i < @$solutionVectorRef ; $i++)
141
  {
142
    $$solutionVectorRef[$i] = Math::BigInt->new($$solutionVectorRef[$i]);
143
    if($$solutionVectorRef[$i]->is_nan())
144
    {
145
      print STDERR "\n\n\"$line\" from the solution vector file holds an invalid";
146
      print STDERR " number at index $i.";
147
      print STDERR " Aborting the program!\n\n";
148
      exit($EX_DATAERR);
149
    }
150
  } # End for
151
  push(@solutionMatrix, $solutionVectorRef);
152
} # End while ($line = <SVF>)
153
#
154
# Read the positions file
155
#
156
$lineNum = 1;
157
while($line=<POF>)
158
{
159
  chomp($line);
160
  # Coalescing all consecutive spaces to one.
161
  $line =~ s/(\s)+/ /g;
162
  # Removing leading and trailling spaces.
163
  $line =~ s/(^\s*)||(\s*$)//g;
164
  # The number of chunks is all by itself on the first line.
165
  if ($lineNum == 1)
166
  {
167
    $numChunks = $line;
168
  }
169
  # All the element are on line 2
170
  if ($lineNum == 2)
171
  {
172
    @positionsVector = split(/ /, $line);
173
    if (@positionsVector != $numChunks)
174
    {
175
      print STDERR "\n\n\ The announced number of chunks($numChunks) ";
176
      print STDERR "does not match the actual number of chuncks ", @positionsVector * 1;
177
      print STDERR ". Aborting the program!\n\n";
178
      exit($EX_DATAERR);
179
    }
180
    last; # Only two data lines in positions file.
181
  }
182
  $lineNum++;
183
} # end scanning POF
184
#
185
# Degree of the polynomial
186
#
187
if ((@{$solutionMatrix[0]} % $numChunks) != 0)
188
{
189
  print STDERR "\n\n\ The number of elements per solution (",(@{$solutionMatrix[0]} * 1);
190
  print STDERR ") and the number of chunks ($numChunks) do not match: ";
191
  print STDERR "NumElemPerSols mod NumChunks should be zero.\n";
192
  print STDERR "Aborting the program!\n\n";
193
  exit($EX_DATAERR);
194
}
195
$polynomialDegree = @{$solutionMatrix[0]} / $numChunks - 1;
196
# TODO: add a check on the degree (>=0)
197
$mainDivisor = "2^$decPointPos";
198
#
199
# Compute the chunk multipliers
200
#
201
for (my $i = 0 ; $i < @positionsVector ; $i++)
202
{
203
  push(@chunkMultipliers, "2^$positionsVector[$i]");
204
} # End for $i
205
#
206
# Shape a polynomial for Maple
207
#
208
for (my $i = 0 ; $i < @solutionMatrix ; $i++)
209
{
210
  $polynomialString = "P:=$variableName->";
211
  for (my $j = 0 ; $j <= $polynomialDegree ; $j++)
212
  {
213
    if ($j != 0)
214
    {
215
      $polynomialString .= " + $variableName * (("
216
    }
217
    for (my $k = 0 ; $k < @chunkMultipliers ; $k++)
218
    {
219
      if (${$solutionMatrix[$i]}[$j] < 0)
220
      {
221
        $polynomialString .= "((${$solutionMatrix[$i]}[$j * $numChunks + $k])";
222
      }
223
      else
224
      {
225
        $polynomialString .= "(${$solutionMatrix[$i]}[$j * $numChunks + $k]";
226
      }
227
      $polynomialString .= "*$chunkMultipliers[$k])";
228
      $polynomialString .= "/$mainDivisor";
229
      if ($k != (@chunkMultipliers - 1))
230
      {
231
        $polynomialString .= "+";
232
      }
233
    } # End for $k $so
234
    if ($j != 0)
235
    {
236
      $polynomialString .= ") ";
237
    }
238
  } # End or $j
239
  $polynomialString .= (')' x $polynomialDegree);
240
  $polynomialString .= ";";
241
  push(@polynomials, $polynomialString);
242
} # End for $i
243
#!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
244
#exit(0);
245
#!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
246
#
247
if ($debug)
248
{
249
  print STDERR "\nNum digits: $numDigits\n";
250
  print STDERR "\nDecimal point position: $decPointPos\n";
251
  print STDERR "\n Solutions:\n";
252
  for (my $i = 0 ; $i < @solutionMatrix ; $i++)
253
  {
254
    for (my $j = 0 ; $j < @{$solutionMatrix[$i]} ; $j++)
255
    {
256
      print STDERR "${$solutionMatrix[$i]}[$j] ";
257
    } # End for $j
258
    print STDERR "\n"
259
  } # End for $i
260
  print STDERR "\n";
261
  print STDERR "$numChunks chunk(s) at position(s): ";
262
  for (my $i = 0 ; $i < @positionsVector ; $i++)
263
  {
264
    print STDERR "$positionsVector[$i] ";
265
  } # End for $i
266
  print STDERR "\n";
267
  print STDERR "Degree of the polynomial: $polynomialDegree";
268
  print STDERR "\n";
269
  print STDERR "Main divisor: $mainDivisor";
270
  print STDERR "\n";
271
  print STDERR "Chunk multipliers:";
272
  for (my $i = 0 ; $i < @chunkMultipliers ; $i++)
273
  {
274
    print STDERR "$chunkMultipliers[$i] ";
275
  } # End for $i
276
  print STDERR "Polynomials:\n";
277
  foreach $elem (@polynomials)
278
  {
279
    print STDERR "$elem\n"
280
  }
281
  print STDERR "\n\n";
282
} # End debug
283
#
284
# Run Maple and analyse the output.
285
#
286
#
287
# The interesting part is made of the three last lines.
288
# The last line is made of a simple prompt at the begining of the line: ">"
289
# 
290
for (my $i = 0 ; $i < @polynomials ; $i++)
291
{
292
  $mapleProgram = "with(numapprox);";
293
  $mapleProgram .= $polynomials[$i];
294
  $mapleProgram .= "infnorm(($function-P(x)), x=$lowerBound..$upperBound, 'xmax');";
295
  $mapleResult = &runMaple($mapleProgram);
296
  for (my $j = 0 ; $j <  @$mapleResult; $j++)
297
  {
298
    if ($$mapleResult[$j] =~ /^> quit/)
299
    {
300
      #
301
      # Removing leading and trailling spacesi for the exponent line. 
302
      $$mapleResult[$j - 3] =~ s/(^\s*)||(\s*$)//g;
303
      if ($$mapleResult[$j - 3] ne "")
304
      {
305
        $infNorm = $$mapleResult[$j - 3];
306
        $powerOfTen = '(' . $$mapleResult[$j - 3] . '))';
307
      }
308
      else
309
      {
310
        $infNorm = "";
311
        $powerOfTen = "";
312
      }
313
      # Removing leading and trailling spaces.
314
      $$mapleResult[$j - 2] =~ s/(^\s*)||(\s*$)//g;
315
      print STDERR '.';
316
      if (length($powerOfTen) > 0)
317
      {
318
        @infNorm = split(/ /, $$mapleResult[$j - 2]);
319
        $infNorm = eval("$infNorm[0] * ($infNorm[1]**".$powerOfTen);
320
      }
321
      else
322
      {
323
        $infNorm = $$mapleResult[$j - 2];
324
      }
325
      #print STDERR $infNorm;
326
      while(defined($sortHash{$infNorm}))
327
      {
328
        $infNorm .= ' ';
329
      }
330
      $sortHash{$infNorm} = $solutionMatrix[$i];
331
    } # End if ($$mapleResult[$j] =~ /^> quit/)
332
  } # End for $j
333
} # End for $i
334
print STDERR "\n";
335
@sortedInfNorms = sort(keys(%sortHash));
336
foreach $elem (@sortedInfNorms)
337
{
338
  for (my $i = 0 ; $i < @{$sortHash{$elem}} ; $i ++)
339
  {
340
    print  "${$sortHash{$elem}}[$i] ";
341
  }
342
  print "-\> ";
343
  print "$elem\n";
344
}
345

    
346
# Uncomment this line and remove the next one if using stdExitErrors.
347
exit($EX_OK);
348

    
349
##############################################################################
350
# PROGRAM SUBS
351
##############################################################################
352
#
353
##############################################################################
354
# sub runMaple
355
##############################################################################
356
# Performed task: runs Maple with a program and returns the answer
357
#
358
# input   : $mapleProgram
359
# output  : $mapleResultRef
360
# globals : none
361
# uses    : none
362
# notes   : none
363
#
364
sub runMaple {
365
  my ($mapleProgram) = @_;
366
  if (! defined($mapleProgram))
367
  {
368
    my @call_info = caller(0);
369
    print "\n\n", $call_info[3], " missing parameter. Aborting program!\n\n";
370
    # Uncomment this line and remove the next one if using stdExitErrors.
371
    #exit($EX_SOFTWARE);
372
    exit(1);
373
  }
374
  my $line        = "";
375
  my $mapleResult = [];
376
  #print STDERR $mapleProgram, "\n";
377
  open(MR, "echo \"$mapleProgram\"  | $mapleApplication |");
378
  while($line = <MR>)
379
  {
380
    chomp($line);
381
    push(@$mapleResult,$line); 
382
  }
383
  return($mapleResult);
384
} # End runMaple
385
##############################################################################
386
# sub mySub
387
##############################################################################
388
# Performed task:
389
#
390
# input   :
391
# output  :
392
# globals :
393
# uses    :
394
# notes   :
395
#
396
sub mySub {
397
  my ($someParameter) = @_;
398
  if (! defined($someParameter))
399
  {
400
    my @call_info = caller(0);
401
    print "\n\n", $call_info[3], " missing parameter. Aborting program!\n\n";
402
    # Uncomment this line and remove the next one if using stdExitErrors.
403
    #exit($EX_SOFTWARE);
404
    exit(1);
405
  }
406
} # End mySub
407
##############################################################################
408
# sub $retv = _privateSubName($argv)
409
##############################################################################
410

    
411
sub _privateSubName
412
{
413
	my $argv = shift();
414
	my $retv = undef;
415
	return($retv);
416
};
417
#
418
##############################################################################
419
# sub usage
420
##############################################################################
421
# Performed task: checks the command line parameters
422
#
423
# input   : none
424
# output  : none
425
# globals : $ARGV,others...
426
# uses    : none
427
# notes   :
428
#
429
sub usage {
430
  if (! defined($ARGV[6]))
431
  {
432
    my $scriptName = `basename $0`;
433
    chomp $scriptName;
434
    print STDERR "\n\nUsage: $scriptName numDigits decPointPos solutionsFile \n";
435
    print STDERR "       " , (" " x length("$scriptName ")), "positionsVectorFile";
436
    print STDERR " function lowerBound upperBound\n\n";
437
    print STDERR "  - numDigits: the number of digits of the coefficients;\n";
438
    print STDERR "  - decPointPos: the positions of the decimal point;\n";
439
    print STDERR "  - solutionFile: the name of the solution vector file;\n";
440
    print STDERR "  - positionsFile: the name of the positions (of chunks in\n"; 
441
    print STDERR "    coefficient) vector file;\n";
442
    print STDERR "  - function (e.g. \"cos(x)\", don't forget quoting);\n";
443
    print STDERR "  - lowerBound (e.g. 0, don't forget quoting for complex expressions;\n";
444
    print STDERR "  - upperBound (e.g. \"Pi/4\" , don't forget quoting for complex expressions;\n";
445
    # Uncomment this line and remove the next one if using stdExitErrors.
446
    #exit($EX_USAGE);
447
    exit(1);
448
  } # End if
449
# Some initializations  
450
$numDigits                = $ARGV[0];
451
$decPointPos              = $ARGV[1];
452
$solutionsFileName        = $ARGV[2];
453
$positionsFileName        = $ARGV[3];
454
$function                 = $ARGV[4];
455
$lowerBound               = $ARGV[5];
456
$upperBound               = $ARGV[6];
457
} # End usage
458

    
459

    
460
__END__
461