Browse Source

Sparse Linear Algebra: fix subtraction bug gh-20 gh-18

la-knuth
Christoph Ruegg 15 years ago
parent
commit
ea4e6a6d80
  1. 69
      src/Numerics/LinearAlgebra/Complex/SparseVector.cs
  2. 69
      src/Numerics/LinearAlgebra/Complex32/SparseVector.cs
  3. 69
      src/Numerics/LinearAlgebra/Double/SparseVector.cs
  4. 69
      src/Numerics/LinearAlgebra/Single/SparseVector.cs

69
src/Numerics/LinearAlgebra/Complex/SparseVector.cs

@ -550,9 +550,74 @@ namespace MathNet.Numerics.LinearAlgebra.Complex
return;
}
for (var index = 0; index < Count; index++)
var otherSparse = other as SparseVector;
if (otherSparse == null)
{
base.DoSubtract(other, result);
return;
}
var resultSparse = result as SparseVector;
if (resultSparse == null)
{
base.DoSubtract(other, result);
return;
}
// TODO (ruegg, 2011-10-11): Options to optimize?
if (ReferenceEquals(this, resultSparse))
{
int i = 0, j = 0;
while (i < NonZerosCount || j < otherSparse.NonZerosCount)
{
if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j])
{
_nonZeroValues[i++] -= otherSparse._nonZeroValues[j++];
}
else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j])
{
_nonZeroValues[i] -= otherSparse.At(_nonZeroIndices[i]);
i++;
}
else
{
var otherValue = otherSparse._nonZeroValues[j];
if (otherValue != Complex.Zero)
{
InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], -otherValue);
}
j++;
}
}
}
else
{
result.At(index, At(index) - other.At(index));
result.Clear();
int i = 0, j = 0, last = -1;
while (i < NonZerosCount || j < otherSparse.NonZerosCount)
{
if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j])
{
var next = _nonZeroIndices[i];
if (next != last)
{
last = next;
result.At(next, _nonZeroValues[i] - otherSparse.At(next));
}
i++;
}
else
{
var next = otherSparse._nonZeroIndices[j];
if (next != last)
{
last = next;
result.At(next, At(next) - otherSparse._nonZeroValues[j]);
}
j++;
}
}
}
}

69
src/Numerics/LinearAlgebra/Complex32/SparseVector.cs

@ -580,9 +580,74 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32
return;
}
for (var index = 0; index < Count; index++)
var otherSparse = other as SparseVector;
if (otherSparse == null)
{
base.DoSubtract(other, result);
return;
}
var resultSparse = result as SparseVector;
if (resultSparse == null)
{
base.DoSubtract(other, result);
return;
}
// TODO (ruegg, 2011-10-11): Options to optimize?
if (ReferenceEquals(this, resultSparse))
{
int i = 0, j = 0;
while (i < NonZerosCount || j < otherSparse.NonZerosCount)
{
if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j])
{
_nonZeroValues[i++] -= otherSparse._nonZeroValues[j++];
}
else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j])
{
_nonZeroValues[i] -= otherSparse.At(_nonZeroIndices[i]);
i++;
}
else
{
var otherValue = otherSparse._nonZeroValues[j];
if (otherValue != Complex32.Zero)
{
InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], -otherValue);
}
j++;
}
}
}
else
{
result.At(index, At(index) - other.At(index));
result.Clear();
int i = 0, j = 0, last = -1;
while (i < NonZerosCount || j < otherSparse.NonZerosCount)
{
if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j])
{
var next = _nonZeroIndices[i];
if (next != last)
{
last = next;
result.At(next, _nonZeroValues[i] - otherSparse.At(next));
}
i++;
}
else
{
var next = otherSparse._nonZeroIndices[j];
if (next != last)
{
last = next;
result.At(next, At(next) - otherSparse._nonZeroValues[j]);
}
j++;
}
}
}
}

69
src/Numerics/LinearAlgebra/Double/SparseVector.cs

@ -498,9 +498,74 @@ namespace MathNet.Numerics.LinearAlgebra.Double
return;
}
for (var index = 0; index < Count; index++)
var otherSparse = other as SparseVector;
if (otherSparse == null)
{
base.DoSubtract(other, result);
return;
}
var resultSparse = result as SparseVector;
if (resultSparse == null)
{
base.DoSubtract(other, result);
return;
}
// TODO (ruegg, 2011-10-11): Options to optimize?
if (ReferenceEquals(this, resultSparse))
{
int i = 0, j = 0;
while (i < NonZerosCount || j < otherSparse.NonZerosCount)
{
if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j])
{
_nonZeroValues[i++] -= otherSparse._nonZeroValues[j++];
}
else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j])
{
_nonZeroValues[i] -= otherSparse.At(_nonZeroIndices[i]);
i++;
}
else
{
var otherValue = otherSparse._nonZeroValues[j];
if (otherValue != 0.0)
{
InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], -otherValue);
}
j++;
}
}
}
else
{
result.At(index, At(index) - other.At(index));
result.Clear();
int i = 0, j = 0, last = -1;
while (i < NonZerosCount || j < otherSparse.NonZerosCount)
{
if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j])
{
var next = _nonZeroIndices[i];
if (next != last)
{
last = next;
result.At(next, _nonZeroValues[i] - otherSparse.At(next));
}
i++;
}
else
{
var next = otherSparse._nonZeroIndices[j];
if (next != last)
{
last = next;
result.At(next, At(next) - otherSparse._nonZeroValues[j]);
}
j++;
}
}
}
}

69
src/Numerics/LinearAlgebra/Single/SparseVector.cs

@ -528,9 +528,74 @@ namespace MathNet.Numerics.LinearAlgebra.Single
return;
}
for (var index = 0; index < Count; index++)
var otherSparse = other as SparseVector;
if (otherSparse == null)
{
base.DoSubtract(other, result);
return;
}
var resultSparse = result as SparseVector;
if (resultSparse == null)
{
base.DoSubtract(other, result);
return;
}
// TODO (ruegg, 2011-10-11): Options to optimize?
if (ReferenceEquals(this, resultSparse))
{
int i = 0, j = 0;
while (i < NonZerosCount || j < otherSparse.NonZerosCount)
{
if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j])
{
_nonZeroValues[i++] -= otherSparse._nonZeroValues[j++];
}
else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j])
{
_nonZeroValues[i] -= otherSparse.At(_nonZeroIndices[i]);
i++;
}
else
{
var otherValue = otherSparse._nonZeroValues[j];
if (otherValue != 0.0)
{
InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], -otherValue);
}
j++;
}
}
}
else
{
result.At(index, At(index) - other.At(index));
result.Clear();
int i = 0, j = 0, last = -1;
while (i < NonZerosCount || j < otherSparse.NonZerosCount)
{
if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j])
{
var next = _nonZeroIndices[i];
if (next != last)
{
last = next;
result.At(next, _nonZeroValues[i] - otherSparse.At(next));
}
i++;
}
else
{
var next = otherSparse._nonZeroIndices[j];
if (next != last)
{
last = next;
result.At(next, At(next) - otherSparse._nonZeroValues[j]);
}
j++;
}
}
}
}

Loading…
Cancel
Save