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)

```