1c4762a1bSJed Brown static char help[] = "Tests relaxation for dense matrices.\n\n"; 2c4762a1bSJed Brown 3c4762a1bSJed Brown #include <petscmat.h> 4c4762a1bSJed Brown 5d71ae5a4SJacob Faibussowitsch int main(int argc, char **args) 6d71ae5a4SJacob Faibussowitsch { 7c4762a1bSJed Brown Mat C; 8c4762a1bSJed Brown Vec u, x, b, e; 9c4762a1bSJed Brown PetscInt i, n = 10, midx[3]; 10c4762a1bSJed Brown PetscScalar v[3]; 11c4762a1bSJed Brown PetscReal omega = 1.0, norm; 12c4762a1bSJed Brown 13327415f7SBarry Smith PetscFunctionBeginUser; 14*c8025a54SPierre Jolivet PetscCall(PetscInitialize(&argc, &args, NULL, help)); 159566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetReal(NULL, NULL, "-omega", &omega, NULL)); 169566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL)); 17c4762a1bSJed Brown 189566063dSJacob Faibussowitsch PetscCall(MatCreate(PETSC_COMM_SELF, &C)); 199566063dSJacob Faibussowitsch PetscCall(MatSetSizes(C, n, n, n, n)); 209566063dSJacob Faibussowitsch PetscCall(MatSetType(C, MATSEQDENSE)); 219566063dSJacob Faibussowitsch PetscCall(MatSetUp(C)); 229566063dSJacob Faibussowitsch PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &b)); 239566063dSJacob Faibussowitsch PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &x)); 249566063dSJacob Faibussowitsch PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &u)); 259566063dSJacob Faibussowitsch PetscCall(VecCreateSeq(PETSC_COMM_SELF, n, &e)); 269566063dSJacob Faibussowitsch PetscCall(VecSet(u, 1.0)); 279566063dSJacob Faibussowitsch PetscCall(VecSet(x, 0.0)); 28c4762a1bSJed Brown 299371c9d4SSatish Balay v[0] = -1.; 309371c9d4SSatish Balay v[1] = 2.; 319371c9d4SSatish Balay v[2] = -1.; 32c4762a1bSJed Brown for (i = 1; i < n - 1; i++) { 339371c9d4SSatish Balay midx[0] = i - 1; 349371c9d4SSatish Balay midx[1] = i; 359371c9d4SSatish Balay midx[2] = i + 1; 369566063dSJacob Faibussowitsch PetscCall(MatSetValues(C, 1, &i, 3, midx, v, INSERT_VALUES)); 37c4762a1bSJed Brown } 389371c9d4SSatish Balay i = 0; 399371c9d4SSatish Balay midx[0] = 0; 409371c9d4SSatish Balay midx[1] = 1; 419371c9d4SSatish Balay v[0] = 2.0; 429371c9d4SSatish Balay v[1] = -1.; 439566063dSJacob Faibussowitsch PetscCall(MatSetValues(C, 1, &i, 2, midx, v, INSERT_VALUES)); 449371c9d4SSatish Balay i = n - 1; 459371c9d4SSatish Balay midx[0] = n - 2; 469371c9d4SSatish Balay midx[1] = n - 1; 479371c9d4SSatish Balay v[0] = -1.0; 489371c9d4SSatish Balay v[1] = 2.; 499566063dSJacob Faibussowitsch PetscCall(MatSetValues(C, 1, &i, 2, midx, v, INSERT_VALUES)); 50c4762a1bSJed Brown 519566063dSJacob Faibussowitsch PetscCall(MatAssemblyBegin(C, MAT_FINAL_ASSEMBLY)); 529566063dSJacob Faibussowitsch PetscCall(MatAssemblyEnd(C, MAT_FINAL_ASSEMBLY)); 53c4762a1bSJed Brown 549566063dSJacob Faibussowitsch PetscCall(MatMult(C, u, b)); 55c4762a1bSJed Brown 56c4762a1bSJed Brown for (i = 0; i < n; i++) { 579566063dSJacob Faibussowitsch PetscCall(MatSOR(C, b, omega, SOR_FORWARD_SWEEP, 0.0, 1, 1, x)); 589566063dSJacob Faibussowitsch PetscCall(VecWAXPY(e, -1.0, x, u)); 599566063dSJacob Faibussowitsch PetscCall(VecNorm(e, NORM_2, &norm)); 609566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_SELF, "2-norm of error %g\n", (double)norm)); 61c4762a1bSJed Brown } 629566063dSJacob Faibussowitsch PetscCall(MatDestroy(&C)); 639566063dSJacob Faibussowitsch PetscCall(VecDestroy(&x)); 649566063dSJacob Faibussowitsch PetscCall(VecDestroy(&b)); 659566063dSJacob Faibussowitsch PetscCall(VecDestroy(&u)); 669566063dSJacob Faibussowitsch PetscCall(VecDestroy(&e)); 679566063dSJacob Faibussowitsch PetscCall(PetscFinalize()); 68b122ec5aSJacob Faibussowitsch return 0; 69c4762a1bSJed Brown } 70c4762a1bSJed Brown 71c4762a1bSJed Brown /*TEST 72c4762a1bSJed Brown 73c4762a1bSJed Brown test: 74c4762a1bSJed Brown 75c4762a1bSJed Brown TEST*/ 76