diff options
author | Ilhan Özgen <ilhan.oezgen@wahyd.tu-berlin.de> | 2017-05-02 14:34:39 +0200 |
---|---|---|
committer | Ilhan Özgen <ilhan.oezgen@wahyd.tu-berlin.de> | 2017-05-02 14:34:39 +0200 |
commit | 7108b261cd088f6287879af3d22d61e0644de4fd (patch) | |
tree | 5e57ec196826406c80ff914f5276a872f5acf105 | |
parent | 83ba0f455ddd64a9042d1af4fccf36decc8d0056 (diff) |
polish solver and test
-rw-r--r-- | solver.c | 8 | ||||
-rw-r--r-- | test_pcg.c | 31 |
2 files changed, 17 insertions, 22 deletions
@@ -70,14 +70,12 @@ int pcg(struct diag a, double* b, double* x0, struct diag c, double tol) } z0 = multiply(c, r0); - + memcpy(p, z0, len_c * sizeof(double)); // loop while (count < 100) { - printf(">>> %d\n", count); - alph = dot(r0, z0, len_c) / dot(p, multiply(a, p), len_c); t = multiply(a, p); @@ -114,8 +112,8 @@ int pcg(struct diag a, double* b, double* x0, struct diag c, double tol) free(r1); free(t); - printf(">>> iteration complete w/ err : %g\n", eps); - printf(">>> iteration complete w/ flag: %d\n", flag); + printf(">>> (pcg) :: iteration complete w/ err : %g, flag: %d, steps: %d\n", + eps, flag, count); return flag; @@ -9,6 +9,7 @@ #include <stdio.h> #include <stdlib.h> +#include <math.h> #include "global.h" #include "linalg.h" #include "solver.h" @@ -21,35 +22,31 @@ int main (void) int ids[3] = { -1, 0, 1 }; struct diag m = { a, ids, 2, 2, 3 }; - - for (int i = 0; i < 2; i ++) { - for (int j = 0; j < 2; j ++) { - printf(" %g ", get(m, i, j)); - } - printf("\n"); - } - printf("\n"); double c[2] = { 1.0 / 4.0, 1.0 / 3.0 }; int cids[1] = { 0 }; struct diag pre = { c, cids, 2, 2, 1 }; // Jacobi preconditioner - - for (int i = 0; i < 2; i ++) { - for (int j = 0; j < 2; j ++) { - printf(" %g ", get(pre, i, j)); - } - printf("\n"); - } - printf("\n"); double b[3] = { 1.0, 2.0 }; double x0[2] = { 10000.0, 1000.0 }; - printf(">>> iteration finished w/ flag: %d\n", pcg(m, b, x0, pre, 1.0e-3)); + pcg(m, b, x0, pre, 1.0e-12); + printf(">>> x0: \n"); for (int i = 0; i < 2; i ++) { printf(">>> %g\n", x0[i]); } + if (fabs(x0[0] - (1.0/11.0)) <= 1.0e-10) { + printf(">>> test passed.\n"); + } else { + printf(">>> test failed. expected: %g, found: %g\n", (1.0/11.0), x0[0]); + } + + if (fabs(x0[1] - (7.0/11.0)) <= 1.0e-10) { + printf(">>> test passed.\n"); + } else { + printf(">>> test failed. expected: %g, found: %g\n", (7.0/11.0), x0[1]); + } } |