LA: member methods instead of extension methods for matrix factorizations

src/Numerics/Distributions/InverseWishart.cs

@@ -129,7 +129,7 @@ namespace MathNet.Numerics.Distributions
 
             _freedom = degreesOfFreedom;
             _scale = scale;
-            _chol = Cholesky<double>.Create(_scale);
+            _chol = _scale.Cholesky();
         }
 
         /// <summary>
@@ -217,7 +217,7 @@ namespace MathNet.Numerics.Distributions
                 throw new ArgumentOutOfRangeException("x", Resources.ArgumentMatrixDimensions);
             }
 
-            var chol = Cholesky<double>.Create(x);
+            var chol = x.Cholesky();
             var dX = chol.Determinant;
             var sXi = chol.Solve(Scale);
 

src/Numerics/Distributions/MatrixNormal.cs

@@ -216,8 +216,8 @@ namespace MathNet.Numerics.Distributions
             }
 
             var a = x - _m;
-            var cholV = Cholesky<double>.Create(_v);
+            var cholV = _v.Cholesky();
-            var cholK = Cholesky<double>.Create(_k);
+            var cholK = _k.Cholesky();
 
             return Math.Exp(-0.5*cholV.Solve(a.Transpose()*cholK.Solve(a)).Trace())
                    /Math.Pow(2.0*Constants.Pi, x.RowCount*x.ColumnCount/2.0)
@@ -281,7 +281,7 @@ namespace MathNet.Numerics.Distributions
         /// <returns>a sequence of samples from defined distribution.</returns>
         static Vector<double> SampleVectorNormal(System.Random rnd, Vector<double> mean, Matrix<double> covariance)
         {
-            var chol = Cholesky<double>.Create(covariance);
+            var chol = covariance.Cholesky();
 
             // Sample a standard normal variable.
             var v = DenseVector.CreateRandom(mean.Count, new Normal(rnd));

src/Numerics/Distributions/Wishart.cs

@@ -105,7 +105,7 @@ namespace MathNet.Numerics.Distributions
 
             _degreesOfFreedom = degreesOfFreedom;
             _scale = scale;
-            _chol = Cholesky<double>.Create(_scale);
+            _chol = _scale.Cholesky();
         }
 
         /// <summary>
@@ -281,7 +281,7 @@ namespace MathNet.Numerics.Distributions
                 throw new ArgumentOutOfRangeException(Resources.InvalidDistributionParameters);
             }
 
-            return DoSample(rnd, degreesOfFreedom, scale, Cholesky<double>.Create(scale));
+            return DoSample(rnd, degreesOfFreedom, scale, scale.Cholesky());
         }
 
         /// <summary>

src/Numerics/LinearAlgebra/Complex/DenseMatrix.cs

@@ -33,6 +33,8 @@ using System.Collections.Generic;
 using System.Diagnostics;
 using System.Linq;
 using MathNet.Numerics.Distributions;
+using MathNet.Numerics.LinearAlgebra.Complex.Factorization;
+using MathNet.Numerics.LinearAlgebra.Factorization;
 using MathNet.Numerics.LinearAlgebra.Storage;
 using MathNet.Numerics.Properties;
 using MathNet.Numerics.Providers.LinearAlgebra;
@@ -1006,5 +1008,35 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
 
             return (DenseMatrix)leftSide.Modulus(rightSide);
         }
+
+        public override Cholesky<Complex> Cholesky()
+        {
+            return new DenseCholesky(this);
+        }
+
+        public override LU<Complex> LU()
+        {
+            return new DenseLU(this);
+        }
+
+        public override QR<Complex> QR(QRMethod method = QRMethod.Thin)
+        {
+            return new DenseQR(this, method);
+        }
+
+        public override GramSchmidt<Complex> GramSchmidt()
+        {
+            return new DenseGramSchmidt(this);
+        }
+
+        public override Svd<Complex> Svd(bool computeVectors)
+        {
+            return new DenseSvd(this, computeVectors);
+        }
+
+        public override Evd<Complex> Evd()
+        {
+            return new DenseEvd(this);
+        }
     }
 }

src/Numerics/LinearAlgebra/Complex/ExtensionMethods.cs

@@ -1,106 +0,0 @@
-// <copyright file="ExtensionMethods.cs" company="Math.NET">
-// Math.NET Numerics, part of the Math.NET Project
-// http://numerics.mathdotnet.com
-// http://github.com/mathnet/mathnet-numerics
-// http://mathnetnumerics.codeplex.com
-// Copyright (c) 2009-2010 Math.NET
-// Permission is hereby granted, free of charge, to any person
-// obtaining a copy of this software and associated documentation
-// files (the "Software"), to deal in the Software without
-// restriction, including without limitation the rights to use,
-// copy, modify, merge, publish, distribute, sublicense, and/or sell
-// copies of the Software, and to permit persons to whom the
-// Software is furnished to do so, subject to the following
-// conditions:
-// The above copyright notice and this permission notice shall be
-// included in all copies or substantial portions of the Software.
-// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
-// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
-// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
-// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
-// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
-// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
-// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
-// OTHER DEALINGS IN THE SOFTWARE.
-// </copyright>
-
-using MathNet.Numerics.LinearAlgebra.Complex.Factorization;
-using MathNet.Numerics.LinearAlgebra.Factorization;
-
-namespace MathNet.Numerics.LinearAlgebra.Complex
-{
-
-#if NOSYSNUMERICS
-    using Numerics;
-#else
-    using System.Numerics;
-#endif
-
-    /// <summary>
-    /// Extension methods which return factorizations for the various matrix classes.
-    /// </summary>
-    public static class ExtensionMethods
-    {
-        /// <summary>
-        /// Computes the Cholesky decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The Cholesky decomposition object.</returns>
-        public static Cholesky Cholesky(this Matrix<Complex> matrix)
-        {
-            return (Cholesky)Cholesky<Complex>.Create(matrix);
-        }
-
-        /// <summary>
-        /// Computes the LU decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The LU decomposition object.</returns>
-        public static LU LU(this Matrix<Complex> matrix)
-        {
-            return (LU)LU<Complex>.Create(matrix);
-        }
-
-        /// <summary>
-        /// Computes the QR decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <param name="method">The type of QR factorization to perform.</param>
-        /// <returns>The QR decomposition object.</returns>
-        public static QR QR(this Matrix<Complex> matrix, QRMethod method = QRMethod.Full)
-        {
-            return (QR)QR<Complex>.Create(matrix, method);
-        }
-
-        /// <summary>
-        /// Computes the QR decomposition for a matrix using Modified Gram-Schmidt Orthogonalization.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The QR decomposition object.</returns>
-        public static GramSchmidt GramSchmidt(this Matrix<Complex> matrix)
-        {
-            return (GramSchmidt)GramSchmidt<Complex>.Create(matrix);
-        }
-
-        /// <summary>
-        /// Computes the SVD decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <param name="computeVectors">Compute the singular U and VT vectors or not.</param>
-        /// <returns>The SVD decomposition object.</returns>
-        public static Svd Svd(this Matrix<Complex> matrix, bool computeVectors)
-        {
-            return (Svd)Svd<Complex>.Create(matrix, computeVectors);
-        }
-
-        /// <summary>
-        /// Computes the EVD decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The EVD decomposition object.</returns>
-        public static Evd Evd(this Matrix<Complex> matrix)
-        {
-            return (Evd)Evd<Complex>.Create(matrix);
-        }
-    }
-}

