Hi,

For some reason, I can't get this program to work. I've had other CS majors look at it and they can't figure it out either.

This program performs the Jacobi algorithm (you can see step-by-step instructions and a MATLAB implementation here). BTW, it's different from the Wikipedia article of the same name.

Since NSArray is one-dimensional, I added a method that makes it act like a two-dimensional C array. After running the Jacobi algorithm many times, the diagonal entries in the NSArray (i[0][0], i[1][1], etc.) are supposed to get bigger and the others approach 0. For some reason though, they all increase exponentially. For instance, i[2][4] should equal 0.0000009, not 9999999, while i[2][2] should be big.

Thanks in advance,

Chris

NSArray+Matrix.m

@implementation NSArray (Matrix)

@dynamic offValue, transposed;

- (double)offValue {
    double sum = 0.0;
    for ( MatrixItem *item in self )
        if ( item.nonDiagonal )
            sum += pow( item.value, 2.0 );
    return sum;
}

- (NSMutableArray *)transposed {
    NSMutableArray *transpose = [[[NSMutableArray alloc] init] autorelease];
    int i, j;

    for ( i = 0; i < 5; i++ ) {
        for ( j = 0; j < 5; j++ ) {
            [transpose addObject:[self objectAtRow:j andColumn:i]];
        }
    }
    return transpose;
}

- (id)objectAtRow:(NSUInteger)row andColumn:(NSUInteger)column {
    NSUInteger index = 5 * row + column;
    return [self objectAtIndex:index];
}

- (NSMutableArray *)multiplyWithMatrix:(NSArray *)array {
    NSMutableArray *result = [[NSMutableArray alloc] init];
    int i = 0, j = 0, k = 0;
    double value;

    for ( i = 0; i < 5; i++ ) {
        value = 0.0;
        for ( j = 0; j < 5; j++ ) {
            for ( k = 0; k < 5; k++ ) {
                MatrixItem *firstItem = [self objectAtRow:i andColumn:k];
                MatrixItem *secondItem = [array objectAtRow:k andColumn:j];
                value += firstItem.value * secondItem.value;
            }
            MatrixItem *item = [[MatrixItem alloc] initWithValue:value];
            item.row = i;
            item.column = j;
            [result addObject:item];
        }
    }
    return result;
}
@end

Jacobi_AlgorithmAppDelegate.m

// ...

- (void)jacobiAlgorithmWithEntry:(MatrixItem *)entry {
    MatrixItem *b11 = [matrix objectAtRow:entry.row andColumn:entry.row];
    MatrixItem *b22 = [matrix objectAtRow:entry.column andColumn:entry.column];

    double muPlus = ( b22.value + b11.value ) / 2.0;
    muPlus += sqrt( pow((b22.value - b11.value), 2.0) + 4.0 * pow(entry.value, 2.0) );

    Vector *u1 = [[[Vector alloc] initWithX:(-1.0 * entry.value) andY:(b11.value - muPlus)] autorelease];
    [u1 normalize];
    Vector *u2 = [[[Vector alloc] initWithX:-u1.y andY:u1.x] autorelease];

    NSMutableArray *g = [[[NSMutableArray alloc] init] autorelease];
    for ( int i = 0; i <= 24; i++ ) {
        MatrixItem *item = [[[MatrixItem alloc] init] autorelease];
        if ( i == 6*entry.row )
            item.value = u1.x;
        else if ( i == 6*entry.column )
            item.value = u2.y;
        else if ( i == ( 5*entry.row + entry.column ) || i == ( 5*entry.column + entry.row ) )
            item.value = u1.y;
        else if ( i % 6 == 0 )
            item.value = 1.0;
        else
            item.value = 0.0;

        [g addObject:item];
    }

    NSMutableArray *firstResult = [[g.transposed multiplyWithMatrix:matrix] autorelease];

    matrix = [firstResult multiplyWithMatrix:g];
}

// ...