[system] / trunk / pg / macros / parserImplicitEquation.pl Repository:
ViewVC logotype

View of /trunk/pg/macros/parserImplicitEquation.pl

Parent Directory Parent Directory | Revision Log Revision Log


Revision 5392 - (download) (as text) (annotate)
Sun Aug 19 21:50:23 2007 UTC (12 years, 6 months ago) by dpvc
File size: 13489 byte(s)
Move the context creation to the _init routine so that modifications
to main:: occur at initialization time, not when the file is loaded.

    1 loadMacros("MathObjects.pl");
    2 
    3 sub _parserImplicitEquation_init {ImplicitEquation::Init()}; # don't reload this file
    4 
    5 =head1 DESCRIPTION
    6 
    7  ######################################################################
    8  #
    9  #  This is a MathObject class that implements an answer checker for
   10  #  implicitly defined equations.  The checker looks for the zeros of
   11  #  the equation and tests that the student and professor equations
   12  #  both have the same solutions.
   13  #
   14  #  This type of check is very subtle, and there are important issues
   15  #  that you may have to take into account.  The solutions to the
   16  #  equations are found numerically, and so they will not be exact;
   17  #  that means that there are tolerances that may need to be adjusted
   18  #  for your particular equation.  Also, it is always possible for the
   19  #  student to represent the function in a form that will exceed those
   20  #  tolerances, and so be marked as incorrect.  The answer checker
   21  #  attempts to set the parameters based on the values of the functions
   22  #  involved, but this may not work perfectly.
   23  #
   24  #  The method used to locate the solutions of A=B is by finding zeros
   25  #  of A-B, and it requires this function to take on both positive and
   26  #  negative values (that is, it can only find transverse intersections
   27  #  of the surface A-B=0 and the plane at height 0).  For example, even
   28  #  though the solutions of (x^2+y^1-1)^2=0 form a circle, the
   29  #  algorithm will fail to find any solutions for this equation.
   30  #
   31  #  In order to locate the zeros, you may need to change the limits so
   32  #  that they include regions where the function is both positive and
   33  #  negative (see below).  The algorithm will avoid discontinuities, so
   34  #  you can specify things like x-y=1/(x+y) rather than x^2-y^2=1.
   35  #
   36  #  Because the solutions are found using a random search, it is
   37  #  possible the randomly chosen starting points don't locate enough
   38  #  zeros for the function, in which case the check will fail.  This
   39  #  can happen for both the professor's function and the student's,
   40  #  since zeros are found for both.  This means that a correct answer
   41  #  can sometimes be marked incorrect depending on the random points
   42  #  chosen initially.  These points also affect the values selected for
   43  #  the tolerances used to determine when a function's value is zero,
   44  #  and so can affect whether the student's function is marked as
   45  #  correct or not.
   46  #
   47  #  If an equation has several components or branches, it is possible
   48  #  that the random location of solutions will not find zeros on some
   49  #  of the branches, and so might incorrectly mark as correct an
   50  #  equation that only is zero on one of the components.  For example,
   51  #  x^2-y^2=0 has solutions along the lines y=x and y=-x, so it is
   52  #  possible that x-y=0 or x+y=0 will be marked as correct if the
   53  #  random points are unluckily chosen.  One way to reduce this problem
   54  #  is to increase the number of solutions that are required (by
   55  #  setting the ImplicitPoints flag in the Context).  Another is to
   56  #  specify the solutions yourself, so that you are sure there are
   57  #  points on each component.
   58  #
   59  #  These problems should be rare, and the values for the various
   60  #  parameters have been set in an attempt to minimize the possibility
   61  #  of these errors, but they can occur, and you should be aware of
   62  #  them, and their possible solutions.
   63  #
   64  #
   65  #  Usage examples:
   66  #
   67  #     Context("ImplicitEquation");
   68  #     $f = ImplicitEquation("x^2 = cos(y)");
   69  #     $f = ImplicitEquation("x^2 - 2y^2 = 5",limits=>[[-3,3],[-2,2]]);
   70  #     $f = ImplicitEquation("x=1/y",tolerance=>.0001);
   71  #
   72  #  Then use
   73  #
   74  #     ANS($f->cmp);
   75  #
   76  #  to get the answer checker for $f.
   77  #
   78  #  There are a number of Context flags that control the answer checker.
   79  #  These include:
   80  #
   81  #    ImplicitPoints => 7                   (the number of solutions to test)
   82  #    ImplicitTolerance => 1E-6             (relative tolerance value for when
   83  #                                           the tested function is zero)
   84  #    ImplicitAbsoluteMinTolerance => 1E-3  (the minimum tolerance allowed)
   85  #    ImplicitAbsoluteMaxTolerance => 1E-3  (the maximum tolerance allowed)
   86  #    ImplicitPointTolerance => 1E-9        (relative tolerance for how close
   87  #                                           the solution point must be to an
   88  #                                           actual solution)
   89  #    BisectionTolerance => .01             (extra factor used for the tolerance
   90  #                                           when finding the solutions)
   91  #    BisectionCutoff => 40                 (maximum number of bisections to
   92  #                                           perform when looking for a solution)
   93  #
   94  #  You may set any of these using Context()->flags->set(...).
   95  #
   96  #  In addition to the Context flags, you can set some values within
   97  #  the ImplicitEquation itself:
   98  #
   99  #    tolerance                (the absolute tolerance for zeros of the function)
  100  #    bisect_tolerance         (the tolerance used when searching for zeros)
  101  #    point_tolerance          (the absolute tolerance for how close to an
  102  #                              actual solution the located solution must be)
  103  #    limits                   (the domain to use for the function; see the
  104  #                              documentation for the Formula object)
  105  #    solutions                (a reference to an array of references to arrays
  106  #                              that contain the coordinates of the points
  107  #                              that are the solutions of the equation)
  108  #
  109  #  These can be set in the in the ImplicitEquation() call that creates
  110  #  the object, as in the examples below:
  111  #
  112  #  For example:
  113  #
  114  #    $f = ImplicitEquation("x^2-y^2=0",
  115  #            solutions => [[0,0],[1,1],[-1,1],[-1,-1],[1,-1]],
  116  #            tolerance => .001
  117  #         );
  118  #
  119  #
  120  #    $f = ImplicitEquation("xy=5",limits=>[-3,3]);
  121  #
  122  #  The limits value can be set globally within the Context, if you wish,
  123  #  and the others can be controlled by the Context flags discussed
  124  #  above.
  125  #
  126  ######################################################################
  127 
  128 =cut
  129 
  130 #
  131 #  Create the ImplicitEquation package
  132 #
  133 package ImplicitEquation;
  134 our @ISA = qw(Value::Formula);
  135 
  136 sub Init {
  137   my $context = $main::context{ImplicitEquation} = Parser::Context->getCopy("Numeric");
  138   $context->variables->are(x=>'Real',y=>'Real');
  139   $context{precedence}{ImplicitEquation} = $context->{precedence}{special};
  140   Parser::BOP::equality->Allow($context);
  141   $context->flags->set(
  142     ImplicitPoints => 10,
  143     ImplicitTolerance => 1E-6,
  144     ImplicitAbsoluteMinTolerance => 1E-3,
  145     ImplicitAbsoluteMaxTolerance => 1,
  146     ImplicitPointTolerance => 1E-9,
  147     BisectionTolerance => .01,
  148     BisectionCutoff => 40,
  149   );
  150 
  151   main::Context("ImplicitEquation");  ### FIXEME:  probably should require author to set this explicitly
  152 
  153   main::PG_restricted_eval('sub ImplicitEquation {ImplicitEquation->new(@_)}');
  154 }
  155 
  156 sub new {
  157   my $self = shift; my $class = ref($self) || $self;
  158   my $context = (Value::isContext($_[0]) ? shift : $self->context);
  159   my $f = shift; return $f if ref($f) eq $class;
  160   $f = main::Formula($f);
  161   Value::Error("Your formula doesn't look like an implicit equation")
  162     unless $f->type eq 'Equality';
  163   my $F = ($context->Package("Formula")->new($context,$f->{tree}{lop}) -
  164      $context->Package("Formula")->new($context,$f->{tree}{rop}))->reduce;
  165   $F = $context->Package("Formula")->new($F) unless Value::isFormula($F);
  166   Value::Error("Your equation must be real-valued") unless $F->isRealNumber;
  167   Value::Error("Your equation should not be constant") if $F->isConstant;
  168   Value::Error("Your equation can not contain adaptive parameters")
  169     if ($F->usesOneOf($context->variables->parameters));
  170   $F = bless $F, $class;
  171   my %options = (@_);  # user can supply limits, tolerance, etc.
  172   foreach my $id (keys %options) {$F->{$id} = $options{$id}}
  173   $F->{F} = $f; $F->{isValue} = $F->{isFormula} = 1;
  174   $F->createPoints unless $F->{solutions};
  175   return $F;
  176 }
  177 
  178 #
  179 #  Override the comparison method.
  180 #
  181 #  Turn the right-hand equation into an ImplicitEquation.  This creates
  182 #  the test points (i.e., finds the solution points).  Then check
  183 #  the professor's function on the student's test points and the
  184 #  student's function on the professor's test points.
  185 #
  186 
  187 sub compare {
  188   my ($l,$r) = @_; my $self = $l; my $tolerance;
  189   my @params; @params = (limits=>$l->{limits}) if $l->{limits};
  190   $r = ImplicitEquation->new($r,@params);
  191   Value::Error("Functions from different contexts can't be compared")
  192     unless $l->{context} == $r->{context};
  193 
  194   #
  195   #  They are not equal if couldn't get solutions for one of them
  196   #
  197   return 1 unless $l->{solutions} && $r->{solutions};
  198 
  199   #
  200   #  Test the right-hand function on the solutions of the left-hand one
  201   #  and vice-versa
  202   #
  203   my $rzeros = $r->createPointValues($l->{solutions});
  204   my $lzeros = $l->createPointValues($r->{solutions});
  205   return 1 unless $lzeros && $rzeros;
  206 
  207   #
  208   #  Check that the values are, in fact, zeros
  209   #
  210   $tolerance = $r->getFlag('tolerance',1E-3);
  211   foreach my $v (@{$rzeros}) {return 1 unless abs($v) < $tolerance}
  212   $tolerance = $l->getFlag('tolerance',1E-3);
  213   foreach my $v (@{$lzeros}) {return 1 unless abs($v) < $tolerance}
  214 
  215   return 0; # equal
  216 }
  217 
  218 #
  219 #  Use the original equation for these (but not for perl(), since we
  220 #  need that to use perlFunction).
  221 #
  222 sub string {shift->{F}->string}
  223 sub TeX    {shift->{F}->TeX}
  224 
  225 sub cmp_class {'an Implicit Equation'}
  226 sub showClass {shift->cmp_class}
  227 
  228 #
  229 #  Locate points that satisfy the equation
  230 #
  231 sub createPoints {
  232   my $self = shift;
  233   my $num_points = int($self->getFlag('ImplicitPoints',10));
  234   $num_points = 1 if $num_points < 1;
  235 
  236   #
  237   #  Get some positive and negative test points (try up to 5 times)
  238   #
  239   my ($OK,$p,$n,$z,$f,@zero); my $k = 5;
  240   while (!$OK && --$k) {($OK,$p,$n,$z,$f) = $self->getPositiveNegativeZero($num_points)}
  241   Value::Error("Can't find any solutions to your equation") unless $OK;
  242   my ($P,@intervals) = $self->getIntervals($p,$n);
  243 
  244   #
  245   #  Get relative tolerance values and make them absolute
  246   #
  247   my $minTolerance = $self->getFlag('ImplicitAbsoluteMinTolerance');
  248   my $maxTolerance = $self->getFlag('ImplicitAbsoluteMaxTolerance');
  249   my $tolerance = $f * $self->getFlag('ImplicitTolerance');
  250   $tolerance = $minTolerance if $tolerance < $minTolerance;
  251   $tolerance = $maxTolerance if $tolerance > $maxTolerance;
  252   $self->{tolerance} = $tolerance unless $self->{tolerance};
  253   $self->{bisect_tolerance} = $self->{tolerance} * $self->getFlag('BisectionTolerance')
  254     unless $self->{bisect_tolerance};
  255   $self->{point_tolerance} = $P * $self->getFlag('ImplicitPointTolerance')
  256     unless $self->{point_tolerance};
  257 
  258   #
  259   #  Locate solutions to be used for comparison test
  260   #
  261   @zero = @{$z}; @zero = $zero[0..$num_points-1] if (scalar(@zero) > $num_points);
  262   for ($i = 0; scalar(@zero) < $num_points && $i < scalar(@intervals); $i++) {
  263     my $Q = $self->Bisect($intervals[$i][0],$intervals[$i][1]);
  264     push(@zero,$Q) if $Q;
  265   }
  266   Value::Error("Can't find enough solutions for an effective test")
  267     unless scalar(@zero) == $num_points;
  268 
  269   #
  270   #  Save the solutions to the equation
  271   #
  272   $self->{solutions} = [@zero];
  273 }
  274 
  275 #
  276 #  Get random points and sort them by sign of the function.
  277 #  Also, get the average function value, and indicate if
  278 #  we actually did find both positive and negative values.
  279 #
  280 sub getPositiveNegativeZero {
  281   my $self = shift; my $n = shift;
  282   my ($p,$v) = $self->SUPER::createRandomPoints(3*$n);
  283   my (@pos,@neg,@zero);
  284   my $f = 0; my $k = 0;
  285   foreach my $i (0..scalar(@{$v})-1) {
  286     if ($v->[$i] == 0) {push(@zero,$p->[$i])} else {
  287       $f += abs($v->[$i]); $k++;
  288       if ($v->[$i] > 0) {push(@pos,$p->[$i])} else {push(@neg,$p->[$i])}
  289     }
  290   }
  291   $f /= $k if $k;
  292   return (scalar(@pos) && scalar(@neg),[@pos],[@neg],[@zero],$f);
  293 }
  294 
  295 #
  296 #  Make a list of positive and negative points sorted by
  297 #  the distance between them.  Also return the average distance
  298 #  between points.
  299 #
  300 sub getIntervals {
  301   my $self = shift; my $pos = shift; my $neg = shift;
  302   my @intervals = (); my $D = 0;
  303   my $point = $self->Package("Point");
  304   my $context = $self->context;
  305   foreach my $p (@{$pos}) {
  306     foreach my $n (@{$neg}) {
  307       my $d = abs($point->make($context,@{$p}) - $point->make($context,@{$n}));
  308       push(@intervals,[$p,$n,$d]); $D += $d;
  309     }
  310   }
  311   @intervals = main::PGsort(sub {$_[0]->[2] < $_[1]->[2]},@intervals);
  312   return($D/scalar(@intervals),@intervals);
  313 }
  314 
  315 #
  316 #  Use bisection algorithm to find a point where the function is zero
  317 #  (i.e., where the original equation is an equality)
  318 #  If we can't find a point (e.g., we are at a discontinuity),
  319 #  return an undefined value.
  320 #
  321 sub Bisect {
  322   my $self = shift;
  323   my $tolerance  = $self->getFlag('bisect_tolerance',1E-5);
  324   my $ptolerance = $self->getFlag('point_tolerance',1E-9);
  325   my $m = $self->getFlag('BisectionCutoff',30); my ($P,$f);
  326   my $point = $self->Package("Point"); my $context = $self->context;
  327   my $P0 = $point->make($context,@{$_[0]}); my $P1 = $point->make($context,@{$_[1]});
  328   my ($f0,$f1) = @{$self->createPointValues([$P0->data,$P1->data],1)};
  329   for (my $i = 0; $i < $m; $i++) {
  330     $P = ($P0+$P1)/2; $f = $self->createPointValues([$P->data]);
  331     return unless ref($f);
  332     $f = $f->[0];
  333     return [$P->value] if abs($f) < $tolerance && abs($P1-$P0) < $ptolerance;
  334     if ($f > 0) {$P0 = $P; $f0 = $f} else {$P1 = $P; $f1 = $f}
  335   }
  336   return;
  337 }
  338 
  339 1;

aubreyja at gmail dot com
ViewVC Help
Powered by ViewVC 1.0.9