src/Numerics/LinearAlgebra/Complex/Matrix.cs

@@ -29,6 +29,8 @@
 // </copyright>
 
 using System;
+using MathNet.Numerics.LinearAlgebra.Complex.Factorization;
+using MathNet.Numerics.LinearAlgebra.Factorization;
 using MathNet.Numerics.LinearAlgebra.Storage;
 using MathNet.Numerics.Properties;
 
@@ -455,5 +457,35 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
 
             return sum;
         }
+
+        public override Cholesky<Complex> Cholesky()
+        {
+            return new UserCholesky(this);
+        }
+
+        public override LU<Complex> LU()
+        {
+            return new UserLU(this);
+        }
+
+        public override QR<Complex> QR(QRMethod method = QRMethod.Thin)
+        {
+            return new UserQR(this, method);
+        }
+
+        public override GramSchmidt<Complex> GramSchmidt()
+        {
+            return new UserGramSchmidt(this);
+        }
+
+        public override Svd<Complex> Svd(bool computeVectors)
+        {
+            return new UserSvd(this, computeVectors);
+        }
+
+        public override Evd<Complex> Evd()
+        {
+            return new UserEvd(this);
+        }
     }
 }

src/Numerics/LinearAlgebra/Complex/Solvers/Iterator.cs

@@ -46,6 +46,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Solvers
     /// </summary>
     public static class Iterator
     {
+
+        // TODO: Refactor
+
         /// <summary>
         /// Creates a default iterator with all the <see cref="IIterationStopCriterium{T}"/> objects.
         /// </summary>

src/Numerics/LinearAlgebra/Complex32/DenseMatrix.cs

@@ -33,6 +33,8 @@ using System.Collections.Generic;
 using System.Diagnostics;
 using System.Linq;
 using MathNet.Numerics.Distributions;
+using MathNet.Numerics.LinearAlgebra.Complex32.Factorization;
+using MathNet.Numerics.LinearAlgebra.Factorization;
 using MathNet.Numerics.LinearAlgebra.Storage;
 using MathNet.Numerics.Properties;
 using MathNet.Numerics.Providers.LinearAlgebra;
@@ -1001,5 +1003,35 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
 
             return (DenseMatrix)leftSide.Modulus(rightSide);
         }
+
+        public override Cholesky<Complex32> Cholesky()
+        {
+            return new DenseCholesky(this);
+        }
+
+        public override LU<Complex32> LU()
+        {
+            return new DenseLU(this);
+        }
+
+        public override QR<Complex32> QR(QRMethod method = QRMethod.Thin)
+        {
+            return new DenseQR(this, method);
+        }
+
+        public override GramSchmidt<Complex32> GramSchmidt()
+        {
+            return new DenseGramSchmidt(this);
+        }
+
+        public override Svd<Complex32> Svd(bool computeVectors)
+        {
+            return new DenseSvd(this, computeVectors);
+        }
+
+        public override Evd<Complex32> Evd()
+        {
+            return new DenseEvd(this);
+        }
     }
 }

src/Numerics/LinearAlgebra/Complex32/ExtensionMethods.cs

@@ -1,101 +0,0 @@
-// <copyright file="ExtensionMethods.cs" company="Math.NET">
-// Math.NET Numerics, part of the Math.NET Project
-// http://numerics.mathdotnet.com
-// http://github.com/mathnet/mathnet-numerics
-// http://mathnetnumerics.codeplex.com
-// Copyright (c) 2009-2010 Math.NET
-// Permission is hereby granted, free of charge, to any person
-// obtaining a copy of this software and associated documentation
-// files (the "Software"), to deal in the Software without
-// restriction, including without limitation the rights to use,
-// copy, modify, merge, publish, distribute, sublicense, and/or sell
-// copies of the Software, and to permit persons to whom the
-// Software is furnished to do so, subject to the following
-// conditions:
-// The above copyright notice and this permission notice shall be
-// included in all copies or substantial portions of the Software.
-// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
-// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
-// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
-// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
-// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
-// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
-// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
-// OTHER DEALINGS IN THE SOFTWARE.
-// </copyright>
-
-using MathNet.Numerics.LinearAlgebra.Complex32.Factorization;
-using MathNet.Numerics.LinearAlgebra.Factorization;
-
-namespace MathNet.Numerics.LinearAlgebra.Complex32
-{
-    using Numerics;
-
-    /// <summary>
-    /// Extension methods which return factorizations for the various matrix classes.
-    /// </summary>
-    public static class ExtensionMethods
-    {
-        /// <summary>
-        /// Computes the Cholesky decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The Cholesky decomposition object.</returns>
-        public static Cholesky Cholesky(this Matrix<Complex32> matrix)
-        {
-            return (Cholesky)Cholesky<Complex32>.Create(matrix);
-        }
-
-        /// <summary>
-        /// Computes the LU decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The LU decomposition object.</returns>
-        public static LU LU(this Matrix<Complex32> matrix)
-        {
-            return (LU)LU<Complex32>.Create(matrix);
-        }
-
-        /// <summary>
-        /// Computes the QR decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <param name="method">The type of QR factorization to perform.</param>
-        /// <returns>The QR decomposition object.</returns>
-        public static QR QR(this Matrix<Complex32> matrix, QRMethod method = QRMethod.Full)
-        {
-            return (QR)QR<Complex32>.Create(matrix, method);
-        }
-
-        /// <summary>
-        /// Computes the QR decomposition for a matrix using Modified Gram-Schmidt Orthogonalization.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The QR decomposition object.</returns>
-        public static GramSchmidt GramSchmidt(this Matrix<Complex32> matrix)
-        {
-            return (GramSchmidt)GramSchmidt<Complex32>.Create(matrix);
-        }
-
-        /// <summary>
-        /// Computes the SVD decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <param name="computeVectors">Compute the singular U and VT vectors or not.</param>
-        /// <returns>The SVD decomposition object.</returns>
-        public static Svd Svd(this Matrix<Complex32> matrix, bool computeVectors)
-        {
-            return (Svd)Svd<Complex32>.Create(matrix, computeVectors);
-        }
-
-        /// <summary>
-        /// Computes the EVD decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The EVD decomposition object.</returns>
-        public static Evd Evd(this Matrix<Complex32> matrix)
-        {
-            return (Evd)Evd<Complex32>.Create(matrix);
-        }
-    }
-}

