Jacbi_J

标签：xk Jacbi r1 r2 ++ float int

#include < stdio.h > 
#include < math.h > 
#define N 3
#define FO "%-10.5g"
void multm(float * a, float * b, int m, int n, int k, float * c) {
    int i,
    j,
    l,
    u;
    for (i = 0; i < m; i++) for (j = 0; j < k; j++) {
        u = i * k + j;
        c = 0.0;
        for (l = 0; l < n; l++) c += a * b[l * k + j];
    }
}
void scanfm(float * a, int m, int n, char aa) {
    int i,
    j;
    if (n == 1) for (i = 0; i < m; i++) {
        printf("%c[%d]=", aa, i + 1);
        scanf("%f", a + i);
    } else {
        for (i = 0; i < m; i++) for (j = 0; j < n; j++) {
            printf("%c[%d][%d]=", aa, i + 1, j + 1);
            scanf("%f", a + i * n + j);
        }
    }
}
void printm(float * a, int m, int n, char aa) {
    int i,
    j,
    k;
    printf("%c=\n", aa);
    for (i = 0; i < m; i++) {
        for (j = 0; j < n; j++) {
            k = i * n + j;
            if (fabs(a[k]) <= 1e-6) a[k] = 0.0;
            printf(FO, a[k]);
        }
        printf("\n");
    }
}
void copym(float * a, float * b, int m, int n) {
    int i,
    j,
    u;
    for (i = 0; i < m; i++) for (j = 0; j < n; j++) {
        u = i * n + j;
        b = a;
    }
}
float fanshu(float * a, int m, int n) {
    int i,
    j,
    u;
    float r = 0.0,
    r1;
    for (i = 0; i < m; i++) {
        r1 = 0.0;
        for (j = 0; j < n; j++) {
            u = i * n + j;
            r1 += fabs(a);
        }
        if (r1 > r) r = r1;
    }
    return r;
}
int Jacobi(float * a, int n, float * f, float * x, float eps, float * xk) {
    float r,
    r1,
    r2;
    int i,
    j,
    u;
    for (i = 0; i < n; i++) {
        u = i * n + i;
        if (fabs(a) <= 1e-5) return 0;
        for (j = 0; j < n; j++) if (j != i) a /= -a;
        f /= a;
        a = 0.0;
    }
    r1 = fanshu(a, n, n);
    if (r1 >= 1) return 0;
    r1 = r1 / (1.0 - r1);
    r = 1.0;
    for (i = 0; i < n; i++) xk = 0.0;
    while (r > eps) {
        for (i = 0; i < n; i++) {
            r2 = 0.0;
            for (j = 0; j < n; j++) if (j != i) r2 += a * xk[j];
            x = r2 + f;
        }
        for (i = 0; i < n; i++) xk = x - xk;
        r = r1 * fanshu(xk, n, 1);
        for (i = 0; i < n; i++) xk = x;
    }
    return 1;
}
int JGS(float * a, int n, float * f, float * x, float eps, float * xk, float * L, float * U) {
    float r,
    r1,
    r2;
    int i,
    j,
    u;
    for (i = 0; i < n; i++) {
        u = i * n + i;
        if (fabs(a) <= 1e-5) return 0;
        for (j = 0; j < n; j++) if (j != i) a /= -a;
        f /= a;
        a = 0.0;
    }
    for (i = 0; i < n; i++) {
        L = 1.0;
        for (j = 0; j < n; j++) {
            u = i * n + j;
            if (j < i) L = a;
            else if (j > i) U = a;
        }
    }
    multm(L, U, n, n, n, a);
    r1 = fanshu(a, n, n);
    if (r1 >= 1) return 0;
    r1 = r1 / (1.0 - r1);
    r = 1.0;
    for (i = 0; i < n; i++) xk = 0.0;
    while (r > eps) {
        for (i = 0; i < n; i++) {
            r2 = 0.0;
            for (j = 0; j < n; j++) {
                u = i * n + j;
                if (j > i) r2 += U * xk[j];
                else if (j < i) r2 += L * x[j];
            }
            x = r2 + f;
        }
        for (i = 0; i < n; i++) xk = x - xk;
        r = r1 * fanshu(xk, n, 1);
        for (i = 0; i < n; i++) xk = x;
    }
    return 1;
}
void main() {
    static float A[N][N],
    A1[N][N],
    L[N][N],
    U[N][N],
    f[N],
    f1[N],
    x[N],
    xk[N];
    float eps = 1e-6;
    int i;
    scanfm(A[0], N, N, 'A');
    scanfm(f, N, 1, 'b');
    printm(A[0], N, N, 'A');
    printm(f, N, 1, 'b');
    copym(A[0], A1[0], N, N);
    copym(f, f1, N, 1);
    i = Jacobi(A[0], N, f, x, eps, xk);
    if (i == 0) printf("Jacobi Failed!\n");
    else {
        printf("Jacobi:\n");
        printm(x, N, 1, 'x');
    }
    i = JGS(A1[0], N, f1, x, eps, xk, L[0], U[0]);
    if (i == 0) printf("JGS Failed!\n");
    else {
        printf("JGS:\n");
        printm(x, N, 1, 'x');
    }
}

标签：xk,Jacbi,r1,r2,++,float,int
From： https://blog.51cto.com/u_16076050/6195958