[system] / trunk / pg / lib / Value / Formula.pm Repository: Repository Listing bbplugincoursesdistsnplrochestersystemwww

# Diff of /trunk/pg/lib/Value/Formula.pm

Revision 5516 Revision 5517
179 # 179 #
181 # Get the tolerances, and check each adapted point relative 181 # Get the tolerances, and check each adapted point relative
182 # to the ORIGINAL correct answer. (This will have to be 182 # to the ORIGINAL correct answer. (This will have to be
183 # fixed if we ever do adaptive parameters for non-real formulas) 183 # fixed if we ever do adaptive parameters for non-real formulas)
184 #
185 # FIXME: it doesn't make sense to apply the ORIGINAL value's
186 # tolerance, and causes problems when the values
187 # differ in magnitude by much. Gavin has found several
188 # situations where this is a problem.
184 # 189 #
187 my \$tolerance = \$self->getFlag('tolerance',1E-4); 192 my \$tolerance = \$self->getFlag('tolerance',1E-4);
188 my \$isRelative = \$self->getFlag('tolType','relative') eq 'relative'; 193 my \$isRelative = \$self->getFlag('tolType','relative') eq 'relative';
410# 415#
411# Find adaptive parameters, if any 416# Find adaptive parameters, if any
412# 417#
414 my \$l = shift; my \$r = shift; 419 my \$l = shift; my \$r = shift;
415 my @params = @_; my \$d = scalar(@params); 420 my @params = @_; my \$d = scalar(@params); my \$D;
416 return 0 if \$d == 0; return 0 unless \$l->usesOneOf(@params); 421 return 0 if \$d == 0; return 0 unless \$l->usesOneOf(@params);
417 \$l->Error("Adaptive parameters can only be used for real-valued formulas") 422 \$l->Error("Adaptive parameters can only be used for real-valued formulas")
418 unless \$l->{tree}->isRealNumber; 423 unless \$l->{tree}->isRealNumber;
424
425 #
426 # Try up to three times (the random points might not work the first time)
427 #
428 foreach my \$attempt (1..3) {
419 # 429 #
420 # Get coefficient matrix of adaptive parameters 430 # Get coefficient matrix of adaptive parameters
421 # and value vector for linear system 431 # and value vector for linear system
422 # 432 #
423 my (\$p,\$v) = \$l->createRandomPoints(\$d); 433 my (\$p,\$v) = \$l->createRandomPoints(\$d);
424 my @P = (0) x \$d; my (\$f,\$F) = (\$l->{f},\$r->{f}); 434 my @P = (0) x \$d; my (\$f,\$F) = (\$l->{f},\$r->{f});
425 my @A = (); my @b = (); 435 my @A = (); my @b = ();
426 foreach my \$i (0..\$d-1) { 436 foreach my \$i (0..\$d-1) {
427 my @a = (); my @p = @{\$p->[\$i]}; 437 my @a = (); my @p = @{\$p->[\$i]};
428 foreach my \$j (0..\$d-1) { 438 foreach my \$j (0..\$d-1) {
429 \$P[\$j] = 1; push(@a,(&\$f(@p,@P)-\$v->[\$i])->value); 439 \$P[\$j] = 1; push(@a,(&\$f(@p,@P)-\$v->[\$i])->value);
430 \$P[\$j] = 0; 440 \$P[\$j] = 0;
431 }
432 push @A, [@a]; push @b, [(&\$F(@p,@P)-\$v->[\$i])->value];
433 }
434 #
435 # Use MatrixReal1.pm to solve system of linear equations
436 #
437 my \$M = MatrixReal1->new(\$d,\$d); \$M->[0] = \@A;
438 my \$B = MatrixReal1->new(\$d,1); \$B->[0] = \@b;
439 (\$M,\$B) = \$M->normalize(\$B);
440 \$M = \$M->decompose_LR;
441 if ((\$d,\$B,\$M) = \$M->solve_LR(\$B)) {
442 if (\$d == 0) {
443 #
444 # Get parameter values and recompute the points using them
445 #
446 my @a; my \$i = 0; my \$max = \$l->getFlag('max_adapt',1E8);
447 foreach my \$row (@{\$B->[0]}) {
448 if (abs(\$row->[0]) > \$max) {
449 \$max = Value::makeValue(\$max); \$row->[0] = Value::makeValue(\$row->[0]);
450 \$l->Error(["Constant of integration is too large: %s\n(maximum allowed is %s)",
451 \$row->[0]->string,\$max->string]) if \$params[\$i] eq 'C0';
452 \$l->Error(["Adaptive constant is too large: %s = %s\n(maximum allowed is %s)",
453 \$params[\$i],\$row->[0]->string,\$max->string]);
454 }
455 push @a, \$row->[0]; \$i++;
456 } 441 }
442 push @A, [@a]; push @b, [(&\$F(@p,@P)-\$v->[\$i])->value];
443 }
444 #
445 # Use MatrixReal1.pm to solve system of linear equations
446 #
447 my \$M = MatrixReal1->new(\$d,\$d); \$M->[0] = \@A;
448 my \$B = MatrixReal1->new(\$d,1); \$B->[0] = \@b;
449 (\$M,\$B) = \$M->normalize(\$B);
450 \$M = \$M->decompose_LR;
451 if ((\$D,\$B,\$M) = \$M->solve_LR(\$B)) {
452 if (\$D == 0) {
453 #
454 # Get parameter values and recompute the points using them
455 #
456 my @a; my \$i = 0; my \$max = \$l->getFlag('max_adapt',1E8);
457 foreach my \$row (@{\$B->[0]}) {
458 if (abs(\$row->[0]) > \$max) {
459 \$max = Value::makeValue(\$max); \$row->[0] = Value::makeValue(\$row->[0]);
460 \$l->Error(["Constant of integration is too large: %s\n(maximum allowed is %s)",
461 \$row->[0]->string,\$max->string]) if \$params[\$i] eq 'C0' or \$params[\$i] eq 'n00';
462 \$l->Error(["Adaptive constant is too large: %s = %s\n(maximum allowed is %s)",
463 \$params[\$i],\$row->[0]->string,\$max->string]);
464 }
465 push @a, \$row->[0]; \$i++;
466 }
457 my \$context = \$l->context; 467 my \$context = \$l->context;
458 foreach my \$i (0..\$#a) {\$context->{variables}{\$params[\$i]}{value} = \$a[\$i]} 468 foreach my \$i (0..\$#a) {\$context->{variables}{\$params[\$i]}{value} = \$a[\$i]}
459 \$l->{parameters} = [@a]; 469 \$l->{parameters} = [@a];