src/Numerics/LinearAlgebra/Complex32/Matrix.cs

@@ -28,6 +28,8 @@
 // OTHER DEALINGS IN THE SOFTWARE.
 // </copyright>
 
+using MathNet.Numerics.LinearAlgebra.Complex32.Factorization;
+using MathNet.Numerics.LinearAlgebra.Factorization;
 using MathNet.Numerics.LinearAlgebra.Storage;
 using MathNet.Numerics.Properties;
 using System;
@@ -450,5 +452,35 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
 
             return sum;
         }
+
+        public override Cholesky<Complex32> Cholesky()
+        {
+            return new UserCholesky(this);
+        }
+
+        public override LU<Complex32> LU()
+        {
+            return new UserLU(this);
+        }
+
+        public override QR<Complex32> QR(QRMethod method = QRMethod.Thin)
+        {
+            return new UserQR(this, method);
+        }
+
+        public override GramSchmidt<Complex32> GramSchmidt()
+        {
+            return new UserGramSchmidt(this);
+        }
+
+        public override Svd<Complex32> Svd(bool computeVectors)
+        {
+            return new UserSvd(this, computeVectors);
+        }
+
+        public override Evd<Complex32> Evd()
+        {
+            return new UserEvd(this);
+        }
     }
 }

src/Numerics/LinearAlgebra/Complex32/Solvers/Iterator.cs

@@ -41,6 +41,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Solvers
     /// </summary>
     public static class Iterator
     {
+
+        // TODO: Refactor
+
         /// <summary>
         /// Creates a default iterator with all the <see cref="IIterationStopCriterium{T}"/> objects.
         /// </summary>

src/Numerics/LinearAlgebra/Double/DenseMatrix.cs

@@ -33,6 +33,8 @@ using System.Collections.Generic;
 using System.Diagnostics;
 using System.Linq;
 using MathNet.Numerics.Distributions;
+using MathNet.Numerics.LinearAlgebra.Double.Factorization;
+using MathNet.Numerics.LinearAlgebra.Factorization;
 using MathNet.Numerics.LinearAlgebra.Storage;
 using MathNet.Numerics.Properties;
 using MathNet.Numerics.Providers.LinearAlgebra;
@@ -1032,5 +1034,35 @@ namespace MathNet.Numerics.LinearAlgebra.Double
 
             return (DenseMatrix)leftSide.Modulus(rightSide);
         }
+
+        public override Cholesky<double> Cholesky()
+        {
+            return new DenseCholesky(this);
+        }
+
+        public override LU<double> LU()
+        {
+            return new DenseLU(this);
+        }
+
+        public override QR<double> QR(QRMethod method = QRMethod.Thin)
+        {
+            return new DenseQR(this, method);
+        }
+
+        public override GramSchmidt<double> GramSchmidt()
+        {
+            return new DenseGramSchmidt(this);
+        }
+
+        public override Svd<double> Svd(bool computeVectors)
+        {
+            return new DenseSvd(this, computeVectors);
+        }
+
+        public override Evd<double> Evd()
+        {
+            return new DenseEvd(this);
+        }
     }
 }

src/Numerics/LinearAlgebra/Double/ExtensionMethods.cs

@@ -1,99 +0,0 @@
-// <copyright file="ExtensionMethods.cs" company="Math.NET">
-// Math.NET Numerics, part of the Math.NET Project
-// http://numerics.mathdotnet.com
-// http://github.com/mathnet/mathnet-numerics
-// http://mathnetnumerics.codeplex.com
-// Copyright (c) 2009-2010 Math.NET
-// Permission is hereby granted, free of charge, to any person
-// obtaining a copy of this software and associated documentation
-// files (the "Software"), to deal in the Software without
-// restriction, including without limitation the rights to use,
-// copy, modify, merge, publish, distribute, sublicense, and/or sell
-// copies of the Software, and to permit persons to whom the
-// Software is furnished to do so, subject to the following
-// conditions:
-// The above copyright notice and this permission notice shall be
-// included in all copies or substantial portions of the Software.
-// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
-// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
-// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
-// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT
-// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY,
-// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
-// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR
-// OTHER DEALINGS IN THE SOFTWARE.
-// </copyright>
-
-using MathNet.Numerics.LinearAlgebra.Double.Factorization;
-using MathNet.Numerics.LinearAlgebra.Factorization;
-
-namespace MathNet.Numerics.LinearAlgebra.Double
-{
-    /// <summary>
-    /// Extension methods which return factorizations for the various matrix classes.
-    /// </summary>
-    public static class ExtensionMethods
-    {
-        /// <summary>
-        /// Computes the Cholesky decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The Cholesky decomposition object.</returns>
-        public static Cholesky Cholesky(this Matrix<double> matrix)
-        {
-            return (Cholesky)Cholesky<double>.Create(matrix);
-        }
-
-        /// <summary>
-        /// Computes the LU decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The LU decomposition object.</returns>
-        public static LU LU(this Matrix<double> matrix)
-        {
-            return (LU)LU<double>.Create(matrix);
-        }
-
-        /// <summary>
-        /// Computes the QR decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <param name="method">The type of QR factorization to perform.</param>
-        /// <returns>The QR decomposition object.</returns>
-        public static QR QR(this Matrix<double> matrix, QRMethod method = QRMethod.Full)
-        {
-            return (QR)QR<double>.Create(matrix, method);
-        }
-
-        /// <summary>
-        /// Computes the QR decomposition for a matrix using Modified Gram-Schmidt Orthogonalization.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The QR decomposition object.</returns>
-        public static GramSchmidt GramSchmidt(this Matrix<double> matrix)
-        {
-            return (GramSchmidt)GramSchmidt<double>.Create(matrix);
-        }
-
-        /// <summary>
-        /// Computes the SVD decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <param name="computeVectors">Compute the singular U and VT vectors or not.</param>
-        /// <returns>The SVD decomposition object.</returns>
-        public static Svd Svd(this Matrix<double> matrix, bool computeVectors)
-        {
-            return (Svd)Svd<double>.Create(matrix, computeVectors);
-        }
-
-        /// <summary>
-        /// Computes the EVD decomposition for a matrix.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>The EVD decomposition object.</returns>
-        public static Evd Evd(this Matrix<double> matrix)
-        {
-            return (Evd)Evd<double>.Create(matrix);
-        }
-    }
-}

src/Numerics/LinearAlgebra/Double/Matrix.cs

@@ -29,6 +29,8 @@
 // </copyright>
 
 using System;
