// This contains routines for integrating a function numerically using // interval techniques. // This function takes an anonymous one-argument function f and numerically // integrates it over n steps from lower to upper. The result is an interval // that should be guaranted to contain the true value of the function. IntervalIntegrate[f, lower, upper, steps] := { // Make sure steps is an integer or this will fail. steps = ceil[steps] stepsize = (upper-lower) / steps sum = 0 f[lower] stepsize // This is necessary to get the units right. min = lower max = lower // By counting for a guaranteed number of steps, we can avoid roundoff // error and "missing" and endpoint when using floating-point boundaries. while min < upper { min = max max = min + stepsize if (max > upper) max = upper width = max-min mid = (min + max)/2 int = new interval[min, mid, max] /* println["int: $int"] println["f[int]: " + f[int]] println["width: $width"] println["sum: $sum"]*/ area = f[int] * width // println["$int\t$area"] sum = sum + area } return sum } // This equation tests the integration with an anonymous function and returns // an interval. //IntervalIntegrate[{|x| sin[x]/x}, 1e-10, 2 pi, 100000] // This represents the light from a long cylindrical lightbulb //IntervalIntegrate[{|h| 1000. lumens/ft / (h^2 + (12 in)^2)}, -1/2 ft, 1/2 ft, 10000] -> lumens/(foot^2)