forum.alglib.net
http://forum.alglib.net/

Matrix inversion produces unexpected results
http://forum.alglib.net/viewtopic.php?f=2&t=4466		 Page 1 of 1
Author:  epsPro [ Sun Jun 05, 2022 9:06 pm ]
Post subject:  Matrix inversion produces unexpected results
Please see the first code excerpt below along with an input matrix, Z, and an output matrix, Zinverse. Zinverse is produced by the alglib.rmatrixinverse; and it is obviously wrong. However, in the second code excerpt, alglib.rmatrixinverse produces a correct result. I do not see reason for the difference in the results. As far as I can tell, the problem is occuring within alglib.rmatrixinverse. Can you please help me understand (1) how this results difference is happening, and (2) how to diagnose issues within alglib methods?

Thanks in advance.





// search for the non-singular scaling matrix with largest termCount
termCount = 6; // (m / 4) + 2;
while (searching && termCount >= MIN_TERM_COUNT)
{
matrix.Multiply(Yt, Y, termCount, ref Z);
if (alglib.rmatrixdet(Z, termCount) != 0)
{
File.WriteAllText(Path.Combine(dataPath, $"Z_{termCount}.txt"), matrix.FormatMatrix(Z, "Z", fieldWidth));
searching = false;
// compute inverse
alglib.rmatrixinverse(ref Z, termCount, out aglibInfo, out aglibReport); // inverse returns in Z matrix
if (aglibInfo == 1) // inversion was successful
{
File.WriteAllText(Path.Combine(dataPath, $"Zinverse_{termCount}.txt"), matrix.FormatMatrix(Z, "Zinverse", fieldWidth));

matrix.Multiply(Z, Yt, ref S);
matrix.Multiply(S, x, ref a);
SetScalingConstants(a);
TermCount = termCount;
rv = true;
break;
}
else
throw new InvalidOperationException($"alglib.rmatrixinverse: aglibInfo = {aglibInfo}");
}
else
{
matrix.Initialize(ref Z, Matrix.InitializationTypes.Null);
termCount--;
}
}

File.WriteAllText(Path.Combine(dataPath, $"Y_{termCount}.txt"), matrix.FormatMatrix(Y, "Y", fieldWidth));
File.WriteAllText(Path.Combine(dataPath, $"Yt_{termCount}.txt"), matrix.FormatMatrix(Yt, "Yt", fieldWidth));
File.WriteAllText(Path.Combine(dataPath, $"S_{termCount}.txt"), matrix.FormatMatrix(S, "S", fieldWidth));


Z = { 6.00 -9.28 3.21 -1.81 -0.01 0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ -9.28 17.65 -6.75 3.21 0.01 -0.04 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 3.21 -6.75 2.69 -1.19 -0.01 0.02 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ -1.81 3.21 -1.19 0.60 0.00 -0.01 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ -0.01 0.01 -0.01 0.00 0.00 -0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.01 -0.04 0.02 -0.01 -0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }

Zinverse = { 10646.89 201458.61 164016.83 -660336.75 -4243385.60 4725122.52 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 201458.61 3918229.40 3163567.51 -12893771.96 -89957384.22 90987096.05 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 164016.83 3163567.51 2560885.83 -10397909.15 -70690760.75 73714030.96 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ -660336.75 -12893771.96 -10397909.15 42453443.94 299596879.81 -298962625.64 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ -4243385.60 -89957384.22 -70690760.75 299596879.82 2674360516.42 -2003109303.81 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 4725122.52 90987096.05 73714030.96 -298962625.64 -2003109303.73 2126085570.33 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }
{ 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 }


public void Test_ZInverse()
{
int info;
alglib.matinvreport rep;

Matrix matrix = new Matrix();
double[,] I;
double[,] P = null;
int termCount = 6;
double[,] X= null;
double[,] Z = null;
double[,] Zinverse = null;

initializeZInputSet(ref Zinverse, ref X);
if (alglib.rmatrixdet(Zinverse, termCount) != 0)
{
//lglib.rmatrixinverse(ref Zinverse, termCount, out info, out rep);
alglib.rmatrixinverse(ref Zinverse, termCount, out info, out rep);
Assert.True(info == 1);
}
else
Assert.True(false);

initializeZInputSet(ref Z, ref X);
matrix.Multiply(Z, Zinverse, termCount, ref P);
I = matrix.NewMatrix(termCount, termCount, Matrix.InitializationTypes.Identity);
Assert.True(matrix.IsEqual(P, I));
}
Page 1 of 1 All times are UTC
Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group
http://www.phpbb.com/