+using MathNet.Numerics.LinearAlgebra.Double.Factorization;
+using MathNet.Numerics.LinearAlgebra.Factorization;
 using MathNet.Numerics.LinearAlgebra.Storage;
 using MathNet.Numerics.Properties;
 
@@ -456,5 +458,35 @@ namespace MathNet.Numerics.LinearAlgebra.Double
 
             return sum;
         }
+
+        public override Cholesky<double> Cholesky()
+        {
+            return new UserCholesky(this);
+        }
+
+        public override LU<double> LU()
+        {
+            return new UserLU(this);
+        }
+
+        public override QR<double> QR(QRMethod method = QRMethod.Thin)
+        {
+            return new UserQR(this, method);
+        }
+
+        public override GramSchmidt<double> GramSchmidt()
+        {
+            return new UserGramSchmidt(this);
+        }
+
+        public override Svd<double> Svd(bool computeVectors)
+        {
+            return new UserSvd(this, computeVectors);
+        }
+
+        public override Evd<double> Evd()
+        {
+            return new UserEvd(this);
+        }
     }
 }

src/Numerics/LinearAlgebra/Double/Solvers/Iterator.cs

@@ -39,6 +39,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Solvers
     /// </summary>
     public static class Iterator
     {
+
+        // TODO: Refactor
+
         /// <summary>
         /// Creates a default iterator with all the <see cref="IIterationStopCriterium{T}"/> objects.
         /// </summary>

src/Numerics/LinearAlgebra/Factorization/Cholesky.cs

@@ -32,12 +32,6 @@ using System;
 
 namespace MathNet.Numerics.LinearAlgebra.Factorization
 {
-    using Numerics;
-
-#if !NOSYSNUMERICS
-    using System.Numerics;
-#endif
-
     /// <summary>
     /// <para>A class which encapsulates the functionality of a Cholesky factorization.</para>
     /// <para>For a symmetric, positive definite matrix A, the Cholesky factorization
@@ -51,68 +45,10 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization
     public abstract class Cholesky<T> : ISolver<T>
     where T : struct, IEquatable<T>, IFormattable
     {
-        /// <summary>
-        /// Internal method which routes the call to perform the Cholesky factorization to the appropriate class.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>A cholesky factorization object.</returns>
-        internal static Cholesky<T> Create(Matrix<T> matrix)
-        {
-            if (typeof(T) == typeof(double))
-            {
-                var dense = matrix as Double.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Double.Factorization.DenseCholesky(dense) as Cholesky<T>;
-                }
-
-                return new Double.Factorization.UserCholesky(matrix as Matrix<double>) as Cholesky<T>;
-            }
-
-            if (typeof(T) == typeof(float))
-            {
-                var dense = matrix as Single.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Single.Factorization.DenseCholesky(dense) as Cholesky<T>;
-                }
-
-                return new Single.Factorization.UserCholesky(matrix as Matrix<float>) as Cholesky<T>;
-            }
-
-            if (typeof(T) == typeof(Complex))
-            {
-                var dense = matrix as LinearAlgebra.Complex.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex.Factorization.DenseCholesky(dense) as Cholesky<T>;
-                }
-
-                return new LinearAlgebra.Complex.Factorization.UserCholesky(matrix as Matrix<Complex>) as Cholesky<T>;
-            }
-
-            if (typeof(T) == typeof(Complex32))
-            {
-                var dense = matrix as LinearAlgebra.Complex32.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex32.Factorization.DenseCholesky(dense) as Cholesky<T>;
-                }
-
-                return new LinearAlgebra.Complex32.Factorization.UserCholesky(matrix as Matrix<Complex32>) as Cholesky<T>;
-            }
-
-            throw new NotSupportedException();
-        }
-
         /// <summary>
         /// Gets or sets the lower triangular form of the Cholesky matrix
         /// </summary>
-        protected Matrix<T> CholeskyFactor
+        protected Matrix<T> CholeskyFactor { get; set; }
-        {
-            get;
-            set;
-        }
 
         /// <summary>
         /// Gets the lower triangular form of the Cholesky matrix.
@@ -128,18 +64,12 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization
         /// <summary>
         /// Gets the determinant of the matrix for which the Cholesky matrix was computed.
         /// </summary>
