Basic Integration in IronPython QuickStart Sample
Illustrates the basic numerical integration classes in IronPython.
View this sample in: C# Visual Basic F#
```Python import numerics from math import * # The numerical integration classes reside in the # Extreme.Mathematics.Calculus namespace. from Extreme.Mathematics import * from Extreme.Mathematics.Calculus import * from Extreme.Mathematics.Algorithms import * # Illustrates the basic use of the numerical integration # classes in the Extreme.Mathematics.Calculus namespace of the Extreme # Optimization Mathematics Library for .NET. # Numerical integration algorithms fall into two # main categories: adaptive and non-adaptive. # This QuickStart Sample illustrates the use of # the non-adaptive numerical integrators. # # All numerical integration classes derive from # NumericalIntegrator. This abstract base class # defines properties and methods that are shared # by all numerical integration classes. # # The integrand # # The function we are integrating must be # provided as a Func<double, double>. f = sin # # SimpsonIntegrator # # The simplest numerical integration algorithm # is Simpson's rule. simpson = SimpsonIntegrator() # You can set the relative or absolute tolerance # to which to evaluate the integral. simpson.RelativeTolerance = 1e-5 # You can select the type of tolerance using the # ConvergenceCriterion property: simpson.ConvergenceCriterion = ConvergenceCriterion.WithinRelativeTolerance # The Integrate method performs the actual # integration: result = simpson.Integrate(sin, 0, 2) print "sin(x) on [0,2]" print "Simpson integrator:" # The result is also available in the Result # property: print " Value:", simpson.Result # To see whether the algorithm ended normally, # inspect the Status property: print " Status:", simpson.Status # You can find out the estimated error of the result # through the EstimatedError property: print " Estimated error:", simpson.EstimatedError # The number of iterations to achieve the result # is available through the IterationsNeeded property. print " Iterations:", simpson.IterationsNeeded # The number of function evaluations is available # through the EvaluationsNeeded property. print " Function evaluations:", simpson.EvaluationsNeeded # # Gauss-Kronrod Integration # # Gauss-Kronrod integrators also use a fixed point # scheme, but with certain optimizations in the # choice of points where the integrand is evaluated. # The NonAdaptiveGaussKronrodIntegrator uses a # succession of 10, 21, 43, and 87 point rules # to approximate the integral. nagk = NonAdaptiveGaussKronrodIntegrator() nagk.Integrate(sin, 0, 2) print "Non-adaptive Gauss-Kronrod rule:" print " Value:", nagk.Result print " Status:", nagk.Status print " Estimated error:", nagk.EstimatedError print " Iterations:", nagk.IterationsNeeded print " Function evaluations:", nagk.EvaluationsNeeded # # Romberg Integration # # Romberg integration combines Simpson's Rule # with a scheme to accelerate convergence. # This algorithm is useful for smooth integrands. romberg = RombergIntegrator() result = romberg.Integrate(sin, 0, 2) print "Romberg integration:" print " Value:", romberg.Result print " Status:", romberg.Status print " Estimated error:", romberg.EstimatedError print " Iterations:", romberg.IterationsNeeded print " Function evaluations:", romberg.EvaluationsNeeded # However, it breaks down if the integration # algorithm contains singularities or # discontinuities. # # The AdaptiveIntegrator can handle this type # of integrand, and many other difficult cases. # See the AdvancedIntegration QuickStart sample # for details. result = romberg.Integrate(lambda x: 0.0 if x <= 0.0 else x ** -0.9 * log(1/x), 0.0, 1.0) print "Romberg on hard integrand:" print " Value:", romberg.Result print " Actual value: 100" print " Status:", romberg.Status print " Estimated error:", romberg.EstimatedError print " Iterations:", romberg.IterationsNeeded print " Function evaluations:", romberg.EvaluationsNeeded # Function that will cause difficulties to the # simplistic integration algorithms. # </summary> def HardIntegrand(x): # This is put in because some integration rules # evaluate the function at x=0. if x <= 0: return 0 return x**-0.9 * log(1/x) ```