ArrayList<Node> nodes = new ArrayList<Node>( nodeSet.size() );
for ( Node node : nodeSet )
{
nodes.add( node );
}
DoubleMatrix hMatrix = new DoubleMatrix();
DoubleMatrix qMatrix = new DoubleMatrix();
for ( int i = 0; i < nodes.size(); ++i )
{
qMatrix.set( 0, i, values.get( nodes.get( i ) ) );
}
int localIterations = 1;
// The main arnoldi iteration loop
while ( true )
{
++totalIterations;
Map<Node,Double> newValues = new HashMap<Node,Double>();
// "matrix multiplication"
for ( Relationship relationship : relationshipSet )
{
if ( relationDirection.equals( Direction.BOTH )
|| relationDirection.equals( Direction.OUTGOING ) )
{
processRelationship( newValues, relationship, false );
}
if ( relationDirection.equals( Direction.BOTH )
|| relationDirection.equals( Direction.INCOMING ) )
{
processRelationship( newValues, relationship, true );
}
}
// Orthogonalize
for ( int j = 0; j < localIterations; ++j )
{
DoubleVector qj = qMatrix.getRow( j );
// vector product
double product = 0;
for ( int i = 0; i < nodes.size(); ++i )
{
Double d1 = newValues.get( nodes.get( i ) );
Double d2 = qj.get( i );
if ( d1 != null && d2 != null )
{
product += d1 * d2;
}
}
hMatrix.set( j, localIterations - 1, product );
if ( product != 0.0 )
{
// vector subtraction
for ( int i = 0; i < nodes.size(); ++i )
{
Node node = nodes.get( i );
Double value = newValues.get( node );
if ( value == null )
{
value = 0.0;
}
Double qValue = qj.get( i );
if ( qValue != null )
{
newValues.put( node, value - product * qValue );
}
}
}
}
double normalizeFactor = normalize( newValues );
values = newValues;
DoubleVector qVector = new DoubleVector();
for ( int i = 0; i < nodes.size(); ++i )
{
qVector.set( i, newValues.get( nodes.get( i ) ) );
}
qMatrix.setRow( localIterations, qVector );
if ( normalizeFactor == 0.0 || localIterations >= nodeSet.size()
|| localIterations >= iterations )
{
break;
}
hMatrix.set( localIterations, localIterations - 1, normalizeFactor );
++localIterations;
}
// employ the power method to find eigenvector to h
Random random = new Random( System.currentTimeMillis() );
DoubleVector vector = new DoubleVector();
for ( int i = 0; i < nodeSet.size(); ++i )
{
vector.set( i, random.nextDouble() );
}
MatrixUtil.normalize( vector );
boolean powerDone = false;
int its = 0;
double powerPrecision = 0.1;
while ( !powerDone )
{
DoubleVector newVector = MatrixUtil.multiply( hMatrix, vector );
MatrixUtil.normalize( newVector );
powerDone = true;
for ( Integer index : vector.getIndices() )
{
if ( newVector.get( index ) == null )
{
continue;
}
double factor = Math.abs( newVector.get( index )
/ vector.get( index ) );
if ( factor - powerPrecision > 1.0
|| factor + powerPrecision < 1.0 )
{
powerDone = false;
break;
}
}
vector = newVector;
++its;
if ( its > 100 )
{
break;
}
}
// multiply q and vector to get a ritz vector
DoubleVector ritzVector = new DoubleVector();
for ( int r = 0; r < nodeSet.size(); ++r )
{
for ( int c = 0; c < localIterations; ++c )
{
ritzVector.incrementValue( r, vector.get( c )
* qMatrix.get( c, r ) );
}
}
for ( int i = 0; i < nodeSet.size(); ++i )
{
values.put( nodes.get( i ), ritzVector.get( i ) );