-        public abstract T Determinant
+        public abstract T Determinant { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Gets the log determinant of the matrix for which the Cholesky matrix was computed.
         /// </summary>
-        public abstract T DeterminantLn
+        public abstract T DeterminantLn { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Solves a system of linear equations, <b>AX = B</b>, with A Cholesky factorized.

src/Numerics/LinearAlgebra/Factorization/Evd.cs

@@ -60,118 +60,39 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization
         /// <summary>
         /// Gets or sets a value indicating whether matrix is symmetric or not
         /// </summary>
-        public bool IsSymmetric
+        public bool IsSymmetric { get; protected set; }
-        {
-            get;
-            protected set;
-        }
 
         /// <summary>
         /// Gets the absolute value of determinant of the square matrix for which the EVD was computed.
         /// </summary>
-        public abstract T Determinant
+        public abstract T Determinant { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Gets the effective numerical matrix rank.
         /// </summary>
         /// <value>The number of non-negligible singular values.</value>
-        public abstract int Rank
+        public abstract int Rank { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Gets a value indicating whether the matrix is full rank or not.
         /// </summary>
         /// <value><c>true</c> if the matrix is full rank; otherwise <c>false</c>.</value>
-        public abstract bool IsFullRank
+        public abstract bool IsFullRank { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Gets or sets the eigen values (λ) of matrix in ascending value.
         /// </summary>
-        protected Vector<Complex> VectorEv
+        protected Vector<Complex> VectorEv { get; set; }
-        {
-            get;
-            set;
-        }
         
         /// <summary>
         /// Gets or sets eigenvectors.
         /// </summary>
-        protected Matrix<T> MatrixEv
+        protected Matrix<T> MatrixEv { get; set; }
-        {
-            get;
-            set;
-        }
 
         /// <summary>
         /// Gets or sets the block diagonal eigenvalue matrix.
         /// </summary>
-        protected Matrix<T> MatrixD
+        protected Matrix<T> MatrixD { get; set; }
-        {
-            get;
-            set;
-        }
-
-        /// <summary>
-        /// Internal method which routes the call to perform the singular value decomposition to the appropriate class.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>An EVD object.</returns>
-        internal static Evd<T> Create(Matrix<T> matrix)
-        {
-            if (typeof(T) == typeof(double))
-            {
-                var dense = matrix as Double.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Double.Factorization.DenseEvd(dense) as Evd<T>;
-                }
-
-                return new Double.Factorization.UserEvd(matrix as Matrix<double>) as Evd<T>;
-            }
-
-            if (typeof(T) == typeof(float))
-            {
-                var dense = matrix as Single.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Single.Factorization.DenseEvd(dense) as Evd<T>;
-                }
-
-                return new Single.Factorization.UserEvd(matrix as Matrix<float>) as Evd<T>;
-            }
-
-            if (typeof(T) == typeof(Complex))
-            {
-                var dense = matrix as LinearAlgebra.Complex.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex.Factorization.DenseEvd(dense) as Evd<T>;
-                }
-
-                return new LinearAlgebra.Complex.Factorization.UserEvd(matrix as Matrix<Complex>) as Evd<T>;
-            }
-
-            if (typeof(T) == typeof(Complex32))
-            {
-                var dense = matrix as LinearAlgebra.Complex32.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex32.Factorization.DenseEvd(dense) as Evd<T>;
-                }
-
-                return new LinearAlgebra.Complex32.Factorization.UserEvd(matrix as Matrix<Complex32>) as Evd<T>;
-            }
-
-            throw new NotSupportedException();
-        }
 
         /// <summary>Returns the eigen values as a <see cref="Vector{T}"/>.</summary>
         /// <returns>The eigen values.</returns>

src/Numerics/LinearAlgebra/Factorization/GramSchmidt.cs

@@ -28,12 +28,6 @@ using System;
 
 namespace MathNet.Numerics.LinearAlgebra.Factorization
 {
-    using Numerics;
-
-#if !NOSYSNUMERICS
-    using System.Numerics;
-#endif
-
     /// <summary>
     /// <para>A class which encapsulates the functionality of the QR decomposition Modified Gram-Schmidt Orthogonalization.</para>
     /// <para>Any real square matrix A may be decomposed as A = QR where Q is an orthogonal mxn matrix and R is an nxn upper triangular matrix.</para>
@@ -45,58 +39,5 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization
     public abstract class GramSchmidt<T> : QR<T>
     where T : struct, IEquatable<T>, IFormattable
     {
-        /// <summary>
-        /// Internal method which routes the call to perform the QR factorization to the appropriate class.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>A QR factorization object.</returns>
-        internal static GramSchmidt<T> Create(Matrix<T> matrix)
-        {
-            if (typeof(T) == typeof(double))
-            {
-                var dense = matrix as Double.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Double.Factorization.DenseGramSchmidt(dense) as GramSchmidt<T>;
-                }
-
-                return new Double.Factorization.UserGramSchmidt(matrix as Matrix<double>) as GramSchmidt<T>;
-            }
-
-            if (typeof(T) == typeof(float))
-            {
-                var dense = matrix as Single.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Single.Factorization.DenseGramSchmidt(dense) as GramSchmidt<T>;
-                }
-
-                return new Single.Factorization.UserGramSchmidt(matrix as Matrix<float>) as GramSchmidt<T>;
-            }
-
-            if (typeof(T) == typeof(Complex))
-            {
-                var dense = matrix as LinearAlgebra.Complex.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex.Factorization.DenseGramSchmidt(dense) as GramSchmidt<T>;
-                }
-
-                return new LinearAlgebra.Complex.Factorization.UserGramSchmidt(matrix as Matrix<Complex>) as GramSchmidt<T>;
-            }
-
-            if (typeof(T) == typeof(Complex32))
-            {
-                var dense = matrix as LinearAlgebra.Complex32.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex32.Factorization.DenseGramSchmidt(dense) as GramSchmidt<T>;
-                }
-
-                return new LinearAlgebra.Complex32.Factorization.UserGramSchmidt(matrix as Matrix<Complex32>) as GramSchmidt<T>;
-            }
-
-            throw new NotSupportedException();
-        }
     }
 }

src/Numerics/LinearAlgebra/Factorization/LU.cs

@@ -28,12 +28,6 @@ using System;
 
 namespace MathNet.Numerics.LinearAlgebra.Factorization
 {
-    using Numerics;
-
-#if !NOSYSNUMERICS
-    using System.Numerics;
-#endif
-
     /// <summary>
     /// <para>A class which encapsulates the functionality of an LU factorization.</para>
     /// <para>For a matrix A, the LU factorization is a pair of lower triangular matrix L and
@@ -56,74 +50,12 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization
         /// <summary>
         /// Gets or sets both the L and U factors in the same matrix.
         /// </summary>
-        protected Matrix<T> Factors
+        protected Matrix<T> Factors { get; set; }
-        {
-            get;
-            set;
-        }
 
         /// <summary>
         /// Gets or sets the pivot indices of the LU factorization.
         /// </summary>
-        protected int[] Pivots
+        protected int[] Pivots { get; set; }
-        {
-            get;
-            set;
-        }
-
-        /// <summary>
-        /// Internal method which routes the call to perform the LU factorization to the appropriate class.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <returns>An LU factorization object.</returns>
-        internal static LU<T> Create(Matrix<T> matrix)
-        {
-            if (typeof(T) == typeof(double))
-            {
-                var dense = matrix as Double.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Double.Factorization.DenseLU(dense) as LU<T>;
-                }
-
-                return new Double.Factorization.UserLU(matrix as Matrix<double>) as LU<T>;
-            }
-
-            if (typeof(T) == typeof(float))
-            {
-                var dense = matrix as Single.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Single.Factorization.DenseLU(dense) as LU<T>;
-                }
-
-                return new Single.Factorization.UserLU(matrix as Matrix<float>) as LU<T>;
-            }
-
-            if (typeof(T) == typeof(Complex))
-            {
-                var dense = matrix as LinearAlgebra.Complex.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex.Factorization.DenseLU(dense) as LU<T>;
-                }
-
-                return new LinearAlgebra.Complex.Factorization.UserLU(matrix as Matrix<Complex>) as LU<T>;
-            }
-
-            if (typeof(T) == typeof(Complex32))
-            {
-                var dense = matrix as LinearAlgebra.Complex32.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex32.Factorization.DenseLU(dense) as LU<T>;
-                }
-
-                return new LinearAlgebra.Complex32.Factorization.UserLU(matrix as Matrix<Complex32>) as LU<T>;
-            }
-
-            throw new NotSupportedException();
-        }
 
         /// <summary>
         /// Gets the lower triangular factor.
@@ -167,10 +99,7 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization
         /// <summary>
         /// Gets the determinant of the matrix for which the LU factorization was computed.
         /// </summary>
