|
[ 2,-4, 2 | 6]
[ 1, 1, -1, | 0]
[ 0,-1, 3, | 8]
[ 0, 2, -6, |-16]
3 inconnues, 4 équations, 1 column augmented
GetSolutionState return Infinite.
The echelonned reduct matrix is:
( 1 0 0 2 )
( 0 1 0 1 )
( 0 0 1 3 )
( 0 0 0 0 )
The system has one solution: x=2, y=1, z=3
Not Infinite.
-- modified 1-Apr-20 12:15pm.
|
|
|
|
|
Hi, you are completely right. GetSolutionState is buggy when working with free variables. I have updated the code. Can you please check it and let me know your feedback.
private static MatrixEliminationResult FindSolution(MatrixEliminationResult result) {
double?[] solution = null;
MatrixSolutionState solutionState = MatrixSolutionState.None;
if (result.AugmentedColumnCount == 1) {
if (result.UnknownsCount > result.TotalRowCount)
solutionState = MatrixSolutionState.Infinite;
else
solutionState = MatrixSolutionState.Unique;
solution = new double?[result.UnknownsCount];
int pivotRow = 0;
for (int col = 0; col < result.UnknownsCount; col++) {
double unknownValue = result.FullMatrix[pivotRow, col];
double solutionValue = result.FullMatrix[pivotRow, result.TotalColumnCount - 1];
if (unknownValue == 1) {
solution[col] = solutionValue;
pivotRow++;
} else {
if (solutionState == MatrixSolutionState.Unique) {
if (solutionValue == 0)
solutionState = MatrixSolutionState.Infinite;
else
solutionState = MatrixSolutionState.None;
}
}
}
}
Mohammad Elsheimy
|
|
|
|
|
Hello,
I don't understand your code without comment.
Please test with the following set of examples.
Attention, for me, A = (A1 | A2), n = number of rows of A1, m = number of columns of A1, and p = number of columns of A2.
case 1 ; n=m=3, p=1 Cramer
Local $Ax=[ _
[ 2,-4, 2, 6], _ ; 1 0 0 | 2
[ 1, 1,-1, 0], _ ; 0 1 0 | 1
[ 0,-1, 3, 8]] ; 0 0 1 | 3 rg=3, det=14, sol UNIQUE (2,1,3)
case 8 ; n=m=3 col2=col0+col1
Local $Ax=[ _
[ 2,-4, -2, 0], _ ; 1 0 1 0
[ 1, 1, 2, 0], _ ; 0 1 1 0
[ 0,-1, -1, 0]] ; 0 0 0 0
;Infinité de solutions (-x2, -x2, x2) où x2 réel quelconque
;Rang(A1)=2
;Déterminant de min(n,m)=0
case 9 ; n=m wikipedia Cramer
Local $Ax=[ _
[ 1,-1, 2, 5], _ ; 1 0 0 1
[ 3, 2, 1, 10], _ ; 0 1 0 2
[ 2,-3,-2, -10]] ; 0 0 1 3 rang=3, det=35, sol UNIQUE (1,2,3)
case 15 ; m=n p=1 3 équations, 3 inconnues col2=col0+col1
Local $Ax=[ _
[ 2,-4, -2, 0], _ ; 1 0 1 0
[ 1, 1, 2, 0], _ ; 0 1 1 0
[ 0,-1, -1, 1]] ; 0 0 0 1 rang=2, det=0 None
; f: R^3-->R^3 dim Im(f) = 2 base V0(2, 1, 0),V1(-4,1,-1)
; Famille liée: V2=V0+V1
case 17 ; n=m=6 p=1
Local $Ax = [ _
[0, 0, 0, 0, 18, 18, 18], _
[0, 0, -3, -3, -3, -3, -3], _
[6, 6, 3, 3, 21, 21, 21], _
[0, 0, 0, 0, 0, 0, 0], _
[0, 0, -9, -9, 9, 9, 9], _
[0, 0, 0, 0, 18, 18, 18]]
;( (1) 1 0 0 0 0 | 0 ) L0
;( 0 0 (1) 1 0 0 | 0 ) L1
;( 0 0 0 0 (1) 1 | 1 ) L2
;( 0 0 0 0 0 0 | 0 ) L3
;( 0 0 0 0 0 0 | 0 ) L4
;( 0 0 0 0 0 0 | 0 ) L5
; sol = ( -x1 , x1, -x3 , x3, 1-x5 , x5)
case 18 ;n=m=4 les pivots ne sont pas tous dans les colonnes 0 et 1 ( <rang)
Local $Ax = [ _
[0, -54, 0, 0, -6], _
[0, 0, 18, 18, 18], _
[0, -108, -18, -18, -27], _
[0, -162, -36, -36, -51]]
;( 0 1 0 0 | 0.093 )
;( 0 0 1 1 | 1 )
;( 0 0 0 0 | 1 )
;( 0 0 0 0 | -1 )
; pas de solution
;Rang(A1)=2
;Déterminant de min(n,m)=0
case 19 ; m=n=4
Local $Ax = [ _
[0, 12, 120, 120, 120], _
[0, 24, 240, 240, 240], _
[0, 72, 720, 720, 720], _
[0, 6, 6, 6, 6]]
;( (0) (1) 0 0 | 0 )
;( 0 0 (1) 1 | 1 )
;( 0 0 0 (0) | 0 )
;( 0 0 0 0 | 0 )
;Infinité de solutions (x0, 0, 1-x3, x3) où x0, x3 réels quelconques
;Rang(A1)=2
;Déterminant de min(n,m)=0
case 20
Local $Ax = [ _
[0, 0, 108, 108, 108], _
[0, 0, 216, 216, 216], _
[0, 0, 648, 648, 648], _
[0, 0, 0, 0, 0]]
#cs
( 0 0 1 1 | 1 )
( 0 0 0 0 | 0 )
( 0 0 0 0 | 0 )
( 0 0 0 0 | 0 )
Infinité de solutions (x0, x1, 1-x3, x3) où x0, x1, x3 réels quelconques
Rang(A1)=1
Déterminant de min(n,m)=0
#ce
case 21
Local $Ax = [ _
[0, 0, 0, 0, 0], _
[0, 0, 0, 0, 0], _
[0, 0, 0, 0, 0], _
[0, 0, 0, 0, 0]]
#cs
Tout élément de |R^4 est solution du système
Rang(A1)=0
Déterminant de min(n,m)=0
#ce
case 2 ; n>m
Local $Ax=[ _
[ 2,-4, 2, 6], _ ; 1 0 0 2
[ 1, 1,-1, 0], _ ; 0 1 0 1
[ 0,-1, 3, 8], _ ; 0 0 1 3
[ 1, 1, 1, 1]] ; 0 0 0 -5 rang=3, det=14, None
; f:R^3 --->R^4 dim Im(f)=3 base: V0, V1, V2
; Famille libre
case 3 ; n>m
Local $Ax=[ _
[ 2,-4, 2, 6], _ ; 1 0 0 2
[ 1, 1,-1, 0], _ ; 0 1 0 1
[ 0,-1, 3, 8], _ ; 0 0 1 3
[ 1, 1, 1, 1], _ ; 0 0 0 -5
[ 0, 0, 0, 0]] ; 0 0 0 0 rang=3, det=14, None, FREE
case 4 ; inverse d'une matrice p=3
Local $Ax=[ _
[ 2,-1, 0, 1, 0, 0], _ ; 1 0 0 .75 .5 .25
[-1, 2,-1, 0, 1, 0], _ ; 0 1 0 .5 1 .5
[ 0,-1, 2, 0, 0, 1]] ; 0 0 1 .25 .5 .75 det=4=2(3)-(-1)(-2)
$p=3 ; 3 colonnes ajoutées
case 5 ; n>m 2 éq identiques, f: |R^3 ----> |R^4
Local $Ax=[ _
[ 2,-4, 2, 6], _ ; 1 0 0 2
[ 1, 1,-1, 0], _ ; 0 1 0 1
[ 0,-1, 3, 8], _ ; 0 0 1 3
[ 0, 2, -6, -16]] ; 0 0 0 0 Unique (2,1,3) rang=3 det=28
; Sol (2,1,3)
; f:R^3 --->R^4 dim Im(f)=3 base: V0, V1, V2
; Famille libre
case 6 ; n>m 4 équations, 3 inconnues, -16 devient -15
Local $Ax=[ _
[ 2,-4, 2, 6], _ ; 1 0 0 1.929
[ 1, 1,-1, 0], _ ; 0 1 0 0.857
[ 0,-1, 3, 8], _ ; 0 0 1 2.786
[ 0, 2, -6, -15]] ; 0 0 0 0.5 Rang=3, det=28, None
case 7 ; n<m, col3=col0+col1
Local $Ax=[ _
[ 2,-4, 2,-2, 6], _ ; 1 0 0 1 0.714
[ 1, 1,-1, 2, 0], _ ; 0 1 0 1 -1.571
[ 0,-1, 3,-1, -1]] ; 0 0 1 0 -0.857 rang=3, det=14, Infinite
; Sol part. (0.7, -1.6, -0.9, 0) sol géné (0.7-x3, -1.6-x3 , -0.9, x3)
; f:R^4 --->R^3 dim Im(f)=3 base: V0(2,1,0), V1(-4,1,-1), V2(2,-1,3)
; Famille liée: V3=V0+V1
case 10 ; n<m wikipedia OK
Local $Ax=[ _
[ 1, 2, 2, -3, 2, 3], _ ; 1 2 0 1 -4 -5
[ 2, 4, 1, 0, -5, -6], _ ; 0 0 1 -2 3 4
[ 4, 8, 5,-6, -1, 0], _ ; 0 0 0 0 0 0
[-1,-2,-1, 1, 1, 1]] ; 0 0 0 0 0 0 infinite, rang=2 det=0
; Sol particulière: (-5,0,4,0,0), générale (-5-2x1-x3+4x4 ,x1, 4+2x3-3x4 ,x3,x4)
; f:R^5 --->R^4 dim Im(f)=2 base: V0(1,2,4,-1), V2(2,1,5,-1)
; Famille liée, V1=2V0, V3=V0-2V2, V4=-4V0+3V2
case 11 ; n<m wikipedia OK
Local $Ax=[ _
[ 1, 2, 2, -3, 2, 3], _ ; 1 2 0 1 -4 0
[ 2, 4, 1, 0, -5, -6], _ ; 0 0 1 -2 3 0
[ 4, 8, 5,-6, -1, 0], _ ; 0 0 0 0 0 1
[-1,-2,-1, 1, 1, 2]] ; 0 0 0 0 0 0 None, rang 2=deux lignes non nulles dans 4x5
; f:R^5 --->R^4 dim Im(f)=2 base: V0(1,2,4,-1), V2(2,1,5,-1)
; Famille liée, V1=2V0, V3=V0-2V2, V4=-4V0+3V2
case 12 ; p=0 Déterminant
Local $Ax=[ _
[2, -1, 0], _
[0, -1, 2], _
[-1, 2,-1]] ; det=-4 , rang=3
case 13 ; n<m
Local $Ax=[ _
[ 1, 2, 0, 1,-4, -5], _ ; 1 2 0 0 0 -1
[ 0, 0, 1,-2, 3, 4], _ ; 0 0 1 0 0 1
[ 0, 0, 0, 1, 0, 0], _ ; 0 0 0 1 0 0
[ 0, 0, 0, 0, 2, 2]] ; 0 0 0 0 1 1
#cs
Infinité de solutions (-1-2 x1, x1, 1, 0, 1) où x1 réel quelconque
Rang(A1)=4
Déterminant de min(n,m)=0
; f:R^5 --->R^4 dim Im(f)=4 base: V0(1,0,0,0), V2(0,1,0,0), V3(1,-2,1,0), V4(-4,3,0,2)
; Famille liée, V1=2V0
#ce
case 14 ; n>m
Local $Ax=[ _
[ 1, 2, 0, 1, -5], _ ;(1)0 -2 0 -13
[ 0, 1, 1, -2, 4], _ ; 0(1) 1 0 4
[ 0, 0, 0, 1, 0], _ ; 0 0 0 (1) 0
[ 0, 0, 0, 0, 2], _ ; 0 0 0 0 2
[ 0, 0, 0, 1, 0]] ; 0 0 0 0 0 None, rang=3 det=0
; f:R^4 --->R^5 dim Im(f)=3 base: V0(1,0,0,0,0), V1(2,1,0,0), V3(1,-2,1,0,1)
; Vect colonnes=Famille liée, V2=-2V0 + 1V3
case 16 ; n<m p=0 wikipedia EN f:R^9 ----> R^6 rang=5=dim Im(f)
Local $Ax=[ _
[1,-1, 2,-3, 4,-5, 6,-7, 8], _ ;(1) -1 0 0 79 -93 0 -121 0
[0, 0, 1,-9, 0,-1, 2,-3, 4], _ ; 0 0 (1) 0 -45 53 0 69 0
[0, 0, 0, 1,-5, 6,-7, 8,-9], _ ; 0 0 0 (1) -5 6 0 8 0
[0, 0, 0, 0, 0, 0, 1, 0,-1], _ ; 0 0 0 0 0 0 (1) 0 0
[0, 0, 0, 0, 0, 0, 0, 0, 1], _ ; 0 0 0 0 0 0 0 0 (1)
[0, 0, 0, 0, 0, 0, 0, 0, 0]] ; 0 0 0 0 0 0 0 0 0 rang=5, det=0, None
$p=0
|
|
|
|
|
Thanks for the interesting piece. Shouldn’t ‘Laprace Expansion method’ be Laplace Expansion method? Also some of your methods require both the matrix row count and the column count as parameters whereas others use int Array.GetLength(int dimension) method to determine the required values within the method body. It may be a good idea to refractor the former methods to reduce the number of parameters.
|
|
|
|
|
1. You are completely right. I have updated the article with the correct spelling.
2. The reason I have added lengths as input arguments is that I found that calling GetLength multiple times incurs much performance overhead. I don't know exactly how GetLength is implemented (it's a natively-implemented function), however, I tested the code with and without GetLength in many scenarios and found that caching values of GetLength saved a minimum of 1 third of the total processing time!
Mohammad Elsheimy
|
|
|
|
|
This article is very useful to me, since I work in a Scientific field. I am excited to try some of your routines. If you have an update, perhaps you could consider some CPU time comparisons, between C# with garbage collection, and C++ using static variables.
|
|
|
|
|
Sounds great. what about Python?
---
Mohammad Elsheimy
http://JustLikeAMagic.com
|
|
|
|
|
Relative to your Projection function, I would prefer something like
public static double project_vec_on_vec(double[] vec0, double[] vec1)
{
return DotProduct(vec0, vec1) / Sqrt(DotProduct(vec1, vec1);
}
It makes clear which vector is being projected and avoids an extra function call.
Louis T Klauder Jr, Marshall, MI, lou@klauder.org
|
|
|
|
|