Seidel_0.8

/*
 * =====================================================================================
 *
 *       Filename:  exp_seidel.c
 *
 *    Description:  Solve example in page 74 - GT GTS LETRONGVINH
 *
 *        Version:  0.8
 *        Created:  06/09/2011 10:02:20 PM
 *       Revision:  none
 *       Compiler:  gcc
 *
 *         Author:  FamiHug (fh), [email protected]
 *        Company:  familug.blogspot.com
 *
 * =====================================================================================
 */

//Todo List
//nhap ma tran
//viet ve dang alphaX + beta
//tach alpha thanh alpha 1 + alpha2
//thuc hien lap seidel

#include <stdio.h>
#include <string.h>
#include <stdlib.h>


//Nhap ma tran

void NhapMaTran(float *matrix, int size)
{
    int i,j;

    printf("Nhap vao ma tran: \n");
    for(i = 0; i < size; i++)
    {
        for(j = 0; j < size + 1; j++)
        {
            scanf("%f",matrix + i * (size + 1)  + j);
        }
    }
    printf("Nhap xong!\n");
    return;
}

void InDeBai(float *matrix, int size)
{
    int i,j;

    printf("Ma Tran: \n");
    for(i = 0; i < size; i++)
    {
        for(j = 0; j < size + 1; j++)
        {
            printf("%3.3f ",*(matrix + i * (size + 1) + j));
        }
        printf("\n");
    }
    return;
}


void InMaTranVuong(float *matrix, int size)
{
    int i,j;

    for(i = 0; i < size; i++)
    {
        for(j = 0; j < size; j++)
        {
            printf("%3.3f ",*(matrix + i * size + j));
        }
        printf("\n");
    }
    return;
}

void InBeta(float *beta, int size)
{
    int t;
    printf("Beta\n");
    for(t = 0; t < size; t++)
        printf("%f\n",*(beta + t));

    return;
}
//Viet lai theo dang alphaX
//Ver 0.8: tach thah 2 matran alpha va beta reing biet. VIet lai theo kieu void
void TaoAlphaBeta(float *matrix, float *alpha, float *beta)
{
    int i,j;
    int size = 3;

    for(i = 0; i < size; i++)
    {
        for(j = 0; j < size; j++)
        {
            if(i != j) *(alpha + i * size + j) = -(*(matrix + i * (size + 1) + j))/(*(matrix + i * (size + 1) + i));
        }
    }

    for(i = 0; i < size; i++)
    {
        //Cac phan tu cua Beta
        (*(beta + i)) = (*(matrix + i * (size + 1) + size)) / (*(matrix + i * (size + 1) + i));
        //Alpha[i][i] = 0
        (*(alpha + i * size + i)) = 0;
    }
    return;
}

//tach alpha thanh alpha1 + alpha 2
void TachAlpha(float *alpha, float *al1, float *al2)
{
    int i,j;
    int size = 3;

    //Copy alpha -> al1, al2
    memcpy(al1, alpha, sizeof(float) * size * size);
    memcpy(al2, alpha, sizeof(float) * size * size);

    for(i = 0; i < size; i++)
    {
        for(j = 0; j < size; j++)
        {
    //xoa nua tren duong cheo cua alpha de dc al1
            if(j >= i) *(al1 + i * size + j) = 0;
    //Xoa nua duoi duong cheo cua matran alpha de dc al2
            else *(al2 + i * size + j) = 0;
        }
    }

    return;
}


//Kiem tra tinh lap don

float absolute(float x)
{
    if(x >= 0) return x;
    else return -x;
}

int isLapDon(float **alpha)
{
    int i,j;
    int size = 3;
    float max = 0;
    float tempSum = 0;
    for(i = 0; i < size; i++)
    {
        tempSum = 0;
        for(j = 0; j < size; j++)
        {
            tempSum += absolute(*(*alpha + i * (size + 1) + j));
        }
        if(max < tempSum) max = tempSum;
    }
    printf("Max: %f\n",max);
    if(max < 1) return 1;
    else return 0;
}


//Chon xap xi dau, so lan lap k
//Thuc hien lap va dua ra ket qua

int main(int argc, char **argv)
{
    int size = 3;
    float *matrix;
    float *alpha, *al1, *al2, *beta;

    //malloc for matrix "size x (size + 1)", the matrix from the system of equations
    matrix = (float *)malloc(sizeof(float) * size * (size + 1));
    //malloc for matrix "size x size"
    alpha  = (float *)malloc(sizeof(float) * size * size);
    al1  = (float *)malloc(sizeof(float) * size * size);
    al2  = (float *)malloc(sizeof(float) * size * size);
    //malloc for matrix "size x 1"
    beta = (float *)malloc(sizeof(float) * size * 1);

    NhapMaTran(matrix,size);
    InDeBai(matrix,size);
    TaoAlphaBeta(matrix,alpha,beta);
    printf("Alpha: \n");
    InMaTranVuong(alpha,size);
    printf("Co tinh lap don: %d\n",isLapDon(&alpha));

    TachAlpha(alpha,al1,al2);

    printf("Alpha1: \n");
    InMaTranVuong(al1,size);
    printf("Alpha2: \n");
    InMaTranVuong(al2,size);

    InBeta(beta,size);

    //Nhap vao xap xi dau
    int k = 4; //k la so lan lap
    int i,g;

    //Ham nay can viet lai cho cac truong hop size khac 3
    float R[3][10];
    R[0][0] = 0, R[1][0] = 0, R[2][0] = 0;

    for(g = 1; g < k + 1; g++)
    {
        for(i = 0; i < size; i++)
        {
            R[i][g] = R[0][g - 1] * (*(al1 + i * size + 0)) + R[1][g - 1] * (*(al1 + i * size + 1)) + R[2][g - 1]  * (*(al1 + i * size + 2)) + (*(al2 + i * size + 0)) + R[1][g - 1] * (*(al2 + i * size + 1)) + R[2][g - 1]  * (*(al2 + i * size + 2)) + *(beta + i + 0);
        }
    }

    printf("Ket Qua: \n");

    for(i = 0; i < size; i++)
    {
        for(g = 0; g < k + 1; g++)
        {
            printf("%f ",R[i][g]);
        }
        printf("\n");
    }

    return 0;
}