-        public abstract T Determinant
+        public abstract T Determinant { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Solves a system of linear equations, <b>AX = B</b>, with A LU factorized.

src/Numerics/LinearAlgebra/Factorization/QR.cs

@@ -28,12 +28,6 @@ using System;
 
 namespace MathNet.Numerics.LinearAlgebra.Factorization
 {
-    using Numerics;
-
-#if !NOSYSNUMERICS
-    using System.Numerics;
-#endif
-
     /// <summary>
     /// The type of QR factorization go perform.
     /// </summary>
@@ -69,85 +63,17 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization
         /// <summary>
         /// Gets or sets orthogonal Q matrix
         /// </summary>
-        protected Matrix<T> MatrixQ
+        protected Matrix<T> MatrixQ { get; set; }
-        {
-            get;
-            set;
-        }
 
         /// <summary>
         /// Gets or sets upper triangular factor R
         /// </summary>
-        protected Matrix<T> MatrixR
+        protected Matrix<T> MatrixR { get; set; }
-        {
-            get;
-            set;
-        }
 
         /// <summary>
         /// The QR factorization method.
         /// </summary>
-        protected QRMethod QrMethod 
+        protected QRMethod QrMethod { get; set; }
-        { 
-            get; 
-            set; 
-        }
-
-        /// <summary>
-        /// Internal method which routes the call to perform the QR factorization to the appropriate class.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <param name="method">The type of QR factorization to perform.</param>
-        /// <returns>A QR factorization object.</returns>
-        internal static QR<T> Create(Matrix<T> matrix, QRMethod method = QRMethod.Full)
-        {
-
-            if (typeof(T) == typeof(double))
-            {
-                var dense = matrix as Double.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Double.Factorization.DenseQR(dense, method) as QR<T>;
-                }
-
-                return new Double.Factorization.UserQR(matrix as Matrix<double>, method) as QR<T>;
-            }
-
-            if (typeof(T) == typeof(float))
-            {
-                var dense = matrix as Single.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Single.Factorization.DenseQR(dense, method) as QR<T>;
-                }
-
-                return new Single.Factorization.UserQR(matrix as Matrix<float>, method) as QR<T>;
-            }
-
-            if (typeof(T) == typeof(Complex))
-            {
-                var dense = matrix as LinearAlgebra.Complex.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex.Factorization.DenseQR(dense, method) as QR<T>;
-                }
-
-                return new LinearAlgebra.Complex.Factorization.UserQR(matrix as Matrix<Complex>, method) as QR<T>;
-            }
-
-            if (typeof(T) == typeof(Complex32))
-            {
-                var dense = matrix as LinearAlgebra.Complex32.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex32.Factorization.DenseQR(dense, method) as QR<T>;
-                }
-
-                return new LinearAlgebra.Complex32.Factorization.UserQR(matrix as Matrix<Complex32>, method) as QR<T>;
-            }
-
-            throw new NotSupportedException();
-        }
 
         /// <summary>
         /// Gets orthogonal Q matrix
@@ -174,19 +100,13 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization
         /// <summary>
         /// Gets the absolute determinant value of the matrix for which the QR matrix was computed.
         /// </summary>
-        public abstract T Determinant
+        public abstract T Determinant { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Gets a value indicating whether the matrix is full rank or not.
         /// </summary>
         /// <value><c>true</c> if the matrix is full rank; otherwise <c>false</c>.</value>
-        public abstract bool IsFullRank
+        public abstract bool IsFullRank { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Solves a system of linear equations, <b>AX = B</b>, with A QR factorized.

src/Numerics/LinearAlgebra/Factorization/Svd.cs

@@ -28,17 +28,11 @@
 // OTHER DEALINGS IN THE SOFTWARE.
 // </copyright>
 
-using MathNet.Numerics.Properties;
 using System;
+using MathNet.Numerics.Properties;
 
 namespace MathNet.Numerics.LinearAlgebra.Factorization
 {
-    using Numerics;
-
-#if !NOSYSNUMERICS
-    using System.Numerics;
-#endif
-
     /// <summary>
     /// <para>A class which encapsulates the functionality of the singular value decomposition (SVD).</para>
     /// <para>Suppose M is an m-by-n matrix whose entries are real numbers. 
@@ -60,128 +54,45 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization
         /// <summary>
         /// Gets or sets a value indicating whether to compute U and VT matrices during SVD factorization or not
         /// </summary>
-        protected bool ComputeVectors
+        protected bool ComputeVectors { get; set; }
-        {
-            get;
-            set;
-        }
 
         /// <summary>
         /// Gets or sets the singular values (Σ) of matrix in ascending value.
         /// </summary>
-        protected Vector<T> VectorS
+        protected Vector<T> VectorS { get; set; }
-        {
-            get;
-            set;
-        }
 
         /// <summary>
         /// Gets or sets left singular vectors (U - m-by-m unitary matrix)
         /// </summary>
-        protected Matrix<T> MatrixU
+        protected Matrix<T> MatrixU { get; set; }
-        {
-            get;
-            set;
-        }
 
         /// <summary>
         /// Gets or sets transpose right singular vectors (transpose of V, an n-by-n unitary matrix
         /// </summary>
-        protected Matrix<T> MatrixVT
+        protected Matrix<T> MatrixVT { get; set; }
-        {
-            get;
-            set;
-        }
 
         /// <summary>
         /// Gets the effective numerical matrix rank.
         /// </summary>
         /// <value>The number of non-negligible singular values.</value>
-        public abstract int Rank
+        public abstract int Rank { get; }
-        {
-            get;
-        }
-
-        /// <summary>
-        /// Internal method which routes the call to perform the singular value decomposition to the appropriate class.
-        /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
-        /// <param name="computeVectors">Compute the singular U and VT vectors or not.</param>
-        /// <returns>An SVD object.</returns>
-        internal static Svd<T> Create(Matrix<T> matrix, bool computeVectors)
-        {
-            if (typeof(T) == typeof(double))
-            {
-                var dense = matrix as Double.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Double.Factorization.DenseSvd(dense, computeVectors) as Svd<T>;
-                }
-
-                return new Double.Factorization.UserSvd(matrix as Matrix<double>, computeVectors) as Svd<T>;
-            }
-
-            if (typeof(T) == typeof(float))
-            {
-                var dense = matrix as Single.DenseMatrix;
-                if (dense != null)
-                {
-                    return new Single.Factorization.DenseSvd(dense, computeVectors) as Svd<T>;
-                }
-
-                return new Single.Factorization.UserSvd(matrix as Matrix<float>, computeVectors) as Svd<T>;
-            }
-
-            if (typeof(T) == typeof(Complex))
-            {
-                var dense = matrix as LinearAlgebra.Complex.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex.Factorization.DenseSvd(dense, computeVectors) as Svd<T>;
-                }
-
-                return new LinearAlgebra.Complex.Factorization.UserSvd(matrix as Matrix<Complex>, computeVectors) as Svd<T>;
-            }
-
-            if (typeof(T) == typeof(Complex32))
-            {
-                var dense = matrix as LinearAlgebra.Complex32.DenseMatrix;
-                if (dense != null)
-                {
-                    return new LinearAlgebra.Complex32.Factorization.DenseSvd(dense, computeVectors) as Svd<T>;
-                }
-
-                return new LinearAlgebra.Complex32.Factorization.UserSvd(matrix as Matrix<Complex32>, computeVectors) as Svd<T>;
-            }
-
-            throw new NotSupportedException();
-        }
 
         /// <summary>
         /// Gets the two norm of the <see cref="Matrix{T}"/>.
         /// </summary>
         /// <returns>The 2-norm of the <see cref="Matrix{T}"/>.</returns>
-        public abstract T Norm2
+        public abstract T Norm2 { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Gets the condition number <b>max(S) / min(S)</b>
         /// </summary>
         /// <returns>The condition number.</returns>
-        public abstract T ConditionNumber
+        public abstract T ConditionNumber { get; }
-        {
-            get;
-        }
 
         /// <summary>
         /// Gets the determinant of the square matrix for which the SVD was computed.
         /// </summary>
-        public abstract T Determinant
+        public abstract T Determinant { get; }
-        {
-            get;
-        }
 
         /// <summary>Returns the left singular vectors as a <see cref="Matrix{T}"/>.</summary>
         /// <returns>The left singular vectors. The matrix will be <c>null</c>, if <b>computeVectors</b> in the constructor is set to <c>false</c>.</returns>

src/Numerics/LinearAlgebra/Matrix.Arithmetic.cs

@@ -1032,7 +1032,7 @@ namespace MathNet.Numerics.LinearAlgebra
         /// <returns>effective numerical rank, obtained from SVD</returns>
         public virtual int Rank()
         {
-            return Svd<T>.Create(this, false).Rank;
+            return Svd(false).Rank;
         }
 
         /// <summary>Calculates the condition number of this matrix.</summary>
@@ -1040,7 +1040,7 @@ namespace MathNet.Numerics.LinearAlgebra
         /// <remarks>The condition number is calculated using singular value decomposition.</remarks>
         public virtual T ConditionNumber()
         {
-            return Svd<T>.Create(this, false).ConditionNumber;
+            return Svd(false).ConditionNumber;
         }
 
         /// <summary>Computes the determinant of this matrix.</summary>
@@ -1052,7 +1052,7 @@ namespace MathNet.Numerics.LinearAlgebra
                 throw new ArgumentException(Resources.ArgumentMatrixSquare);
             }
 
-            return LU<T>.Create(this).Determinant;
+            return LU().Determinant;
         }
 
         /// <summary>Computes the inverse of this matrix.</summary>
