# Nonlinear Systems in IronPython QuickStart Sample

Illustrates the use of the NewtonRaphsonSystemSolver and DoglegSystemSolver classes for solving systems of nonlinear equations in IronPython.

View this sample in: C# Visual Basic F#

```Python import numerics from math import exp, sin, cos from System import Array, Func # The equation solver classes reside in the # Extreme.Mathematics.EquationSolvers namespace. from Extreme.Mathematics.EquationSolvers import * # Function delegates reside in the Extreme.Mathematics # namespace. from Extreme.Mathematics import * # Vectors reside in the Extreme.Mathematics.LinearAlgebra # namespace. from Extreme.Mathematics.LinearAlgebra import * # Illustrates solving systems of non-linear equations using # classes in the Extreme.Mathematics.EquationSolvers namespace # of Extreme Numerics.NET. # # Target function # # The function we are trying to solve can be provided # on one of two ways. The first is as an array of # Func<Vector, double> delegates. See the end of this # sample for definitions of the methods that are referenced here. f = Array[Func[Vector, float]]([ \ lambda x: exp(x[0])*cos(x[1]) - x[0]*x[0] + x[1]*x[1], \ lambda x: exp(x[0])*sin(x[1]) - 2*x[0]*x[1] ]) # We can also supply the Jacobian, which is the matrix of partial # derivatives. We do so by providing the gradient of each target # function as a delegate that takes a second argument: # the vector that is to hold the return value. This avoids unnecessary # creation of new Vector instances. def df1(x,y): y[0] = exp(x[0])*cos(x[1]) - 2*x[0] y[1] = -exp(x[0])*sin(x[1]) + 2*x[1] return y def df2(x,y): y[0] = exp(x[0])*sin(x[1]) - 2*x[1] y[1] = exp(x[0])*cos(x[1]) - 2*x[0] return y df = Array[Func[Vector,Vector,Vector]]([ df1, df2 ]) # The initial values are supplied as a vector: initialGuess = Vector([0.5, 0.5]) # # Newton-Raphson Method # # The Newton-Raphson method is implemented by # the NewtonRaphsonSystemSolver class. solver = NewtonRaphsonSystemSolver(f, df, initialGuess) # and call the Solve method to obtain the solution: solution = solver.Solve() print "N-dimensional Newton-Raphson Solver:" print "exp(x)*cos(y) - x^2 + y^2 = 0" print "exp(x)*sin(y) - 2xy = 0" print " Initial guess: {0:F2}", initialGuess # The Status property indicates # the result of running the algorithm. print " Status:", solver.Status # The result is also available through the # Result property. print " Solution:", solver.Result print " Function value:", solver.ValueTest.Error # You can find out the estimated error of the result # through the EstimatedError property: print " Estimated error:", solver.EstimatedError print " # iterations:", solver.IterationsNeeded print " # evaluations:", solver.EvaluationsNeeded ### Symbolic derivatives are currently not supported ### ## In .NET 4.0, you can use Automatic Differentiation to compute the Jacobian. ## To do so, set the target functions using the SetSymbolicTargetFunctions ## method: #solver = NewtonRaphsonSystemSolver() #solver.InitialGuess = initialGuess # #solver.SetSymbolicTargetFunctions( \ # lambda x: exp(x[0]) * cos(x[1]) - x[0] * x[0] + x[1] * x[1], \ # lambda x: exp(x[0]) * sin(x[1]) - 2 * x[0] * x[1]) # #solution = solver.Solve() #print "Using Automatic Differentiation:" #print " Solution:", solver.Result #print " Function value:", solver.ValueTest.Error #print " # iterations:", solver.IterationsNeeded # When you don't have the Jacobian for the target functions # and you don't use Automatic Differentiation, the equation solver # will use a numerical approximation. # # Controlling the process # print "Same with modified parameters:" # You can set the maximum # of iterations: # If the solution cannot be found in time, the # Status will return a value of # IterationStatus.IterationLimitExceeded solver.MaxIterations = 10 # The ValueTest property returns the convergence # test based on the function value. We can set # its tolerance property: solver.ValueTest.Tolerance = 1e-10 # Its Norm property determines how the error # is calculated. Here, we choose the maximum # of the function values: solver.ValueTest.Norm = Algorithms.VectorConvergenceNorm.Maximum # The SolutionTest property returns the test # on the change in location of the solution. solver.SolutionTest.Tolerance = 1e-8 # You can specify how convergence is to be tested # through the ConvergenceCriterion property: solver.SolutionTest.ConvergenceCriterion = ConvergenceCriterion.WithinRelativeTolerance solver.InitialGuess = initialGuess solution = solver.Solve() print " Status:", solver.Status print " Solution:", solver.Result # The estimated error will be less than 5e-14 print " Estimated error:", solver.SolutionTest.Error print " # iterations:", solver.IterationsNeeded print " # evaluations:", solver.EvaluationsNeeded # # Powell's dogleg method # # The dogleg method is more robust than Newton's method. # It will converge often when Newton's method fails. dogleg = DoglegSystemSolver(f, df, initialGuess) # Unique to the dogleg method is the TrustRegionRadius property. # Any step of the algorithm is not larger than this value. # It is adjusted at each iteration. dogleg.TrustRegionRadius = 0.5 # Call the Solve method to obtain the solution: solution = dogleg.Solve() print "Powell's Dogleg Solver:" print " Initial guess: {0:F2}", initialGuess print " Status:", dogleg.Status print " Solution:", dogleg.Result print " Estimated error:", dogleg.EstimatedError print " # iterations:", dogleg.IterationsNeeded print " # evaluations:", dogleg.EvaluationsNeeded # The dogleg method can work without derivatives. Care is taken # to keep the number of evaluations down to a minimum. dogleg.JacobianFunction = None # Call the Solve method to obtain the solution: solution = dogleg.Solve() print "Powell's Dogleg Solver (no derivatives):" print " Initial guess: {0:F2}", initialGuess print " Status:", dogleg.Status print " Solution:", dogleg.Result print " Estimated error:", dogleg.EstimatedError print " # iterations:", dogleg.IterationsNeeded print " # evaluations:", dogleg.EvaluationsNeeded ```