Actual source code: ex71.c
1: static char help[] = "This example illustrates the use of PCBDDC/FETI-DP with 2D/3D DMDA.\n\
2: It solves the constant coefficient Poisson problem or the Elasticity problem \n\
3: on a uniform grid of [0,cells_x] x [0,cells_y] x [0,cells_z]\n\n";
5: /* Contributed by Wim Vanroose <wim@vanroo.se> */
7: #include <petscksp.h>
8: #include <petscpc.h>
9: #include <petscdm.h>
10: #include <petscdmda.h>
11: #include <petscdmplex.h>
13: static PetscScalar poiss_1D_emat[] = {
14: 1.0000000000000000e+00, -1.0000000000000000e+00,
15: -1.0000000000000000e+00, 1.0000000000000000e+00
16: };
17: static PetscScalar poiss_2D_emat[] = {
18: 6.6666666666666674e-01, -1.6666666666666666e-01, -1.6666666666666666e-01, -3.3333333333333337e-01,
19: -1.6666666666666666e-01, 6.6666666666666674e-01, -3.3333333333333337e-01, -1.6666666666666666e-01,
20: -1.6666666666666666e-01, -3.3333333333333337e-01, 6.6666666666666674e-01, -1.6666666666666666e-01,
21: -3.3333333333333337e-01, -1.6666666666666666e-01, -1.6666666666666666e-01, 6.6666666666666674e-01
22: };
23: static PetscScalar poiss_3D_emat[] = {
24: 3.3333333333333348e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -8.3333333333333343e-02, 0.0000000000000000e+00, -8.3333333333333343e-02, -8.3333333333333343e-02, -8.3333333333333356e-02,
25: 0.0000000000000000e+00, 3.3333333333333337e-01, -8.3333333333333343e-02, 0.0000000000000000e+00, -8.3333333333333343e-02, 0.0000000000000000e+00, -8.3333333333333356e-02, -8.3333333333333343e-02,
26: 0.0000000000000000e+00, -8.3333333333333343e-02, 3.3333333333333337e-01, 0.0000000000000000e+00, -8.3333333333333343e-02, -8.3333333333333356e-02, 0.0000000000000000e+00, -8.3333333333333343e-02,
27: -8.3333333333333343e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, 3.3333333333333348e-01, -8.3333333333333356e-02, -8.3333333333333343e-02, -8.3333333333333343e-02, 0.0000000000000000e+00,
28: 0.0000000000000000e+00, -8.3333333333333343e-02, -8.3333333333333343e-02, -8.3333333333333356e-02, 3.3333333333333337e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -8.3333333333333343e-02,
29: -8.3333333333333343e-02, 0.0000000000000000e+00, -8.3333333333333356e-02, -8.3333333333333343e-02, 0.0000000000000000e+00, 3.3333333333333337e-01, -8.3333333333333343e-02, 0.0000000000000000e+00,
30: -8.3333333333333343e-02, -8.3333333333333356e-02, 0.0000000000000000e+00, -8.3333333333333343e-02, 0.0000000000000000e+00, -8.3333333333333343e-02, 3.3333333333333337e-01, 0.0000000000000000e+00,
31: -8.3333333333333356e-02, -8.3333333333333343e-02, -8.3333333333333343e-02, 0.0000000000000000e+00, -8.3333333333333343e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, 3.3333333333333337e-01
32: };
33: static PetscScalar elast_1D_emat[] = {
34: 3.0000000000000000e+00, -3.0000000000000000e+00,
35: -3.0000000000000000e+00, 3.0000000000000000e+00
36: };
37: static PetscScalar elast_2D_emat[] = {
38: 1.3333333333333335e+00, 5.0000000000000000e-01, -8.3333333333333337e-01, 0.0000000000000000e+00, 1.6666666666666671e-01, 0.0000000000000000e+00, -6.6666666666666674e-01, -5.0000000000000000e-01,
39: 5.0000000000000000e-01, 1.3333333333333335e+00, 0.0000000000000000e+00, 1.6666666666666671e-01, 0.0000000000000000e+00, -8.3333333333333337e-01, -5.0000000000000000e-01, -6.6666666666666674e-01,
40: -8.3333333333333337e-01, 0.0000000000000000e+00, 1.3333333333333335e+00, -5.0000000000000000e-01, -6.6666666666666674e-01, 5.0000000000000000e-01, 1.6666666666666674e-01, 0.0000000000000000e+00,
41: 0.0000000000000000e+00, 1.6666666666666671e-01, -5.0000000000000000e-01, 1.3333333333333335e+00, 5.0000000000000000e-01, -6.6666666666666674e-01, 0.0000000000000000e+00, -8.3333333333333337e-01,
42: 1.6666666666666671e-01, 0.0000000000000000e+00, -6.6666666666666674e-01, 5.0000000000000000e-01, 1.3333333333333335e+00, -5.0000000000000000e-01, -8.3333333333333337e-01, 0.0000000000000000e+00,
43: 0.0000000000000000e+00, -8.3333333333333337e-01, 5.0000000000000000e-01, -6.6666666666666674e-01, -5.0000000000000000e-01, 1.3333333333333335e+00, 0.0000000000000000e+00, 1.6666666666666674e-01,
44: -6.6666666666666674e-01, -5.0000000000000000e-01, 1.6666666666666674e-01, 0.0000000000000000e+00, -8.3333333333333337e-01, 0.0000000000000000e+00, 1.3333333333333335e+00, 5.0000000000000000e-01,
45: -5.0000000000000000e-01, -6.6666666666666674e-01, 0.0000000000000000e+00, -8.3333333333333337e-01, 0.0000000000000000e+00, 1.6666666666666674e-01, 5.0000000000000000e-01, 1.3333333333333335e+00
46: };
47: static PetscScalar elast_3D_emat[] = {
48: 5.5555555555555558e-01, 1.6666666666666666e-01, 1.6666666666666666e-01, -2.2222222222222232e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, 1.1111111111111113e-01, 0.0000000000000000e+00, 8.3333333333333356e-02, -1.9444444444444442e-01, -1.6666666666666669e-01, 0.0000000000000000e+00, 1.1111111111111112e-01, 8.3333333333333356e-02, 0.0000000000000000e+00, -1.9444444444444445e-01, 0.0000000000000000e+00, -1.6666666666666669e-01, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -1.3888888888888887e-01, -8.3333333333333356e-02, -8.3333333333333356e-02,
49: 1.6666666666666666e-01, 5.5555555555555558e-01, 1.6666666666666666e-01, 0.0000000000000000e+00, 1.1111111111111113e-01, 8.3333333333333356e-02, 0.0000000000000000e+00, -2.2222222222222232e-01, 0.0000000000000000e+00, -1.6666666666666669e-01, -1.9444444444444442e-01, 0.0000000000000000e+00, 8.3333333333333356e-02, 1.1111111111111112e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -1.9444444444444445e-01, -1.6666666666666669e-01, -8.3333333333333356e-02, -1.3888888888888887e-01, -8.3333333333333356e-02,
50: 1.6666666666666666e-01, 1.6666666666666666e-01, 5.5555555555555558e-01, 0.0000000000000000e+00, 8.3333333333333356e-02, 1.1111111111111112e-01, 8.3333333333333356e-02, 0.0000000000000000e+00, 1.1111111111111112e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222229e-01, -1.6666666666666669e-01, 0.0000000000000000e+00, -1.9444444444444445e-01, 0.0000000000000000e+00, -1.6666666666666669e-01, -1.9444444444444445e-01, -8.3333333333333356e-02, -8.3333333333333356e-02, -1.3888888888888887e-01,
51: -2.2222222222222232e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, 5.5555555555555558e-01, -1.6666666666666666e-01, -1.6666666666666666e-01, -1.9444444444444442e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, 1.1111111111111113e-01, 0.0000000000000000e+00, -8.3333333333333356e-02, -1.9444444444444445e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, 1.1111111111111113e-01, -8.3333333333333356e-02, 0.0000000000000000e+00, -1.3888888888888887e-01, 8.3333333333333356e-02, 8.3333333333333356e-02, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00,
52: 0.0000000000000000e+00, 1.1111111111111113e-01, 8.3333333333333356e-02, -1.6666666666666666e-01, 5.5555555555555558e-01, 1.6666666666666669e-01, 1.6666666666666669e-01, -1.9444444444444442e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222229e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00, -8.3333333333333356e-02, 1.1111111111111112e-01, 0.0000000000000000e+00, 8.3333333333333356e-02, -1.3888888888888887e-01, -8.3333333333333356e-02, 0.0000000000000000e+00, -1.9444444444444448e-01, -1.6666666666666666e-01,
53: 0.0000000000000000e+00, 8.3333333333333356e-02, 1.1111111111111112e-01, -1.6666666666666666e-01, 1.6666666666666669e-01, 5.5555555555555558e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, -8.3333333333333356e-02, 0.0000000000000000e+00, 1.1111111111111112e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, -1.9444444444444445e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01, 8.3333333333333356e-02, -8.3333333333333356e-02, -1.3888888888888887e-01, 0.0000000000000000e+00, -1.6666666666666666e-01, -1.9444444444444448e-01,
54: 1.1111111111111113e-01, 0.0000000000000000e+00, 8.3333333333333356e-02, -1.9444444444444442e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, 5.5555555555555569e-01, -1.6666666666666666e-01, 1.6666666666666669e-01, -2.2222222222222229e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -1.3888888888888887e-01, 8.3333333333333356e-02, -8.3333333333333356e-02, 1.1111111111111112e-01, -8.3333333333333343e-02, 0.0000000000000000e+00, -1.9444444444444448e-01, 0.0000000000000000e+00, -1.6666666666666669e-01,
55: 0.0000000000000000e+00, -2.2222222222222232e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, -1.9444444444444442e-01, 0.0000000000000000e+00, -1.6666666666666666e-01, 5.5555555555555558e-01, -1.6666666666666669e-01, 0.0000000000000000e+00, 1.1111111111111113e-01, -8.3333333333333343e-02, 0.0000000000000000e+00, -1.9444444444444445e-01, 1.6666666666666669e-01, 8.3333333333333356e-02, -1.3888888888888887e-01, 8.3333333333333356e-02, -8.3333333333333343e-02, 1.1111111111111113e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00,
56: 8.3333333333333356e-02, 0.0000000000000000e+00, 1.1111111111111112e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 1.6666666666666669e-01, -1.6666666666666669e-01, 5.5555555555555558e-01, 0.0000000000000000e+00, -8.3333333333333343e-02, 1.1111111111111112e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, -1.9444444444444445e-01, -8.3333333333333356e-02, 8.3333333333333356e-02, -1.3888888888888887e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01, -1.6666666666666669e-01, 0.0000000000000000e+00, -1.9444444444444448e-01,
57: -1.9444444444444442e-01, -1.6666666666666669e-01, 0.0000000000000000e+00, 1.1111111111111113e-01, 0.0000000000000000e+00, -8.3333333333333356e-02, -2.2222222222222229e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, 5.5555555555555558e-01, 1.6666666666666669e-01, -1.6666666666666666e-01, -1.3888888888888887e-01, -8.3333333333333356e-02, 8.3333333333333356e-02, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -1.9444444444444448e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, 1.1111111111111112e-01, 8.3333333333333343e-02, 0.0000000000000000e+00,
58: -1.6666666666666669e-01, -1.9444444444444442e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222229e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, 1.1111111111111113e-01, -8.3333333333333343e-02, 1.6666666666666669e-01, 5.5555555555555558e-01, -1.6666666666666669e-01, -8.3333333333333356e-02, -1.3888888888888887e-01, 8.3333333333333356e-02, 0.0000000000000000e+00, -1.9444444444444448e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00, 8.3333333333333343e-02, 1.1111111111111112e-01, 0.0000000000000000e+00,
59: 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, -8.3333333333333356e-02, 0.0000000000000000e+00, 1.1111111111111112e-01, 0.0000000000000000e+00, -8.3333333333333343e-02, 1.1111111111111112e-01, -1.6666666666666666e-01, -1.6666666666666669e-01, 5.5555555555555558e-01, 8.3333333333333356e-02, 8.3333333333333356e-02, -1.3888888888888887e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, -1.9444444444444448e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, -1.9444444444444448e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01,
60: 1.1111111111111112e-01, 8.3333333333333356e-02, 0.0000000000000000e+00, -1.9444444444444445e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -1.3888888888888887e-01, -8.3333333333333356e-02, 8.3333333333333356e-02, 5.5555555555555569e-01, 1.6666666666666669e-01, -1.6666666666666669e-01, -2.2222222222222227e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, 1.1111111111111112e-01, 0.0000000000000000e+00, -8.3333333333333343e-02, -1.9444444444444448e-01, -1.6666666666666669e-01, 0.0000000000000000e+00,
61: 8.3333333333333356e-02, 1.1111111111111112e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -1.9444444444444445e-01, 1.6666666666666669e-01, -8.3333333333333356e-02, -1.3888888888888887e-01, 8.3333333333333356e-02, 1.6666666666666669e-01, 5.5555555555555558e-01, -1.6666666666666669e-01, 0.0000000000000000e+00, 1.1111111111111112e-01, -8.3333333333333343e-02, 0.0000000000000000e+00, -2.2222222222222227e-01, 0.0000000000000000e+00, -1.6666666666666669e-01, -1.9444444444444448e-01, 0.0000000000000000e+00,
62: 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222229e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, -1.9444444444444445e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, -1.9444444444444445e-01, 8.3333333333333356e-02, 8.3333333333333356e-02, -1.3888888888888887e-01, -1.6666666666666669e-01, -1.6666666666666669e-01, 5.5555555555555558e-01, 0.0000000000000000e+00, -8.3333333333333343e-02, 1.1111111111111113e-01, -8.3333333333333343e-02, 0.0000000000000000e+00, 1.1111111111111113e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02,
63: -1.9444444444444445e-01, 0.0000000000000000e+00, -1.6666666666666669e-01, 1.1111111111111113e-01, -8.3333333333333356e-02, 0.0000000000000000e+00, -1.3888888888888887e-01, 8.3333333333333356e-02, -8.3333333333333356e-02, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, 5.5555555555555558e-01, -1.6666666666666669e-01, 1.6666666666666669e-01, -1.9444444444444448e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, 1.1111111111111112e-01, 0.0000000000000000e+00, 8.3333333333333343e-02,
64: 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00, -8.3333333333333356e-02, 1.1111111111111112e-01, 0.0000000000000000e+00, 8.3333333333333356e-02, -1.3888888888888887e-01, 8.3333333333333356e-02, 0.0000000000000000e+00, -1.9444444444444448e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, 1.1111111111111112e-01, -8.3333333333333343e-02, -1.6666666666666669e-01, 5.5555555555555558e-01, -1.6666666666666666e-01, 1.6666666666666669e-01, -1.9444444444444448e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01, 0.0000000000000000e+00,
65: -1.6666666666666669e-01, 0.0000000000000000e+00, -1.9444444444444445e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01, -8.3333333333333356e-02, 8.3333333333333356e-02, -1.3888888888888887e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, -1.9444444444444448e-01, 0.0000000000000000e+00, -8.3333333333333343e-02, 1.1111111111111113e-01, 1.6666666666666669e-01, -1.6666666666666666e-01, 5.5555555555555558e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 8.3333333333333343e-02, 0.0000000000000000e+00, 1.1111111111111113e-01,
66: -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -1.3888888888888887e-01, 8.3333333333333356e-02, 8.3333333333333356e-02, 1.1111111111111112e-01, -8.3333333333333343e-02, 0.0000000000000000e+00, -1.9444444444444448e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, 1.1111111111111112e-01, 0.0000000000000000e+00, -8.3333333333333343e-02, -1.9444444444444448e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, 5.5555555555555558e-01, -1.6666666666666669e-01, -1.6666666666666669e-01, -2.2222222222222227e-01, 0.0000000000000000e+00, 0.0000000000000000e+00,
67: 0.0000000000000000e+00, -1.9444444444444445e-01, -1.6666666666666669e-01, 8.3333333333333356e-02, -1.3888888888888887e-01, -8.3333333333333356e-02, -8.3333333333333343e-02, 1.1111111111111113e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01, 0.0000000000000000e+00, 1.6666666666666669e-01, -1.9444444444444448e-01, 0.0000000000000000e+00, -1.6666666666666669e-01, 5.5555555555555558e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, 1.1111111111111112e-01, 8.3333333333333343e-02,
68: 0.0000000000000000e+00, -1.6666666666666669e-01, -1.9444444444444445e-01, 8.3333333333333356e-02, -8.3333333333333356e-02, -1.3888888888888887e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01, 1.6666666666666669e-01, 0.0000000000000000e+00, -1.9444444444444448e-01, -8.3333333333333343e-02, 0.0000000000000000e+00, 1.1111111111111113e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, -1.6666666666666669e-01, 1.6666666666666669e-01, 5.5555555555555558e-01, 0.0000000000000000e+00, 8.3333333333333343e-02, 1.1111111111111113e-01,
69: -1.3888888888888887e-01, -8.3333333333333356e-02, -8.3333333333333356e-02, -2.7777777777777769e-02, 0.0000000000000000e+00, 0.0000000000000000e+00, -1.9444444444444448e-01, 0.0000000000000000e+00, -1.6666666666666669e-01, 1.1111111111111112e-01, 8.3333333333333343e-02, 0.0000000000000000e+00, -1.9444444444444448e-01, -1.6666666666666669e-01, 0.0000000000000000e+00, 1.1111111111111112e-01, 0.0000000000000000e+00, 8.3333333333333343e-02, -2.2222222222222227e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, 5.5555555555555558e-01, 1.6666666666666669e-01, 1.6666666666666669e-01,
70: -8.3333333333333356e-02, -1.3888888888888887e-01, -8.3333333333333356e-02, 0.0000000000000000e+00, -1.9444444444444448e-01, -1.6666666666666666e-01, 0.0000000000000000e+00, -2.7777777777777769e-02, 0.0000000000000000e+00, 8.3333333333333343e-02, 1.1111111111111112e-01, 0.0000000000000000e+00, -1.6666666666666669e-01, -1.9444444444444448e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, 1.1111111111111112e-01, 8.3333333333333343e-02, 1.6666666666666669e-01, 5.5555555555555558e-01, 1.6666666666666669e-01,
71: -8.3333333333333356e-02, -8.3333333333333356e-02, -1.3888888888888887e-01, 0.0000000000000000e+00, -1.6666666666666666e-01, -1.9444444444444448e-01, -1.6666666666666669e-01, 0.0000000000000000e+00, -1.9444444444444448e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.2222222222222227e-01, 0.0000000000000000e+00, 0.0000000000000000e+00, -2.7777777777777769e-02, 8.3333333333333343e-02, 0.0000000000000000e+00, 1.1111111111111113e-01, 0.0000000000000000e+00, 8.3333333333333343e-02, 1.1111111111111113e-01, 1.6666666666666669e-01, 1.6666666666666669e-01, 5.5555555555555558e-01
72: };
74: typedef enum {PDE_POISSON, PDE_ELASTICITY} PDEType;
76: typedef struct {
77: PDEType pde;
78: PetscInt dim;
79: PetscInt dof;
80: PetscInt cells[3];
81: PetscBool useglobal;
82: PetscBool dirbc;
83: PetscBool per[3];
84: PetscBool test;
85: PetscScalar *elemMat;
86: PetscBool use_composite_pc;
87: PetscBool random_initial_guess;
88: PetscBool random_real;
89: } AppCtx;
91: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
92: {
93: const char *pdeTypes[2] = {"Poisson", "Elasticity"};
94: PetscInt n,pde;
98: options->pde = PDE_POISSON;
99: options->elemMat = NULL;
100: options->dim = 1;
101: options->cells[0] = 8;
102: options->cells[1] = 6;
103: options->cells[2] = 4;
104: options->useglobal = PETSC_FALSE;
105: options->dirbc = PETSC_TRUE;
106: options->test = PETSC_FALSE;
107: options->per[0] = PETSC_FALSE;
108: options->per[1] = PETSC_FALSE;
109: options->per[2] = PETSC_FALSE;
110: options->use_composite_pc = PETSC_FALSE;
111: options->random_initial_guess = PETSC_FALSE;
112: options->random_real = PETSC_FALSE;
114: PetscOptionsBegin(comm,NULL,"Problem Options",NULL);
115: pde = options->pde;
116: PetscOptionsEList("-pde_type","The PDE type",__FILE__,pdeTypes,2,pdeTypes[options->pde],&pde,NULL);
117: options->pde = (PDEType)pde;
118: PetscOptionsInt("-dim","The topological mesh dimension",__FILE__,options->dim,&options->dim,NULL);
119: PetscOptionsIntArray("-cells","The mesh division",__FILE__,options->cells,(n=3,&n),NULL);
120: PetscOptionsBoolArray("-periodicity","The mesh periodicity",__FILE__,options->per,(n=3,&n),NULL);
121: PetscOptionsBool("-use_global","Test MatSetValues",__FILE__,options->useglobal,&options->useglobal,NULL);
122: PetscOptionsBool("-dirichlet","Use dirichlet BC",__FILE__,options->dirbc,&options->dirbc,NULL);
123: PetscOptionsBool("-test_assembly","Test MATIS assembly",__FILE__,options->test,&options->test,NULL);
124: PetscOptionsBool("-use_composite_pc","Multiplicative composite with BDDC + Richardson/Jacobi",__FILE__,options->use_composite_pc,&options->use_composite_pc,NULL);
125: PetscOptionsBool("-random_initial_guess","Solve A x = 0 with random initial guess, instead of A x = b with random b",__FILE__,options->random_initial_guess,&options->random_initial_guess,NULL);
126: PetscOptionsBool("-random_real","Use real-valued b (or x, if -random_initial_guess) instead of default scalar type",__FILE__,options->random_real,&options->random_real,NULL);
127: PetscOptionsEnd();
129: for (n=options->dim;n<3;n++) options->cells[n] = 0;
130: if (options->per[0]) options->dirbc = PETSC_FALSE;
132: /* element matrices */
133: switch (options->pde) {
134: case PDE_ELASTICITY:
135: options->dof = options->dim;
136: switch (options->dim) {
137: case 1:
138: options->elemMat = elast_1D_emat;
139: break;
140: case 2:
141: options->elemMat = elast_2D_emat;
142: break;
143: case 3:
144: options->elemMat = elast_3D_emat;
145: break;
146: default: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Unsupported dimension %D",options->dim);
147: }
148: break;
149: case PDE_POISSON:
150: options->dof = 1;
151: switch (options->dim) {
152: case 1:
153: options->elemMat = poiss_1D_emat;
154: break;
155: case 2:
156: options->elemMat = poiss_2D_emat;
157: break;
158: case 3:
159: options->elemMat = poiss_3D_emat;
160: break;
161: default: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Unsupported dimension %D",options->dim);
162: }
163: break;
164: default: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Unsupported PDE %D",options->pde);
165: }
166: return 0;
167: }
169: int main(int argc,char **args)
170: {
171: AppCtx user;
172: KSP ksp;
173: PC pc;
174: Mat A;
175: DM da;
176: Vec x,b,xcoor,xcoorl;
177: IS zero;
178: ISLocalToGlobalMapping map;
179: MatNullSpace nullsp = NULL;
180: PetscInt i;
181: PetscInt nel,nen; /* Number of elements & element nodes */
182: const PetscInt *e_loc; /* Local indices of element nodes (in local element order) */
183: PetscInt *e_glo = NULL; /* Global indices of element nodes (in local element order) */
184: PetscInt nodes[3];
185: PetscBool ismatis;
186: #if defined(PETSC_USE_LOG)
187: PetscLogStage stages[2];
188: #endif
189: PetscErrorCode ierr;
191: PetscInitialize(&argc,&args,(char*)0,help);
192: ProcessOptions(PETSC_COMM_WORLD,&user);
193: for (i=0; i<3; i++) nodes[i] = user.cells[i] + !user.per[i];
194: switch (user.dim) {
195: case 3:
196: DMDACreate3d(PETSC_COMM_WORLD,user.per[0] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE,
197: user.per[1] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE,
198: user.per[2] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE,
199: DMDA_STENCIL_BOX,nodes[0],nodes[1],nodes[2],
200: PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,user.dof,
201: 1,PETSC_NULL,PETSC_NULL,PETSC_NULL,&da);
202: break;
203: case 2:
204: DMDACreate2d(PETSC_COMM_WORLD,user.per[0] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE,
205: user.per[1] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE,
206: DMDA_STENCIL_BOX,nodes[0],nodes[1],
207: PETSC_DECIDE,PETSC_DECIDE,user.dof,
208: 1,PETSC_NULL,PETSC_NULL,&da);
209: break;
210: case 1:
211: DMDACreate1d(PETSC_COMM_WORLD,user.per[0] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE,
212: nodes[0],user.dof,1,PETSC_NULL,&da);
213: break;
214: default: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Unsupported dimension %D",user.dim);
215: }
217: PetscLogStageRegister("KSPSetUp",&stages[0]);
218: PetscLogStageRegister("KSPSolve",&stages[1]);
220: DMSetMatType(da,MATIS);
221: DMSetFromOptions(da);
222: DMDASetElementType(da,DMDA_ELEMENT_Q1);
223: DMSetUp(da);
224: {
225: PetscInt M,N,P;
226: DMDAGetInfo(da,0,&M,&N,&P,0,0,0,0,0,0,0,0,0);
227: switch (user.dim) {
228: case 3:
229: user.cells[2] = P - !user.per[2];
230: case 2:
231: user.cells[1] = N - !user.per[1];
232: case 1:
233: user.cells[0] = M - !user.per[0];
234: break;
235: default: SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Unsupported dimension %D",user.dim);
236: }
237: }
238: DMDASetUniformCoordinates(da,0.0,1.0*user.cells[0],0.0,1.0*user.cells[1],0.0,1.0*user.cells[2]);
239: DMGetCoordinates(da,&xcoor);
241: DMCreateMatrix(da,&A);
242: MatSetFromOptions(A);
243: DMGetLocalToGlobalMapping(da,&map);
244: DMDAGetElements(da,&nel,&nen,&e_loc);
245: if (user.useglobal) {
246: PetscMalloc1(nel*nen,&e_glo);
247: ISLocalToGlobalMappingApplyBlock(map,nen*nel,e_loc,e_glo);
248: }
250: /* we reorder the indices since the element matrices are given in lexicographic order,
251: whereas the elements indices returned by DMDAGetElements follow the usual FEM ordering
252: i.e., element matrices DMDA ordering
253: 2---3 3---2
254: / / / /
255: 0---1 0---1
256: */
257: for (i = 0; i < nel; ++i) {
258: PetscInt ord[8] = {0,1,3,2,4,5,7,6};
259: PetscInt j,idxs[8];
262: if (!e_glo) {
263: for (j=0;j<nen;j++) idxs[j] = e_loc[i*nen+ord[j]];
264: MatSetValuesBlockedLocal(A,nen,idxs,nen,idxs,user.elemMat,ADD_VALUES);
265: } else {
266: for (j=0;j<nen;j++) idxs[j] = e_glo[i*nen+ord[j]];
267: MatSetValuesBlocked(A,nen,idxs,nen,idxs,user.elemMat,ADD_VALUES);
268: }
269: }
270: DMDARestoreElements(da,&nel,&nen,&e_loc);
271: MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);
272: MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);
273: MatSetOption(A,MAT_SPD,PETSC_TRUE);
275: /* Boundary conditions */
276: zero = NULL;
277: if (user.dirbc) { /* fix one side of DMDA */
278: Vec nat,glob;
279: PetscScalar *vals;
280: PetscInt n,*idx,j,st;
282: n = PetscGlobalRank ? 0 : (user.cells[1]+1)*(user.cells[2]+1);
283: ISCreateStride(PETSC_COMM_WORLD,n,0,user.cells[0]+1,&zero);
284: if (user.dof > 1) { /* zero all components */
285: const PetscInt *idx;
286: IS bzero;
288: ISGetIndices(zero,(const PetscInt**)&idx);
289: ISCreateBlock(PETSC_COMM_WORLD,user.dof,n,idx,PETSC_COPY_VALUES,&bzero);
290: ISRestoreIndices(zero,(const PetscInt**)&idx);
291: ISDestroy(&zero);
292: zero = bzero;
293: }
294: /* map indices from natural to global */
295: DMDACreateNaturalVector(da,&nat);
296: ISGetLocalSize(zero,&n);
297: PetscMalloc1(n,&vals);
298: for (i=0;i<n;i++) vals[i] = 1.0;
299: ISGetIndices(zero,(const PetscInt**)&idx);
300: VecSetValues(nat,n,idx,vals,INSERT_VALUES);
301: ISRestoreIndices(zero,(const PetscInt**)&idx);
302: PetscFree(vals);
303: VecAssemblyBegin(nat);
304: VecAssemblyEnd(nat);
305: DMCreateGlobalVector(da,&glob);
306: DMDANaturalToGlobalBegin(da,nat,INSERT_VALUES,glob);
307: DMDANaturalToGlobalEnd(da,nat,INSERT_VALUES,glob);
308: VecDestroy(&nat);
309: ISDestroy(&zero);
310: VecGetLocalSize(glob,&n);
311: PetscMalloc1(n,&idx);
312: VecGetOwnershipRange(glob,&st,NULL);
313: VecGetArray(glob,&vals);
314: for (i=0,j=0;i<n;i++) if (PetscRealPart(vals[i]) == 1.0) idx[j++] = i + st;
315: VecRestoreArray(glob,&vals);
316: VecDestroy(&glob);
317: ISCreateGeneral(PETSC_COMM_WORLD,j,idx,PETSC_OWN_POINTER,&zero);
318: MatZeroRowsColumnsIS(A,zero,1.0,NULL,NULL);
319: } else {
320: switch (user.pde) {
321: case PDE_POISSON:
322: MatNullSpaceCreate(PETSC_COMM_WORLD,PETSC_TRUE,0,NULL,&nullsp);
323: break;
324: case PDE_ELASTICITY:
325: MatNullSpaceCreateRigidBody(xcoor,&nullsp);
326: break;
327: }
328: /* with periodic BC and Elasticity, just the displacements are in the nullspace
329: this is no harm since we eliminate all the components of the rhs */
330: MatSetNullSpace(A,nullsp);
331: }
333: if (user.test) {
334: Mat AA;
336: MatConvert(A,MATAIJ,MAT_INITIAL_MATRIX,&AA);
337: MatViewFromOptions(AA,NULL,"-assembled_view");
338: MatDestroy(&AA);
339: }
341: /* Attach near null space for elasticity */
342: if (user.pde == PDE_ELASTICITY) {
343: MatNullSpace nearnullsp;
345: MatNullSpaceCreateRigidBody(xcoor,&nearnullsp);
346: MatSetNearNullSpace(A,nearnullsp);
347: MatNullSpaceDestroy(&nearnullsp);
348: }
350: /* we may want to use MG for the local solvers: attach local nearnullspace to the local matrices */
351: DMGetCoordinatesLocal(da,&xcoorl);
352: PetscObjectTypeCompare((PetscObject)A,MATIS,&ismatis);
353: if (ismatis) {
354: MatNullSpace lnullsp = NULL;
355: Mat lA;
357: MatISGetLocalMat(A,&lA);
358: if (user.pde == PDE_ELASTICITY) {
359: Vec lc;
360: ISLocalToGlobalMapping l2l;
361: IS is;
362: const PetscScalar *a;
363: const PetscInt *idxs;
364: PetscInt n,bs;
366: /* when using a DMDA, the local matrices have an additional local-to-local map
367: that maps from the DA local ordering to the ordering induced by the elements */
368: MatCreateVecs(lA,&lc,NULL);
369: MatGetLocalToGlobalMapping(lA,&l2l,NULL);
370: VecSetLocalToGlobalMapping(lc,l2l);
371: VecSetOption(lc,VEC_IGNORE_NEGATIVE_INDICES,PETSC_TRUE);
372: VecGetLocalSize(xcoorl,&n);
373: VecGetBlockSize(xcoorl,&bs);
374: ISCreateStride(PETSC_COMM_SELF,n/bs,0,1,&is);
375: ISGetIndices(is,&idxs);
376: VecGetArrayRead(xcoorl,&a);
377: VecSetValuesBlockedLocal(lc,n/bs,idxs,a,INSERT_VALUES);
378: VecAssemblyBegin(lc);
379: VecAssemblyEnd(lc);
380: VecRestoreArrayRead(xcoorl,&a);
381: ISRestoreIndices(is,&idxs);
382: ISDestroy(&is);
383: MatNullSpaceCreateRigidBody(lc,&lnullsp);
384: VecDestroy(&lc);
385: } else if (user.pde == PDE_POISSON) {
386: MatNullSpaceCreate(PETSC_COMM_SELF,PETSC_TRUE,0,NULL,&lnullsp);
387: }
388: MatSetNearNullSpace(lA,lnullsp);
389: MatNullSpaceDestroy(&lnullsp);
390: MatISRestoreLocalMat(A,&lA);
391: }
393: KSPCreate(PETSC_COMM_WORLD,&ksp);
394: KSPSetOperators(ksp,A,A);
395: KSPSetType(ksp,KSPCG);
396: KSPGetPC(ksp,&pc);
397: if (user.use_composite_pc) {
398: PC pcksp,pcjacobi;
399: KSP ksprich;
400: PCSetType(pc,PCCOMPOSITE);
401: PCCompositeSetType(pc,PC_COMPOSITE_MULTIPLICATIVE);
402: PCCompositeAddPCType(pc,PCBDDC);
403: PCCompositeAddPCType(pc,PCKSP);
404: PCCompositeGetPC(pc,1,&pcksp);
405: PCKSPGetKSP(pcksp,&ksprich);
406: KSPSetType(ksprich,KSPRICHARDSON);
407: KSPSetTolerances(ksprich,PETSC_DEFAULT,PETSC_DEFAULT,PETSC_DEFAULT,1);
408: KSPSetNormType(ksprich,KSP_NORM_NONE);
409: KSPSetConvergenceTest(ksprich,KSPConvergedSkip,NULL,NULL);
410: KSPGetPC(ksprich,&pcjacobi);
411: PCSetType(pcjacobi,PCJACOBI);
412: } else {
413: PCSetType(pc,PCBDDC);
414: }
415: /* PCBDDCSetDirichletBoundaries(pc,zero); */
416: KSPSetFromOptions(ksp);
417: PetscLogStagePush(stages[0]);
418: KSPSetUp(ksp);
419: PetscLogStagePop();
421: MatCreateVecs(A,&x,&b);
422: if (user.random_initial_guess) {
423: /* Solving A x = 0 with random initial guess allows Arnoldi to run for more iterations, thereby yielding a more
424: * complete Hessenberg matrix and more accurate eigenvalues. */
425: VecZeroEntries(b);
426: VecSetRandom(x,NULL);
427: if (user.random_real) VecRealPart(x);
428: if (nullsp) {
429: MatNullSpaceRemove(nullsp,x);
430: }
431: KSPSetInitialGuessNonzero(ksp,PETSC_TRUE);
432: KSPSetComputeEigenvalues(ksp,PETSC_TRUE);
433: KSPGMRESSetRestart(ksp,100);
434: } else {
435: VecSetRandom(b,NULL);
436: if (user.random_real) VecRealPart(x);
437: if (nullsp) {
438: MatNullSpaceRemove(nullsp,b);
439: }
440: }
441: PetscLogStagePush(stages[1]);
442: KSPSolve(ksp,b,x);
443: PetscLogStagePop();
445: /* cleanup */
446: VecDestroy(&x);
447: VecDestroy(&b);
448: ISDestroy(&zero);
449: PetscFree(e_glo);
450: MatNullSpaceDestroy(&nullsp);
451: KSPDestroy(&ksp);
452: MatDestroy(&A);
453: DMDestroy(&da);
454: PetscFinalize();
455: return 0;
456: }
458: /*TEST
460: test:
461: nsize: 8
462: filter: grep -v "variant HERMITIAN"
463: suffix: bddc_1
464: args: -pde_type Poisson -dim 3 -dirichlet 0 -ksp_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged
465: test:
466: nsize: 8
467: filter: grep -v "variant HERMITIAN"
468: suffix: bddc_2
469: args: -pde_type Poisson -dim 3 -dirichlet 0 -ksp_view -use_global -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged
470: test:
471: nsize: 8
472: filter: grep -v "variant HERMITIAN"
473: suffix: bddc_elast
474: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic
475: test:
476: nsize: 8
477: filter: grep -v "variant HERMITIAN"
478: suffix: bddc_elast_3lev
479: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_levels 1 -pc_bddc_coarsening_ratio 1 -ksp_error_if_not_converged -pc_bddc_monolithic -pc_bddc_use_faces -pc_bddc_coarse_pc_bddc_corner_selection
480: testset:
481: nsize: 8
482: requires: hpddm slepc defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
483: # on some architectures, this test will converge in 19 or 21 iterations
484: filter: grep -v "variant HERMITIAN" | grep -v " tolerance" | sed -e "s/CONVERGED_RTOL iterations [1-2][91]\{0,1\}$/CONVERGED_RTOL iterations 20/g"
485: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_levels 1 -pc_bddc_coarsening_ratio 1 -ksp_error_if_not_converged -pc_bddc_monolithic -pc_bddc_use_faces -pc_bddc_coarse_pc_type hpddm -prefix_push pc_bddc_coarse_ -pc_hpddm_levels_1_sub_pc_type cholesky -pc_hpddm_levels_1_eps_nev 6 -pc_hpddm_levels_1_st_pc_factor_shift_type INBLOCKS -prefix_pop -ksp_type fgmres -ksp_max_it 50 -ksp_converged_reason
486: test:
487: args: -pc_bddc_coarse_pc_hpddm_coarse_mat_type baij -options_left no
488: suffix: bddc_elast_3lev_hpddm_baij
489: test:
490: requires: !complex
491: suffix: bddc_elast_3lev_hpddm
492: test:
493: nsize: 8
494: requires: !single
495: filter: grep -v "variant HERMITIAN"
496: suffix: bddc_elast_4lev
497: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_levels 2 -pc_bddc_coarsening_ratio 2 -ksp_error_if_not_converged -pc_bddc_monolithic -pc_bddc_use_faces -pc_bddc_coarse_pc_bddc_corner_selection -pc_bddc_coarse_l1_pc_bddc_corner_selection -mat_partitioning_type average -options_left 0
498: test:
499: nsize: 8
500: filter: grep -v "variant HERMITIAN"
501: suffix: bddc_elast_deluxe_layers
502: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic -pc_bddc_use_deluxe_scaling -pc_bddc_schur_layers 1
503: test:
504: nsize: 8
505: filter: grep -v "variant HERMITIAN" | sed -e "s/iterations 1[0-9]/iterations 10/g"
506: suffix: bddc_elast_dir_approx
507: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic -pc_bddc_dirichlet_pc_type gamg -pc_bddc_dirichlet_pc_gamg_esteig_ksp_max_it 10 -ksp_converged_reason -pc_bddc_dirichlet_approximate
508: test:
509: nsize: 8
510: filter: grep -v "variant HERMITIAN" | sed -e "s/iterations 1[0-9]/iterations 10/g"
511: suffix: bddc_elast_neu_approx
512: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic -pc_bddc_neumann_pc_type gamg -pc_bddc_neumann_pc_gamg_esteig_ksp_max_it 10 -ksp_converged_reason -pc_bddc_neumann_approximate
513: test:
514: nsize: 8
515: filter: grep -v "variant HERMITIAN" | sed -e "s/iterations 1[0-9]/iterations 10/g"
516: suffix: bddc_elast_both_approx
517: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic -pc_bddc_dirichlet_pc_type gamg -pc_bddc_dirichlet_pc_gamg_esteig_ksp_max_it 10 -pc_bddc_neumann_pc_type gamg -pc_bddc_neumann_pc_gamg_esteig_ksp_max_it 10 -ksp_converged_reason -pc_bddc_neumann_approximate -pc_bddc_dirichlet_approximate
518: test:
519: nsize: 8
520: filter: grep -v "variant HERMITIAN"
521: suffix: fetidp_1
522: args: -pde_type Poisson -dim 3 -dirichlet 0 -ksp_view -ksp_type fetidp -fetidp_ksp_type cg -fetidp_bddc_pc_bddc_coarse_redundant_pc_type svd -ksp_fetidp_fullyredundant -ksp_error_if_not_converged
523: test:
524: nsize: 8
525: filter: grep -v "variant HERMITIAN"
526: suffix: fetidp_2
527: args: -pde_type Poisson -dim 3 -dirichlet 0 -ksp_view -use_global -ksp_type fetidp -fetidp_ksp_type cg -fetidp_bddc_pc_bddc_coarse_redundant_pc_type svd -ksp_fetidp_fullyredundant -ksp_error_if_not_converged
528: test:
529: nsize: 8
530: filter: grep -v "variant HERMITIAN"
531: suffix: fetidp_elast
532: args: -pde_type Elasticity -cells 9,7,8 -dim 3 -ksp_view -ksp_type fetidp -fetidp_ksp_type cg -fetidp_bddc_pc_bddc_coarse_redundant_pc_type svd -ksp_fetidp_fullyredundant -ksp_error_if_not_converged -fetidp_bddc_pc_bddc_monolithic
533: testset:
534: nsize: 8
535: requires: hpddm slepc defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
536: args: -pde_type Elasticity -cells 12,12 -dim 2 -ksp_converged_reason -pc_type hpddm -pc_hpddm_coarse_correction balanced -pc_hpddm_levels_1_pc_type asm -pc_hpddm_levels_1_pc_asm_overlap 1 -pc_hpddm_levels_1_pc_asm_type basic -pc_hpddm_levels_1_sub_pc_type cholesky -pc_hpddm_levels_1_eps_nev 10 -matis_localmat_type {{aij baij sbaij}shared output} -pc_hpddm_coarse_mat_type {{baij sbaij}shared output} -pc_hpddm_levels_1_st_pc_factor_shift_type INBLOCKS
537: test:
538: suffix: hpddm
539: output_file: output/ex71_hpddm.out
540: test:
541: args: -pc_hpddm_levels_1_eps_type lapack -pc_hpddm_levels_1_eps_smallest_magnitude -pc_hpddm_levels_1_st_type shift
542: suffix: hpddm_lapack
543: output_file: output/ex71_hpddm.out
544: testset:
545: nsize: 9
546: args: -test_assembly -assembled_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged
547: test:
548: args: -dim 1 -cells 12 -pde_type Poisson
549: suffix: dmda_matis_poiss_1d_loc
550: output_file: output/ex71_dmda_matis_poiss_1d.out
551: test:
552: args: -dim 1 -cells 12 -pde_type Poisson -use_global
553: suffix: dmda_matis_poiss_1d_glob
554: output_file: output/ex71_dmda_matis_poiss_1d.out
555: test:
556: args: -dim 1 -cells 12 -pde_type Elasticity
557: suffix: dmda_matis_elast_1d_loc
558: output_file: output/ex71_dmda_matis_elast_1d.out
559: test:
560: args: -dim 1 -cells 12 -pde_type Elasticity -use_global
561: suffix: dmda_matis_elast_1d_glob
562: output_file: output/ex71_dmda_matis_elast_1d.out
563: test:
564: args: -dim 2 -cells 5,7 -pde_type Poisson
565: suffix: dmda_matis_poiss_2d_loc
566: output_file: output/ex71_dmda_matis_poiss_2d.out
567: test:
568: args: -dim 2 -cells 5,7 -pde_type Poisson -use_global
569: suffix: dmda_matis_poiss_2d_glob
570: output_file: output/ex71_dmda_matis_poiss_2d.out
571: test:
572: args: -dim 2 -cells 5,7 -pde_type Elasticity
573: suffix: dmda_matis_elast_2d_loc
574: output_file: output/ex71_dmda_matis_elast_2d.out
575: test:
576: args: -dim 2 -cells 5,7 -pde_type Elasticity -use_global
577: suffix: dmda_matis_elast_2d_glob
578: output_file: output/ex71_dmda_matis_elast_2d.out
579: test:
580: args: -dim 3 -cells 3,3,3 -pde_type Poisson
581: suffix: dmda_matis_poiss_3d_loc
582: output_file: output/ex71_dmda_matis_poiss_3d.out
583: test:
584: args: -dim 3 -cells 3,3,3 -pde_type Poisson -use_global
585: suffix: dmda_matis_poiss_3d_glob
586: output_file: output/ex71_dmda_matis_poiss_3d.out
587: test:
588: args: -dim 3 -cells 3,3,3 -pde_type Elasticity
589: suffix: dmda_matis_elast_3d_loc
590: output_file: output/ex71_dmda_matis_elast_3d.out
591: test:
592: args: -dim 3 -cells 3,3,3 -pde_type Elasticity -use_global
593: suffix: dmda_matis_elast_3d_glob
594: output_file: output/ex71_dmda_matis_elast_3d.out
595: test:
596: nsize: 8
597: filter: grep -v "variant HERMITIAN"
598: suffix: bddc_elast_deluxe_layers_adapt
599: requires: mumps !complex
600: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_converged_reason -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic -sub_schurs_mat_solver_type mumps -pc_bddc_use_deluxe_scaling -pc_bddc_adaptive_threshold 2.0 -pc_bddc_schur_layers {{1 10}separate_output} -pc_bddc_adaptive_userdefined {{0 1}separate output} -sub_schurs_schur_mat_type seqdense
601: # gitlab runners have a quite old MKL (2016) which interacts badly with AMD machines (not Intel-based ones!)
602: # this is the reason behind the filtering rule
603: test:
604: nsize: 8
605: suffix: bddc_elast_deluxe_layers_adapt_mkl_pardiso
606: filter: sed -e "s/CONVERGED_RTOL iterations [1-2][0-9]/CONVERGED_RTOL iterations 13/g" | sed -e "s/CONVERGED_RTOL iterations 6/CONVERGED_RTOL iterations 5/g"
607: requires: mkl_pardiso !complex
608: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_converged_reason -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic -sub_schurs_mat_solver_type mkl_pardiso -sub_schurs_mat_mkl_pardiso_65 1 -pc_bddc_use_deluxe_scaling -pc_bddc_adaptive_threshold 2.0 -pc_bddc_schur_layers {{1 10}separate_output} -pc_bddc_adaptive_userdefined {{0 1}separate output} -sub_schurs_schur_mat_type seqdense
609: test:
610: nsize: 8
611: filter: grep -v "variant HERMITIAN"
612: suffix: bddc_cusparse
613: # no kokkos since it seems kokkos's resource demand is too much with 8 ranks and the test will fail on cuda related initialization.
614: requires: cuda !kokkos
615: args: -pde_type Poisson -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_dirichlet_pc_type cholesky -pc_bddc_dirichlet_pc_factor_mat_solver_type cusparse -pc_bddc_dirichlet_pc_factor_mat_ordering_type nd -pc_bddc_neumann_pc_type cholesky -pc_bddc_neumann_pc_factor_mat_solver_type cusparse -pc_bddc_neumann_pc_factor_mat_ordering_type nd -matis_localmat_type aijcusparse
616: test:
617: nsize: 8
618: filter: grep -v "variant HERMITIAN"
619: suffix: bddc_elast_deluxe_layers_adapt_cuda
620: requires: !complex mumps cuda defined(PETSC_HAVE_CUSOLVERDNDPOTRI)
621: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_converged_reason -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic -sub_schurs_mat_solver_type mumps -pc_bddc_use_deluxe_scaling -pc_bddc_adaptive_threshold 2.0 -pc_bddc_schur_layers {{1 10}separate_output} -pc_bddc_adaptive_userdefined {{0 1}separate output} -matis_localmat_type seqaijcusparse -sub_schurs_schur_mat_type {{seqdensecuda seqdense}}
622: test:
623: nsize: 8
624: filter: grep -v "variant HERMITIAN" | grep -v "I-node routines" | sed -e "s/seqaijcusparse/seqaij/g"
625: suffix: bddc_elast_deluxe_layers_adapt_cuda_approx
626: requires: !complex mumps cuda defined(PETSC_HAVE_CUSOLVERDNDPOTRI)
627: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic -sub_schurs_mat_solver_type mumps -pc_bddc_use_deluxe_scaling -pc_bddc_adaptive_threshold 2.0 -pc_bddc_schur_layers 1 -matis_localmat_type {{seqaij seqaijcusparse}separate output} -sub_schurs_schur_mat_type {{seqdensecuda seqdense}} -pc_bddc_dirichlet_pc_type gamg -pc_bddc_dirichlet_approximate -pc_bddc_neumann_pc_type gamg -pc_bddc_neumann_approximate -pc_bddc_dirichlet_pc_gamg_esteig_ksp_max_it 10 -pc_bddc_neumann_pc_gamg_esteig_ksp_max_it 10
628: test:
629: nsize: 8
630: suffix: bddc_elast_deluxe_layers_adapt_mkl_pardiso_cuda
631: requires: !complex mkl_pardiso cuda defined(PETSC_HAVE_CUSOLVERDNDPOTRI)
632: filter: sed -e "s/CONVERGED_RTOL iterations 6/CONVERGED_RTOL iterations 5/g"
633: args: -pde_type Elasticity -cells 7,9,8 -dim 3 -ksp_converged_reason -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged -pc_bddc_monolithic -sub_schurs_mat_solver_type mkl_pardiso -sub_schurs_mat_mkl_pardiso_65 1 -pc_bddc_use_deluxe_scaling -pc_bddc_adaptive_threshold 2.0 -pc_bddc_schur_layers {{1 10}separate_output} -pc_bddc_adaptive_userdefined {{0 1}separate output} -matis_localmat_type seqaijcusparse -sub_schurs_schur_mat_type {{seqdensecuda seqdense}}
635: testset:
636: nsize: 2
637: output_file: output/ex71_aij_dmda_preall.out
638: filter: sed -e "s/CONVERGED_RTOL iterations 7/CONVERGED_RTOL iterations 6/g"
639: args: -pde_type Poisson -dim 1 -cells 6 -pc_type none -ksp_converged_reason
640: test:
641: suffix: aijviennacl_dmda_preall
642: requires: viennacl
643: args: -dm_mat_type aijviennacl -dm_preallocate_only {{0 1}} -dirichlet {{0 1}}
644: # -dm_preallocate_only 0 is broken
645: test:
646: suffix: aijcusparse_dmda_preall
647: requires: cuda
648: args: -dm_mat_type aijcusparse -dm_preallocate_only -dirichlet {{0 1}}
649: test:
650: suffix: aij_dmda_preall
651: args: -dm_mat_type aij -dm_preallocate_only {{0 1}} -dirichlet {{0 1}}
652: testset:
653: nsize: 4
654: args: -dim 2 -cells 16,16 -periodicity 1,1 -random_initial_guess -random_real -sub_0_pc_bddc_switch_static -use_composite_pc -ksp_monitor -ksp_converged_reason -ksp_type gmres -ksp_view_singularvalues -ksp_view_eigenvalues -sub_0_pc_bddc_use_edges 0 -sub_0_pc_bddc_coarse_pc_type svd -sub_1_ksp_ksp_max_it 1 -sub_1_ksp_ksp_richardson_scale 2.3
655: test:
656: args: -sub_0_pc_bddc_interface_ext_type lump
657: suffix: composite_bddc_lumped
658: test:
659: requires: !single
660: args: -sub_0_pc_bddc_interface_ext_type dirichlet
661: suffix: composite_bddc_dirichlet
663: TEST*/