@@ -1064,7 +1064,7 @@ namespace MathNet.Numerics.LinearAlgebra
                 throw new ArgumentException(Resources.ArgumentMatrixSquare);
             }
 
-            return LU<T>.Create(this).Inverse();
+            return LU().Inverse();
         }
 
         /// <summary>

src/Numerics/LinearAlgebra/Single/ExtensionMethods.cs → src/Numerics/LinearAlgebra/Matrix.Factorization.cs

@@ -1,9 +1,11 @@
-// <copyright file="ExtensionMethods.cs" company="Math.NET">
+// <copyright file="Matrix.Factorization.cs" company="Math.NET">
 // Math.NET Numerics, part of the Math.NET Project
 // http://numerics.mathdotnet.com
 // http://github.com/mathnet/mathnet-numerics
 // http://mathnetnumerics.codeplex.com
-// Copyright (c) 2009-2010 Math.NET
+//
+// Copyright (c) 2009-2013 Math.NET
+//
 // Permission is hereby granted, free of charge, to any person
 // obtaining a copy of this software and associated documentation
 // files (the "Software"), to deal in the Software without
@@ -12,8 +14,10 @@
 // copies of the Software, and to permit persons to whom the
 // Software is furnished to do so, subject to the following
 // conditions:
+//
 // The above copyright notice and this permission notice shall be
 // included in all copies or substantial portions of the Software.
+//
 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
 // EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES
 // OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
@@ -25,75 +29,50 @@
 // </copyright>
 
 using MathNet.Numerics.LinearAlgebra.Factorization;
-using MathNet.Numerics.LinearAlgebra.Single.Factorization;
 
-namespace MathNet.Numerics.LinearAlgebra.Single
+namespace MathNet.Numerics.LinearAlgebra
 {
     /// <summary>
-    /// Extension methods which return factorizations for the various matrix classes.
+    /// Defines the base class for <c>Matrix</c> classes.
     /// </summary>
-    public static class ExtensionMethods
+    public abstract partial class Matrix<T>
     {
         /// <summary>
         /// Computes the Cholesky decomposition for a matrix.
         /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
         /// <returns>The Cholesky decomposition object.</returns>
-        public static Cholesky Cholesky(this Matrix<float> matrix)
+        public abstract Cholesky<T> Cholesky();
-        {
-            return (Cholesky)Cholesky<float>.Create(matrix);
-        }
 
         /// <summary>
         /// Computes the LU decomposition for a matrix.
         /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
         /// <returns>The LU decomposition object.</returns>
-        public static LU LU(this Matrix<float> matrix)
+        public abstract LU<T> LU();
-        {
-            return (LU)LU<float>.Create(matrix);
-        }
 
         /// <summary>
         /// Computes the QR decomposition for a matrix.
         /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
         /// <param name="method">The type of QR factorization to perform.</param>
         /// <returns>The QR decomposition object.</returns>
-        public static QR QR(this Matrix<float> matrix, QRMethod method = QRMethod.Full)
+        public abstract QR<T> QR(QRMethod method = QRMethod.Thin);
-        {
-            return (QR)QR<float>.Create(matrix, method);
-        }
 
         /// <summary>
         /// Computes the QR decomposition for a matrix using Modified Gram-Schmidt Orthogonalization.
         /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
         /// <returns>The QR decomposition object.</returns>
-        public static GramSchmidt GramSchmidt(this Matrix<float> matrix)
+        public abstract GramSchmidt<T> GramSchmidt();
-        {
-            return (GramSchmidt)GramSchmidt<float>.Create(matrix);
-        }
 
         /// <summary>
         /// Computes the SVD decomposition for a matrix.
         /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
         /// <param name="computeVectors">Compute the singular U and VT vectors or not.</param>
         /// <returns>The SVD decomposition object.</returns>
-        public static Svd Svd(this Matrix<float> matrix, bool computeVectors)
+        public abstract Svd<T> Svd(bool computeVectors);
-        {
-            return (Svd)Svd<float>.Create(matrix, computeVectors);
-        }
 
         /// <summary>
         /// Computes the EVD decomposition for a matrix.
         /// </summary>
-        /// <param name="matrix">The matrix to factor.</param>
         /// <returns>The EVD decomposition object.</returns>
-        public static Evd Evd(this Matrix<float> matrix)
+        public abstract Evd<T> Evd();
-        {
-            return (Evd)Evd<float>.Create(matrix);
-        }
     }
 }

