Parent Directory
| 
Revision Log
Added in functMaxConstantOfIntegration support in the adaptive
parameter answer checking. Also added a flag to the formula answer
checker that allows checking a formula up to a constant (for
integration problems). This really just adds a new parameter to the
context and adds that to the formula, so there is nothing deep about
this.
ANS(Formula("2x")->cmp(upToConstant=>1));
Finally, don't look for adaptive parameters if they aren't actually
used in the professor's formula (even if they are defined).
1 ########################################################################### 2 # 3 # Implements the Formula class. 4 # 5 package Value::Formula; 6 my $pkg = 'Value::Formula'; 7 8 use strict; 9 use vars qw(@ISA); 10 @ISA = qw(Parser Value); 11 12 use overload 13 '+' => \&add, 14 '-' => \&sub, 15 '*' => \&mult, 16 '/' => \&div, 17 '**' => \&power, 18 '.' => \&dot, 19 'x' => \&cross, 20 '<=>' => \&compare, 21 'cmp' => \&Value::cmp, 22 '~' => sub {Parser::Function->call('conj',$_[0])}, 23 'neg' => sub {$_[0]->neg}, 24 'sin' => sub {Parser::Function->call('sin',$_[0])}, 25 'cos' => sub {Parser::Function->call('cos',$_[0])}, 26 'exp' => sub {Parser::Function->call('exp',$_[0])}, 27 'abs' => sub {Parser::Function->call('abs',$_[0])}, 28 'log' => sub {Parser::Function->call('log',$_[0])}, 29 'sqrt' => sub {Parser::Function->call('sqrt',$_[0])}, 30 'atan2' => \&atan2, 31 'nomethod' => \&Value::nomethod, 32 '""' => \&Value::stringify; 33 34 # 35 # Call Parser to make the new item 36 # 37 sub new {shift; $pkg->SUPER::new(@_)} 38 39 # 40 # Create the new parser with no string 41 # (we'll fill in its tree by hand) 42 # 43 sub blank {$pkg->SUPER::new('')} 44 45 # 46 # Get the type from the tree 47 # 48 sub typeRef {(shift)->{tree}->typeRef} 49 50 ############################################ 51 # 52 # Create a BOP from two operands 53 # 54 # Get the context and variables from the left and right operands 55 # if they are formulas 56 # Make them into Value objects if they aren't already. 57 # Convert '+' to union for intervals or unions. 58 # Make a new BOP with the two operands. 59 # Record the variables. 60 # Evaluate the formula if it is constant. 61 # 62 sub bop { 63 my ($l,$r,$flag,$bop) = @_; 64 if ($l->promotePrecedence($r)) {return $r->add($l,!$flag)} 65 if ($flag) {my $tmp = $l; $l = $r; $r = $tmp} 66 my $formula = $pkg->blank; my $parser = $formula->{context}{parser}; 67 my $vars = {}; 68 if (ref($r) eq $pkg) { 69 $formula->{context} = $r->{context}; 70 $vars = {%{$vars},%{$r->{variables}}}; 71 $r = $r->{tree}->copy($formula); 72 } 73 if (ref($l) eq $pkg) { 74 $formula->{context} = $l->{context}; 75 $vars = {%{$vars},%{$l->{variables}}}; 76 $l = $l->{tree}->copy($formula); 77 } 78 $l = $pkg->new($l) if (!ref($l) && Value::getType($formula,$l) eq "unknown"); 79 $r = $pkg->new($r) if (!ref($r) && Value::getType($formula,$r) eq "unknown"); 80 $l = $parser->{Value}->new($formula,$l) unless ref($l) =~ m/^Parser::/; 81 $r = $parser->{Value}->new($formula,$r) unless ref($r) =~ m/^Parser::/; 82 $bop = 'U' if $bop eq '+' && 83 ($l->type =~ m/Interval|Union/ || $r->type =~ m/Interval|Union/); 84 $formula->{tree} = $parser->{BOP}->new($formula,$bop,$l,$r); 85 $formula->{variables} = {%{$vars}}; 86 return $formula->eval if scalar(%{$vars}) == 0; 87 return $formula; 88 } 89 90 sub add {bop(@_,'+')} 91 sub sub {bop(@_,'-')} 92 sub mult {bop(@_,'*')} 93 sub div {bop(@_,'/')} 94 sub power {bop(@_,'**')} 95 sub cross {bop(@_,'><')} 96 97 # 98 # Make dot work for vector operands 99 # 100 sub dot { 101 my ($l,$r,$flag) = @_; 102 if ($l->promotePrecedence($r)) {return $r->compare($l,!$flag)} 103 return bop(@_,'.') if $l->type eq 'Vector' && 104 Value::isValue($r) && $r->type eq 'Vector'; 105 Value::_dot(@_); 106 } 107 108 ############################################ 109 # 110 # Form the negation of a formula 111 # 112 sub neg { 113 my $self = shift; 114 my $formula = $self->blank; 115 $formula->{context} = $self->{context}; 116 $formula->{variables} = $self->{variables}; 117 $formula->{tree} = $formula->{context}{parser}{UOP}->new($formula,'u-',$self->{tree}->copy($formula)); 118 return $formula->eval if scalar(%{$formula->{variables}}) == 0; 119 return $formula; 120 } 121 122 # 123 # Form the function atan2 function call on two operands 124 # 125 sub atan2 { 126 my ($l,$r,$flag) = @_; 127 if ($l->promotePrecedence($r)) {return $r->compare($l,!$flag)} 128 if ($flag) {my $tmp = $l; $l = $r; $r = $tmp} 129 Parser::Function->call('atan2',$l,$r); 130 } 131 132 ############################################ 133 # 134 # Compare two functions for equality 135 # 136 sub compare { 137 my ($l,$r,$flag) = @_; my $self = $l; 138 if ($l->promotePrecedence($r)) {return $r->compare($l,!$flag)} 139 $r = Value::Formula->new($r) unless Value::isFormula($r); 140 Value::Error("Functions from different contexts can't be compared") 141 unless $l->{context} == $r->{context}; 142 143 # 144 # Get the test points and evaluate the functions at those points 145 # 146 ## FIXME: Check given points for consistency 147 my $points = $l->{test_points} || $r->{test_points} || $l->createRandomPoints; 148 my $lvalues = $l->{test_values} || $l->createPointValues($points,1); 149 my $rvalues = $r->createPointValues($points); 150 # 151 # Note: $l is bigger if $r can't be evaluated at one of the points 152 return 1 unless $rvalues; 153 154 # 155 # Handle parameters 156 # 157 $lvalues = $l->{test_values} 158 if $l->AdaptParameters($r,$self->{context}->variables->parameters); 159 160 # 161 # Look through the two lists to see if they are equal. 162 # If not, return the comparison of the first unequal value 163 # (not good for < and >, but OK for ==). 164 # 165 my ($i, $cmp); 166 foreach $i (0..scalar(@{$lvalues})-1) { 167 $cmp = $lvalues->[$i] <=> $rvalues->[$i]; 168 return $cmp if $cmp; 169 } 170 return 0; 171 } 172 173 # 174 # Create the value list from a given set of test points 175 # 176 sub createPointValues { 177 my $self = shift; 178 my $points = shift || $self->{test_points} || $self->createRandomPoints; 179 my $showError = shift; 180 my @vars = $self->{context}->variables->variables; 181 my @params = $self->{context}->variables->parameters; 182 my @zeros = @{$self->{parameters} || [split('',"0" x scalar(@params))]}; 183 my $f = $self->{f}; $f = $self->{f} = $self->perlFunction(undef,[@vars,@params]) unless $f; 184 185 my $values = []; my $v; 186 foreach my $p (@{$points}) { 187 $v = eval {&$f(@{$p},@zeros)}; 188 if (!defined($v)) { 189 return unless $showError; 190 Value::Error("Can't evaluate formula on test point (".join(',',@{$p}).")"); 191 } 192 push @{$values}, Value::makeValue($v); 193 } 194 195 $self->{test_points} = $points; 196 $self->{test_values} = $values; 197 } 198 199 # 200 # Create a list of random points, making sure that the function 201 # is defined at the given points. Error if we can't find enough. 202 # 203 sub createRandomPoints { 204 my $self = shift; 205 my $num_points = @_[0]; 206 $num_points = int($self->getFlag('num_points',5)) unless defined($num_points); 207 $num_points = 1 if $num_points < 1; 208 209 my @vars = $self->{context}->variables->variables; 210 my @params = $self->{context}->variables->parameters; 211 my @limits = $self->getVariableLimits(@vars); 212 my @make = $self->getVariableTypes(@vars); 213 my @zeros = split('',"0" x scalar(@params)); 214 my $f = $self->{f}; $f = $self->{f} = $self->perlFunction(undef,[@vars,@params]) unless $f; 215 my $seedRandom = $self->{context}->flag('random_seed')? 'PGseedRandom' : 'seedRandom'; 216 my $getRandom = $self->{context}->flag('random_seed')? 'PGgetRandom' : 'getRandom'; 217 218 $self->$seedRandom; 219 my $points = []; my $values = []; 220 my (@P,@p,$v,$i); my $k = 0; 221 while (scalar(@{$points}) < $num_points && $k < 10) { 222 @P = (); $i = 0; 223 foreach my $limit (@limits) { 224 @p = (); foreach my $I (@{$limit}) {push @p, $self->$getRandom(@{$I})} 225 push @P, $make[$i++]->make(@p); 226 } 227 $v = eval {&$f(@P,@zeros)}; 228 if (!defined($v)) {$k++} else { 229 push @{$points}, [@P]; 230 push @{$values}, Value::makeValue($v); 231 $k = 0; # reset count when we find a point 232 } 233 } 234 235 Value::Error("Can't generate enough valid points for comparison") if $k; 236 return ($points,$values) if defined(@_[0]); 237 $self->{test_values} = $values; 238 $self->{test_points} = $points; 239 } 240 241 # 242 # Get the array of variable limits 243 # 244 sub getVariableLimits { 245 my $self = shift; 246 my $userlimits = $self->{limits}; 247 if (defined($userlimits)) { 248 $userlimits = [[[-$userlimits,$userlimits]]] unless ref($userlimits) eq 'ARRAY'; 249 $userlimits = [$userlimits] unless ref($userlimits->[0]) eq 'ARRAY'; 250 $userlimits = [$userlimits] if scalar(@_) == 1 && ref($userlimits->[0][0]) ne 'ARRAY'; 251 foreach my $I (@{$userlimits}) {$I = [$I] unless ref($I->[0]) eq 'ARRAY'}; 252 } 253 $userlimits = [] unless $userlimits; my @limits; 254 my $default; $default = $userlimits->[0][0] if defined($userlimits->[0]); 255 my $default = $default || $self->{context}{flags}{limits} || [-2,2]; 256 my $granularity = $self->getFlag('granularity',1000); 257 my $resolution = $self->getFlag('resolution'); 258 my $i = 0; 259 foreach my $x (@_) { 260 my $def = $self->{context}{variables}{$x}; 261 my $limit = $userlimits->[$i++] || $def->{limits} || []; 262 $limit = [$limit] if defined($limit->[0]) && ref($limit->[0]) ne 'ARRAY'; 263 push(@{$limit},$limit->[0] || $default) while (scalar(@{$limit}) < $def->{type}{length}); 264 pop(@{$limit}) while (scalar(@{$limit}) > $def->{type}{length}); 265 push @limits, $self->addGranularity($limit,$def,$granularity,$resolution); 266 } 267 return @limits; 268 } 269 270 # 271 # Add the granularity to the limit intervals 272 # 273 sub addGranularity { 274 my $self = shift; my $limit = shift; my $def = shift; 275 my $granularity = shift; my $resolution = shift; 276 $granularity = $def->{granularity} || $granularity; 277 $resolution = $def->{resolution} || $resolution; 278 foreach my $I (@{$limit}) { 279 my ($a,$b,$n) = @{$I}; $b = -$a unless defined $b; 280 $I = [$a,$b,($n || $resolution || abs($b-$a)/$granularity)]; 281 } 282 return $limit; 283 } 284 285 # 286 # Get the routines to make the coordinates of the points 287 # 288 sub getVariableTypes { 289 my $self = shift; 290 my @make; 291 foreach my $x (@_) { 292 my $type = $self->{context}{variables}{$x}{type}; 293 if ($type->{name} eq 'Number') { 294 push @make,($type->{length} == 1)? 'Value::Formula::number': 'Value::Complex'; 295 } else { 296 push @make, "Value::$type->{name}"; 297 } 298 } 299 return @make; 300 } 301 302 # 303 # Fake object for making reals (rather than use overhead of Value::Real) 304 # 305 sub Value::Formula::number::make {shift; shift} 306 307 # 308 # Find adaptive parameters, if any 309 # 310 sub AdaptParameters { 311 my $l = shift; my $r = shift; 312 my @params = @_; my $d = scalar(@params); 313 return 0 if $d == 0; return 0 unless $l->usesOneOf(@params); 314 $l->Error("Adaptive parameters can only be used for real-valued functions") 315 unless $l->{tree}->isRealNumber; 316 # 317 # Get coefficient matrix of adaptive parameters 318 # and value vector for linear system 319 # 320 my ($p,$v) = $l->createRandomPoints($d,1); 321 my @P = split('',"0" x $d); my ($f,$F) = ($l->{f},$r->{f}); 322 my @A = (); my @b = (); 323 foreach my $i (0..$d-1) { 324 my @a = (); my @p = @{$p->[$i]}; 325 foreach my $j (0..$d-1) { 326 $P[$j] = 1; push(@a,&$f(@p,@P)-$v->[$i]); 327 $P[$j] = 0; 328 } 329 push @A, [@a]; push @b, [&$F(@p,@P)-$v->[$i]]; 330 } 331 # 332 # Use MatrixReal1.pm to solve system of linear equations 333 # 334 my $M = MatrixReal1->new($d,$d); $M->[0] = \@A; 335 my $B = MatrixReal1->new($d,1); $B->[0] = \@b; 336 ($M,$B) = $M->normalize($B); 337 $M = $M->decompose_LR; 338 if (($d,$B,$M) = $M->solve_LR($B)) { 339 if ($d == 0) { 340 # 341 # Get parameter values and recompute the points using them 342 # 343 my @a; my $i = 0; my $max = $l->getFlag('max_adapt',1E8); 344 foreach my $row (@{$B->[0]}) { 345 if (abs($row->[0]) > $max) { 346 $l->Error("Constant of integration is too large: $row->[0]") 347 if ($params[$i] eq 'C0'); 348 $l->Error("Adaptive constant is too large: $params[$i] = $row->[0]"); 349 } 350 push @a, $row->[0]; $i++; 351 } 352 $l->{parameters} = [@a]; 353 $l->createPointValues; 354 return 1; 355 } 356 } 357 $l->Error("Can't solve for adaptive parameters"); 358 } 359 360 sub usesOneOf { 361 my $self = shift; 362 foreach my $x (@_) {return 1 if $self->{variables}{$x}} 363 return 0; 364 } 365 366 ## 367 ## debugging routine 368 ## 369 #sub main::Format { 370 # my $v = scalar(@_) > 1? [@_]: shift; 371 # $v = [%{$v}] if ref($v) eq 'HASH'; 372 # return $v unless ref($v) eq 'ARRAY'; 373 # my @V; foreach my $x (@{$v}) {push @V, main::Format($x)} 374 # return '['.join(",",@V).']'; 375 #} 376 377 # 378 # Random number generator (replaced by Value::WeBWorK.pm) 379 # 380 sub seedRandom {srand} 381 sub getRandom { 382 my $self = shift; 383 my ($m,$M,$n) = @_; $n = 1 unless $n; 384 return $m + $n*int(rand()*(int(($M-$m)/$n)+1)); 385 } 386 387 # 388 # Get the value of a flag from the object itself, 389 # or from the context, or from the default context 390 # or from the given default, whichever is found first. 391 # 392 sub getFlag { 393 my $self = shift; my $name = shift; 394 return $self->{$name} if defined($self->{$name}); 395 return $self->{context}{flags}{$name} if defined($self->{context}{flags}{$name}); 396 return $$Value::context->{flags}{$name} if defined($$Value::context->{flags}{$name}); 397 return shift; 398 } 399 400 ############################################ 401 # 402 # Check if the value of a formula is constant 403 # (could use shift->{tree}{isConstant}, but I don't trust it) 404 # 405 sub isConstant {scalar(%{shift->{variables}}) == 0} 406 407 ########################################################################### 408 409 1;
| aubreyja at gmail dot com | ViewVC Help |
| Powered by ViewVC 1.0.9 |