Does math.net numerics solve non-linear systems?

[slbdeveloper]

Does math.net numerics have any methods to solve a non-linear system of equations F(x) = 0, namely Ax=b where A=A(x) and/or b=b(x)? The system of course must be solved iteratively and I see there are iterative solvers in math.net. However, do these only work for linear problems. I want to do something similar to the following but for a non-linear system: https://mathnetnumerics.codeplex.com/discussions/474432 Matlab has an fsolve function that does this: http://www.mathworks.com/help/optim/ug/fsolve.html Is it possible to solve example 1 in the above link using math.net numeric? If so can someone show me how.

See also section 6 in following link: http://web.cecs.pdx.edu/~gerry/class/ME448/notes_2012/pdf/stoppingCriteria.pdf

Thanks very much. Looking forward to your feedback.

---

The following is my second attempt of solution which is far more elegant since I use the residual as the input function to the Broyden solver to find the solution of the nonlinear equation. Broyden solves the problem in 14 iterations when using the same default accuracy as fsolve in Matlab. Do you know why the L2-norm so high? Seems like it should be easy enough for you to wrap the residual computation followed by the Broyden solve into a

function which computes the solution of a non-linear matrix equation given an initial guess from the user. Let me no what you think? Look forward to hearing from you.

    /// <summary>
    /// Run example
    /// </summary>
    public void Run()
    {
        // Format matrix output to console
        var formatProvider = (CultureInfo)CultureInfo.InvariantCulture.Clone();
        formatProvider.TextInfo.ListSeparator = " ";

        // Solve next system of non-linear equations (Ax=b):
        // 2*x1 -   x2  = exp(-x1)
        //  -x1 + 2*x2  = exp(-x2)

        // Initial guess
        var vectorX0 = Vector.Build.Dense((new[]
        {
            -5.0, -5.0
        }));
        var initialGuess = vectorX0.ToArray();
        Console.WriteLine(@"Initial vector 'x0'");
        Console.WriteLine(vectorX0.ToVectorString("#0.0000\t", formatProvider));
        Console.WriteLine();

        // Create matrix "A" with coefficients
        var matrixA = DenseMatrix.OfArray(new[,]
        {
            {2.0, -1.0},
            {-1.0, 2.0}
        });
        Console.WriteLine(@"Matrix 'A'");
        Console.WriteLine(matrixA.ToMatrixString("#0.0000\t", formatProvider));
        Console.WriteLine();

        // find root using Broyden
        var xSolution = Broyden.FindRoot(Function, initialGuess);
        var vectorX = new DenseVector(xSolution);
        Console.WriteLine(@"Result vector 'x'");
        Console.WriteLine(vectorX.ToVectorString("#0.0000\t", formatProvider));
    }

    public Double[] Function(Double[] xValues)
    {
        var matrixA = DenseMatrix.OfArray(new[,]
        {
            {2.0, -1.0},
            {-1.0, 2.0}
        });
        var vectorB = new DenseVector(new[]
        {
            Math.Exp(-xValues[0]), Math.Exp(-xValues[1])
        });
        var vectorX = new DenseVector(xValues);

        // Format matrix output to console
        var formatProvider = (CultureInfo)CultureInfo.InvariantCulture.Clone();
        formatProvider.TextInfo.ListSeparator = " ";

        iterationCount ++;         

        var residual = matrixA.Multiply(vectorX).Subtract(vectorB).ToArray();
        var residualVector = new DenseVector(xValues);

        Console.WriteLine(@"Iteration # = {0}", iterationCount);
        Console.WriteLine(@"L2-Norm  = {0}", residualVector.L2Norm());
        Console.WriteLine(@"Current vector 'x'");
        Console.WriteLine(vectorX.ToVectorString("#0.0000\t", formatProvider));
        Console.WriteLine(@"vector 'b'");
        Console.WriteLine(vectorB.ToVectorString("#0.0000\t", formatProvider));

        return residual;
    }
}

---

Written the following very "rough" code which solves the problem using a simplified Newton solver.

It of course does far worse than fsolve and probably anything that math.net numerics could do as it takes 27 iterations versus 10. This is exactly the sort of code I would like to avoid having to write.

using System; using System.Globalization; using MathNet.Numerics.LinearAlgebra.Double;

namespace Examples.LinearAlgebra.IterativeSolversExamples { public class MySolver : IExample { ///

    /// Gets the name of this example
    /// </summary>

    public string Name
    {
        get { return "Solve non-linear Ax=b"; }
    }

    /// <summary>
    /// Gets the description of this example
    /// </summary>
    public string Description
    {
        get { return "Solve non-linear using a Simplified Newton solver"; }
    }

    /// <summary>
    /// Run example
    /// </summary>
    public void Run()
    {
        // Format matrix output to console
        var formatProvider = (CultureInfo) CultureInfo.InvariantCulture.Clone();
        formatProvider.TextInfo.ListSeparator = " ";

        // Solve next system of non-linear equations (Ax=b):
        // 2*x1 -   x2  = exp(-x1)
        //  -x1 + 2*x2  = exp(-x2)

        // Initialize
        var vectorX = Vector.Build.Dense((new[]
        {
            -5.0, -5.0
        }));
        Console.WriteLine(@"Initial vector 'x0'");
        Console.WriteLine(vectorX.ToVectorString("#0.0000\t", formatProvider));
        Console.WriteLine();

        // Linearize A (already linear)
        // Create matrix "A" with coefficients
        var matrixA = DenseMatrix.OfArray(new[,]
        {
            {2.0, -1.0},
            {-1.0, 2.0}
        });
        Console.WriteLine(@"Matrix 'A'");
        Console.WriteLine(matrixA.ToMatrixString("#0.0000\t", formatProvider));
        Console.WriteLine();

        // Loop over k
        for (var k = 0; k < 50; k++)
        {
            // Linearize b
            // Create vector "b" with the constant terms.
            var vectorB = new DenseVector(new[]
            {
                Math.Exp(-vectorX.At(0)), Math.Exp(-vectorX.At(1))
            });

            // Compute residual
            var res = matrixA.Multiply(vectorX).Subtract(vectorB);
            var norm = res.L2Norm();
            if (norm < 1e-6)
            {
                break; // converged
            }

            var delta = -res;
            vectorX = vectorX.Add(delta);

            Console.WriteLine(@"Iteration # = {0}", k);
            Console.WriteLine(@"Residual = {0}", res);
            Console.WriteLine(@"L2-Norm  = {0}", norm);
            Console.WriteLine(@"Current vector 'x'");
            Console.WriteLine(vectorX.ToVectorString("#0.0000\t", formatProvider));
            Console.WriteLine();
        }

        Console.WriteLine(@"Result vector 'x'");
        Console.WriteLine(vectorX.ToVectorString("#0.0000\t", formatProvider));
        Console.WriteLine();
    }
}

}

Create a new C# console app with the following as the main function:

using System; using Examples.LinearAlgebra.IterativeSolversExamples;

namespace ExamplesRunner { public class Program { public static void Main(string[] args) { var mySolver = new MySolver(); mySolver.Run(); Console.ReadKey(); } }

}

[NilesDavis]

Any progress here? I found some pretty promising but basic level approach here: http://www.imagingshop.com/linear-and-nonlinear-least-squares-with-math-net/