src/Numerics/LinearAlgebra/Matrix.cs

@@ -1215,7 +1215,7 @@ namespace MathNet.Numerics.LinearAlgebra
         /// In a later release, it will be replaced with a sparse implementation.</remarks>
         public virtual T L2Norm()
         {
-            return Svd<T>.Create(this, false).Norm2;
+            return Svd(false).Norm2;
         }
 
         /// <summary>Calculates the Frobenius norm of this matrix.</summary>

src/Numerics/LinearAlgebra/Single/DenseMatrix.cs

@@ -33,6 +33,8 @@ using System.Collections.Generic;
 using System.Diagnostics;
 using System.Linq;
 using MathNet.Numerics.Distributions;
+using MathNet.Numerics.LinearAlgebra.Factorization;
+using MathNet.Numerics.LinearAlgebra.Single.Factorization;
 using MathNet.Numerics.LinearAlgebra.Storage;
 using MathNet.Numerics.Properties;
 using MathNet.Numerics.Providers.LinearAlgebra;
@@ -1032,5 +1034,35 @@ namespace MathNet.Numerics.LinearAlgebra.Single
 
             return (DenseMatrix)leftSide.Modulus(rightSide);
         }
+
+        public override Cholesky<float> Cholesky()
+        {
+            return new DenseCholesky(this);
+        }
+
+        public override LU<float> LU()
+        {
+            return new DenseLU(this);
+        }
+
+        public override QR<float> QR(QRMethod method = QRMethod.Thin)
+        {
+            return new DenseQR(this, method);
+        }
+
+        public override GramSchmidt<float> GramSchmidt()
+        {
+            return new DenseGramSchmidt(this);
+        }
+
+        public override Svd<float> Svd(bool computeVectors)
+        {
+            return new DenseSvd(this, computeVectors);
+        }
+
+        public override Evd<float> Evd()
+        {
+            return new DenseEvd(this);
+        }
     }
 }

src/Numerics/LinearAlgebra/Single/Matrix.cs

@@ -29,6 +29,8 @@
 // </copyright>
 
 using System;
+using MathNet.Numerics.LinearAlgebra.Single.Factorization;
+using MathNet.Numerics.LinearAlgebra.Factorization;
 using MathNet.Numerics.LinearAlgebra.Storage;
 using MathNet.Numerics.Properties;
 
@@ -456,5 +458,35 @@ namespace MathNet.Numerics.LinearAlgebra.Single
 
             return sum;
         }
+
+        public override Cholesky<float> Cholesky()
+        {
+            return new UserCholesky(this);
+        }
+
+        public override LU<float> LU()
+        {
+            return new UserLU(this);
+        }
+
+        public override QR<float> QR(QRMethod method = QRMethod.Thin)
+        {
+            return new UserQR(this, method);
+        }
+
+        public override GramSchmidt<float> GramSchmidt()
+        {
+            return new UserGramSchmidt(this);
+        }
+
+        public override Svd<float> Svd(bool computeVectors)
+        {
+            return new UserSvd(this, computeVectors);
+        }
+
+        public override Evd<float> Evd()
+        {
+            return new UserEvd(this);
+        }
     }
 }

src/Numerics/LinearAlgebra/Single/Solvers/Iterator.cs

@@ -39,6 +39,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Solvers
     /// </summary>
     public static class Iterator
     {
+
+        // TODO: Refactor
+
         /// <summary>
         /// Creates a default iterator with all the <see cref="IIterationStopCriterium{T}"/> objects.
         /// </summary>

src/Numerics/Numerics.csproj

@@ -122,6 +122,7 @@
     <Compile Include="Distributions\Wishart.cs" />
     <Compile Include="Distributions\Zipf.cs" />
     <Compile Include="LinearAlgebra\Builder.cs" />
+    <Compile Include="LinearAlgebra\Matrix.Factorization.cs" />
     <Compile Include="LinearAlgebra\Solvers\Iterator.cs" />
     <Compile Include="LinearAlgebra\Solvers\StopCriterium\IterationCountStopCriterium.cs" />
     <Compile Include="Providers\LinearAlgebra\Acml\AcmlLinearAlgebraProvider.Complex.cs" />
@@ -181,7 +182,6 @@
     <Compile Include="LinearAlgebra\Storage\SparseVectorStorage.cs" />
     <Compile Include="LinearAlgebra\Storage\DenseVectorStorage.cs" />
     <Compile Include="LinearAlgebra\Storage\MatrixStorage.Validation.cs" />
-    <Compile Include="LinearAlgebra\Complex32\ExtensionMethods.cs" />
     <Compile Include="LinearAlgebra\Complex32\Factorization\Cholesky.cs" />
     <Compile Include="LinearAlgebra\Complex32\DenseMatrix.cs" />
     <Compile Include="LinearAlgebra\Complex32\DiagonalMatrix.cs" />
@@ -199,7 +199,6 @@
       <SubType>Code</SubType>
     </Compile>
     <Compile Include="LinearAlgebra\Complex\DiagonalMatrix.cs" />
-    <Compile Include="LinearAlgebra\Complex\ExtensionMethods.cs" />
     <Compile Include="LinearAlgebra\Complex\Factorization\Cholesky.cs" />
     <Compile Include="LinearAlgebra\Complex\Factorization\Evd.cs" />
     <Compile Include="LinearAlgebra\Complex\Factorization\GramSchmidt.cs" />
@@ -210,7 +209,6 @@
     <Compile Include="LinearAlgebra\Complex\Solvers\Preconditioners\IPreConditioner.cs" />
     <Compile Include="LinearAlgebra\Complex\SparseMatrix.cs" />
     <Compile Include="LinearAlgebra\Complex\Vector.cs" />
-    <Compile Include="LinearAlgebra\Double\ExtensionMethods.cs" />
     <Compile Include="LinearAlgebra\Double\Factorization\Cholesky.cs" />
     <Compile Include="LinearAlgebra\Double\DenseMatrix.cs" />
     <Compile Include="LinearAlgebra\Double\DiagonalMatrix.cs" />
@@ -291,7 +289,6 @@
     <Compile Include="LinearAlgebra\Single\DenseMatrix.cs" />
     <Compile Include="LinearAlgebra\Single\DenseVector.cs" />
     <Compile Include="LinearAlgebra\Single\DiagonalMatrix.cs" />
-    <Compile Include="LinearAlgebra\Single\ExtensionMethods.cs" />
     <Compile Include="LinearAlgebra\Single\Factorization\Cholesky.cs" />
     <Compile Include="LinearAlgebra\Single\Factorization\DenseGramSchmidt.cs" />
     <Compile Include="LinearAlgebra\Single\Factorization\DenseEvd.cs" />