[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 6185 - (download) (as text) (annotate)
Sat Jan 23 16:39:35 2010 UTC (9 years, 11 months ago) by dpvc
File size: 13752 byte(s)
Fixed wrong default value in comments

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

aubreyja at gmail dot com
ViewVC Help
Powered by ViewVC 1.0.9