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[] = {1.0000000000000000e+00, -1.0000000000000000e+00, -1.0000000000000000e+00, 1.0000000000000000e+00};
 14: static PetscScalar poiss_2D_emat[] = {6.6666666666666674e-01,  -1.6666666666666666e-01, -1.6666666666666666e-01, -3.3333333333333337e-01, -1.6666666666666666e-01, 6.6666666666666674e-01,  -3.3333333333333337e-01, -1.6666666666666666e-01,
 15:                                       -1.6666666666666666e-01, -3.3333333333333337e-01, 6.6666666666666674e-01,  -1.6666666666666666e-01, -3.3333333333333337e-01, -1.6666666666666666e-01, -1.6666666666666666e-01, 6.6666666666666674e-01};
 16: static PetscScalar poiss_3D_emat[] = {3.3333333333333348e-01,  0.0000000000000000e+00,  0.0000000000000000e+00,  -8.3333333333333343e-02, 0.0000000000000000e+00,  -8.3333333333333343e-02, -8.3333333333333343e-02, -8.3333333333333356e-02,
 17:                                       0.0000000000000000e+00,  3.3333333333333337e-01,  -8.3333333333333343e-02, 0.0000000000000000e+00,  -8.3333333333333343e-02, 0.0000000000000000e+00,  -8.3333333333333356e-02, -8.3333333333333343e-02,
 18:                                       0.0000000000000000e+00,  -8.3333333333333343e-02, 3.3333333333333337e-01,  0.0000000000000000e+00,  -8.3333333333333343e-02, -8.3333333333333356e-02, 0.0000000000000000e+00,  -8.3333333333333343e-02,
 19:                                       -8.3333333333333343e-02, 0.0000000000000000e+00,  0.0000000000000000e+00,  3.3333333333333348e-01,  -8.3333333333333356e-02, -8.3333333333333343e-02, -8.3333333333333343e-02, 0.0000000000000000e+00,
 20:                                       0.0000000000000000e+00,  -8.3333333333333343e-02, -8.3333333333333343e-02, -8.3333333333333356e-02, 3.3333333333333337e-01,  0.0000000000000000e+00,  0.0000000000000000e+00,  -8.3333333333333343e-02,
 21:                                       -8.3333333333333343e-02, 0.0000000000000000e+00,  -8.3333333333333356e-02, -8.3333333333333343e-02, 0.0000000000000000e+00,  3.3333333333333337e-01,  -8.3333333333333343e-02, 0.0000000000000000e+00,
 22:                                       -8.3333333333333343e-02, -8.3333333333333356e-02, 0.0000000000000000e+00,  -8.3333333333333343e-02, 0.0000000000000000e+00,  -8.3333333333333343e-02, 3.3333333333333337e-01,  0.0000000000000000e+00,
 23:                                       -8.3333333333333356e-02, -8.3333333333333343e-02, -8.3333333333333343e-02, 0.0000000000000000e+00,  -8.3333333333333343e-02, 0.0000000000000000e+00,  0.0000000000000000e+00,  3.3333333333333337e-01};
 24: static PetscScalar elast_1D_emat[] = {3.0000000000000000e+00, -3.0000000000000000e+00, -3.0000000000000000e+00, 3.0000000000000000e+00};
 25: static PetscScalar elast_2D_emat[] = {1.3333333333333335e+00,  5.0000000000000000e-01,  -8.3333333333333337e-01, 0.0000000000000000e+00,  1.6666666666666671e-01,  0.0000000000000000e+00,  -6.6666666666666674e-01, -5.0000000000000000e-01,
 26:                                       5.0000000000000000e-01,  1.3333333333333335e+00,  0.0000000000000000e+00,  1.6666666666666671e-01,  0.0000000000000000e+00,  -8.3333333333333337e-01, -5.0000000000000000e-01, -6.6666666666666674e-01,
 27:                                       -8.3333333333333337e-01, 0.0000000000000000e+00,  1.3333333333333335e+00,  -5.0000000000000000e-01, -6.6666666666666674e-01, 5.0000000000000000e-01,  1.6666666666666674e-01,  0.0000000000000000e+00,
 28:                                       0.0000000000000000e+00,  1.6666666666666671e-01,  -5.0000000000000000e-01, 1.3333333333333335e+00,  5.0000000000000000e-01,  -6.6666666666666674e-01, 0.0000000000000000e+00,  -8.3333333333333337e-01,
 29:                                       1.6666666666666671e-01,  0.0000000000000000e+00,  -6.6666666666666674e-01, 5.0000000000000000e-01,  1.3333333333333335e+00,  -5.0000000000000000e-01, -8.3333333333333337e-01, 0.0000000000000000e+00,
 30:                                       0.0000000000000000e+00,  -8.3333333333333337e-01, 5.0000000000000000e-01,  -6.6666666666666674e-01, -5.0000000000000000e-01, 1.3333333333333335e+00,  0.0000000000000000e+00,  1.6666666666666674e-01,
 31:                                       -6.6666666666666674e-01, -5.0000000000000000e-01, 1.6666666666666674e-01,  0.0000000000000000e+00,  -8.3333333333333337e-01, 0.0000000000000000e+00,  1.3333333333333335e+00,  5.0000000000000000e-01,
 32:                                       -5.0000000000000000e-01, -6.6666666666666674e-01, 0.0000000000000000e+00,  -8.3333333333333337e-01, 0.0000000000000000e+00,  1.6666666666666674e-01,  5.0000000000000000e-01,  1.3333333333333335e+00};
 33: static PetscScalar elast_3D_emat[] =
 34:   {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,
 35:    -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,
 36:    -2.7777777777777769e-02, 0.0000000000000000e+00,  0.0000000000000000e+00,  -1.3888888888888887e-01, -8.3333333333333356e-02, -8.3333333333333356e-02, 1.6666666666666666e-01,  5.5555555555555558e-01,  1.6666666666666666e-01,
 37:    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,
 38:    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,
 39:    -8.3333333333333356e-02, -1.3888888888888887e-01, -8.3333333333333356e-02, 1.6666666666666666e-01,  1.6666666666666666e-01,  5.5555555555555558e-01,  0.0000000000000000e+00,  8.3333333333333356e-02,  1.1111111111111112e-01,
 40:    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,
 41:    -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,
 42:    -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,
 43:    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,
 44:    -1.3888888888888887e-01, 8.3333333333333356e-02,  8.3333333333333356e-02,  -2.7777777777777769e-02, 0.0000000000000000e+00,  0.0000000000000000e+00,  0.0000000000000000e+00,  1.1111111111111113e-01,  8.3333333333333356e-02,
 45:    -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,
 46:    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,
 47:    0.0000000000000000e+00,  -1.9444444444444448e-01, -1.6666666666666666e-01, 0.0000000000000000e+00,  8.3333333333333356e-02,  1.1111111111111112e-01,  -1.6666666666666666e-01, 1.6666666666666669e-01,  5.5555555555555558e-01,
 48:    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,
 49:    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,
 50:    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,
 51:    -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,
 52:    1.1111111111111112e-01,  -8.3333333333333343e-02, 0.0000000000000000e+00,  -1.9444444444444448e-01, 0.0000000000000000e+00,  -1.6666666666666669e-01, 0.0000000000000000e+00,  -2.2222222222222232e-01, 0.0000000000000000e+00,
 53:    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,
 54:    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,
 55:    0.0000000000000000e+00,  -2.7777777777777769e-02, 0.0000000000000000e+00,  8.3333333333333356e-02,  0.0000000000000000e+00,  1.1111111111111112e-01,  0.0000000000000000e+00,  0.0000000000000000e+00,  -2.7777777777777769e-02,
 56:    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,
 57:    -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,
 58:    -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,
 59:    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,
 60:    -1.9444444444444448e-01, 0.0000000000000000e+00,  1.6666666666666669e-01,  1.1111111111111112e-01,  8.3333333333333343e-02,  0.0000000000000000e+00,  -1.6666666666666669e-01, -1.9444444444444442e-01, 0.0000000000000000e+00,
 61:    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,
 62:    -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,
 63:    8.3333333333333343e-02,  1.1111111111111112e-01,  0.0000000000000000e+00,  0.0000000000000000e+00,  0.0000000000000000e+00,  -2.7777777777777769e-02, -8.3333333333333356e-02, 0.0000000000000000e+00,  1.1111111111111112e-01,
 64:    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,
 65:    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,
 66:    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,
 67:    -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,
 68:    1.1111111111111112e-01,  0.0000000000000000e+00,  -8.3333333333333343e-02, -1.9444444444444448e-01, -1.6666666666666669e-01, 0.0000000000000000e+00,  8.3333333333333356e-02,  1.1111111111111112e-01,  0.0000000000000000e+00,
 69:    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,
 70:    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,
 71:    -1.6666666666666669e-01, -1.9444444444444448e-01, 0.0000000000000000e+00,  0.0000000000000000e+00,  0.0000000000000000e+00,  -2.2222222222222229e-01, 1.6666666666666669e-01,  0.0000000000000000e+00,  -1.9444444444444445e-01,
 72:    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,
 73:    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,
 74:    -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,
 75:    -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,
 76:    -1.9444444444444448e-01, 1.6666666666666669e-01,  0.0000000000000000e+00,  1.1111111111111112e-01,  0.0000000000000000e+00,  8.3333333333333343e-02,  0.0000000000000000e+00,  -2.7777777777777769e-02, 0.0000000000000000e+00,
 77:    -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,
 78:    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,
 79:    0.0000000000000000e+00,  -2.2222222222222227e-01, 0.0000000000000000e+00,  -1.6666666666666669e-01, 0.0000000000000000e+00,  -1.9444444444444445e-01, 0.0000000000000000e+00,  0.0000000000000000e+00,  -2.2222222222222227e-01,
 80:    -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,
 81:    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,
 82:    -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,
 83:    -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,
 84:    5.5555555555555558e-01,  -1.6666666666666669e-01, -1.6666666666666669e-01, -2.2222222222222227e-01, 0.0000000000000000e+00,  0.0000000000000000e+00,  0.0000000000000000e+00,  -1.9444444444444445e-01, -1.6666666666666669e-01,
 85:    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,
 86:    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,
 87:    0.0000000000000000e+00,  1.1111111111111112e-01,  8.3333333333333343e-02,  0.0000000000000000e+00,  -1.6666666666666669e-01, -1.9444444444444445e-01, 8.3333333333333356e-02,  -8.3333333333333356e-02, -1.3888888888888887e-01,
 88:    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,
 89:    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,
 90:    -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,
 91:    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,
 92:    -2.2222222222222227e-01, 0.0000000000000000e+00,  0.0000000000000000e+00,  5.5555555555555558e-01,  1.6666666666666669e-01,  1.6666666666666669e-01,  -8.3333333333333356e-02, -1.3888888888888887e-01, -8.3333333333333356e-02,
 93:    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,
 94:    -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,
 95:    1.6666666666666669e-01,  5.5555555555555558e-01,  1.6666666666666669e-01,  -8.3333333333333356e-02, -8.3333333333333356e-02, -1.3888888888888887e-01, 0.0000000000000000e+00,  -1.6666666666666666e-01, -1.9444444444444448e-01,
 96:    -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,
 97:    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};

 99: typedef enum {
100:   PDE_POISSON,
101:   PDE_ELASTICITY
102: } PDEType;

104: typedef struct {
105:   PDEType      pde;
106:   PetscInt     dim;
107:   PetscInt     dof;
108:   PetscInt     cells[3];
109:   PetscBool    useglobal;
110:   PetscBool    dirbc;
111:   PetscBool    per[3];
112:   PetscBool    test;
113:   PetscScalar *elemMat;
114:   PetscBool    use_composite_pc;
115:   PetscBool    random_initial_guess;
116:   PetscBool    random_real;
117: } AppCtx;

119: static PetscErrorCode ProcessOptions(MPI_Comm comm, AppCtx *options)
120: {
121:   const char *pdeTypes[2] = {"Poisson", "Elasticity"};
122:   PetscInt    n, pde;

124:   PetscFunctionBeginUser;
125:   options->pde                  = PDE_POISSON;
126:   options->elemMat              = NULL;
127:   options->dim                  = 1;
128:   options->cells[0]             = 8;
129:   options->cells[1]             = 6;
130:   options->cells[2]             = 4;
131:   options->useglobal            = PETSC_FALSE;
132:   options->dirbc                = PETSC_TRUE;
133:   options->test                 = PETSC_FALSE;
134:   options->per[0]               = PETSC_FALSE;
135:   options->per[1]               = PETSC_FALSE;
136:   options->per[2]               = PETSC_FALSE;
137:   options->use_composite_pc     = PETSC_FALSE;
138:   options->random_initial_guess = PETSC_FALSE;
139:   options->random_real          = PETSC_FALSE;

141:   PetscOptionsBegin(comm, NULL, "Problem Options", NULL);
142:   pde = options->pde;
143:   PetscCall(PetscOptionsEList("-pde_type", "The PDE type", __FILE__, pdeTypes, 2, pdeTypes[options->pde], &pde, NULL));
144:   options->pde = (PDEType)pde;
145:   PetscCall(PetscOptionsInt("-dim", "The topological mesh dimension", __FILE__, options->dim, &options->dim, NULL));
146:   PetscCall(PetscOptionsIntArray("-cells", "The mesh division", __FILE__, options->cells, (n = 3, &n), NULL));
147:   PetscCall(PetscOptionsBoolArray("-periodicity", "The mesh periodicity", __FILE__, options->per, (n = 3, &n), NULL));
148:   PetscCall(PetscOptionsBool("-use_global", "Test MatSetValues", __FILE__, options->useglobal, &options->useglobal, NULL));
149:   PetscCall(PetscOptionsBool("-dirichlet", "Use dirichlet BC", __FILE__, options->dirbc, &options->dirbc, NULL));
150:   PetscCall(PetscOptionsBool("-test_assembly", "Test MATIS assembly", __FILE__, options->test, &options->test, NULL));
151:   PetscCall(PetscOptionsBool("-use_composite_pc", "Multiplicative composite with BDDC + Richardson/Jacobi", __FILE__, options->use_composite_pc, &options->use_composite_pc, NULL));
152:   PetscCall(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));
153:   PetscCall(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));
154:   PetscOptionsEnd();

156:   for (n = options->dim; n < 3; n++) options->cells[n] = 0;
157:   if (options->per[0]) options->dirbc = PETSC_FALSE;

159:   /* element matrices */
160:   switch (options->pde) {
161:   case PDE_ELASTICITY:
162:     options->dof = options->dim;
163:     switch (options->dim) {
164:     case 1:
165:       options->elemMat = elast_1D_emat;
166:       break;
167:     case 2:
168:       options->elemMat = elast_2D_emat;
169:       break;
170:     case 3:
171:       options->elemMat = elast_3D_emat;
172:       break;
173:     default:
174:       SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Unsupported dimension %" PetscInt_FMT, options->dim);
175:     }
176:     break;
177:   case PDE_POISSON:
178:     options->dof = 1;
179:     switch (options->dim) {
180:     case 1:
181:       options->elemMat = poiss_1D_emat;
182:       break;
183:     case 2:
184:       options->elemMat = poiss_2D_emat;
185:       break;
186:     case 3:
187:       options->elemMat = poiss_3D_emat;
188:       break;
189:     default:
190:       SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Unsupported dimension %" PetscInt_FMT, options->dim);
191:     }
192:     break;
193:   default:
194:     SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Unsupported PDE %d", options->pde);
195:   }
196:   PetscFunctionReturn(PETSC_SUCCESS);
197: }

199: int main(int argc, char **args)
200: {
201:   AppCtx                 user;
202:   KSP                    ksp;
203:   PC                     pc;
204:   Mat                    A;
205:   DM                     da;
206:   Vec                    x, b, xcoor, xcoorl;
207:   IS                     zero;
208:   ISLocalToGlobalMapping map;
209:   MatNullSpace           nullsp = NULL;
210:   PetscInt               i;
211:   PetscInt               nel, nen;     /* Number of elements & element nodes */
212:   const PetscInt        *e_loc;        /* Local indices of element nodes (in local element order) */
213:   PetscInt              *e_glo = NULL; /* Global indices of element nodes (in local element order) */
214:   PetscInt               nodes[3];
215:   PetscBool              ismatis;
216:   PetscLogStage          stages[2];

218:   PetscFunctionBeginUser;
219:   PetscCall(PetscInitialize(&argc, &args, (char *)0, help));
220:   PetscCall(ProcessOptions(PETSC_COMM_WORLD, &user));
221:   for (i = 0; i < 3; i++) nodes[i] = user.cells[i] + !user.per[i];
222:   switch (user.dim) {
223:   case 3:
224:     PetscCall(DMDACreate3d(PETSC_COMM_WORLD, user.per[0] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE, user.per[1] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE, user.per[2] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, nodes[0], nodes[1], nodes[2], PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE,
225:                            user.dof, 1, NULL, NULL, NULL, &da));
226:     break;
227:   case 2:
228:     PetscCall(DMDACreate2d(PETSC_COMM_WORLD, user.per[0] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE, user.per[1] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE, DMDA_STENCIL_BOX, nodes[0], nodes[1], PETSC_DECIDE, PETSC_DECIDE, user.dof, 1, NULL, NULL, &da));
229:     break;
230:   case 1:
231:     PetscCall(DMDACreate1d(PETSC_COMM_WORLD, user.per[0] ? DM_BOUNDARY_PERIODIC : DM_BOUNDARY_NONE, nodes[0], user.dof, 1, NULL, &da));
232:     break;
233:   default:
234:     SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Unsupported dimension %" PetscInt_FMT, user.dim);
235:   }

237:   PetscCall(PetscLogStageRegister("KSPSetUp", &stages[0]));
238:   PetscCall(PetscLogStageRegister("KSPSolve", &stages[1]));

240:   PetscCall(DMSetMatType(da, MATIS));
241:   PetscCall(DMSetFromOptions(da));
242:   PetscCall(DMDASetElementType(da, DMDA_ELEMENT_Q1));
243:   PetscCall(DMSetUp(da));
244:   {
245:     PetscInt M, N, P;
246:     PetscCall(DMDAGetInfo(da, 0, &M, &N, &P, 0, 0, 0, 0, 0, 0, 0, 0, 0));
247:     switch (user.dim) {
248:     case 3:
249:       user.cells[2] = P - !user.per[2]; /* fall through */
250:     case 2:
251:       user.cells[1] = N - !user.per[1]; /* fall through */
252:     case 1:
253:       user.cells[0] = M - !user.per[0];
254:       break;
255:     default:
256:       SETERRQ(PETSC_COMM_WORLD, PETSC_ERR_SUP, "Unsupported dimension %" PetscInt_FMT, user.dim);
257:     }
258:   }
259:   PetscCall(DMDASetUniformCoordinates(da, 0.0, 1.0 * user.cells[0], 0.0, 1.0 * user.cells[1], 0.0, 1.0 * user.cells[2]));
260:   PetscCall(DMGetCoordinates(da, &xcoor));

262:   PetscCall(DMCreateMatrix(da, &A));
263:   PetscCall(MatSetFromOptions(A));
264:   PetscCall(DMGetLocalToGlobalMapping(da, &map));
265:   PetscCall(DMDAGetElements(da, &nel, &nen, &e_loc));
266:   if (user.useglobal) {
267:     PetscCall(PetscMalloc1(nel * nen, &e_glo));
268:     PetscCall(ISLocalToGlobalMappingApplyBlock(map, nen * nel, e_loc, e_glo));
269:   }

271:   /* we reorder the indices since the element matrices are given in lexicographic order,
272:      whereas the elements indices returned by DMDAGetElements follow the usual FEM ordering
273:      i.e., element matrices     DMDA ordering
274:                2---3              3---2
275:               /   /              /   /
276:              0---1              0---1
277:   */
278:   for (i = 0; i < nel; ++i) {
279:     PetscInt ord[8] = {0, 1, 3, 2, 4, 5, 7, 6};
280:     PetscInt j, idxs[8];

282:     PetscCheck(nen <= 8, PETSC_COMM_WORLD, PETSC_ERR_SUP, "Not coded");
283:     if (!e_glo) {
284:       for (j = 0; j < nen; j++) idxs[j] = e_loc[i * nen + ord[j]];
285:       PetscCall(MatSetValuesBlockedLocal(A, nen, idxs, nen, idxs, user.elemMat, ADD_VALUES));
286:     } else {
287:       for (j = 0; j < nen; j++) idxs[j] = e_glo[i * nen + ord[j]];
288:       PetscCall(MatSetValuesBlocked(A, nen, idxs, nen, idxs, user.elemMat, ADD_VALUES));
289:     }
290:   }
291:   PetscCall(DMDARestoreElements(da, &nel, &nen, &e_loc));
292:   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
293:   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
294:   PetscCall(MatSetOption(A, MAT_SPD, PETSC_TRUE));
295:   PetscCall(MatSetOption(A, MAT_SPD_ETERNAL, PETSC_TRUE));

297:   /* Boundary conditions */
298:   zero = NULL;
299:   if (user.dirbc) { /* fix one side of DMDA */
300:     Vec          nat, glob;
301:     PetscScalar *vals;
302:     PetscInt     n, *idx, j, st;

304:     n = PetscGlobalRank ? 0 : (user.cells[1] + 1) * (user.cells[2] + 1);
305:     PetscCall(ISCreateStride(PETSC_COMM_WORLD, n, 0, user.cells[0] + 1, &zero));
306:     if (user.dof > 1) { /* zero all components */
307:       const PetscInt *idx;
308:       IS              bzero;

310:       PetscCall(ISGetIndices(zero, (const PetscInt **)&idx));
311:       PetscCall(ISCreateBlock(PETSC_COMM_WORLD, user.dof, n, idx, PETSC_COPY_VALUES, &bzero));
312:       PetscCall(ISRestoreIndices(zero, (const PetscInt **)&idx));
313:       PetscCall(ISDestroy(&zero));
314:       zero = bzero;
315:     }
316:     /* map indices from natural to global */
317:     PetscCall(DMDACreateNaturalVector(da, &nat));
318:     PetscCall(ISGetLocalSize(zero, &n));
319:     PetscCall(PetscMalloc1(n, &vals));
320:     for (i = 0; i < n; i++) vals[i] = 1.0;
321:     PetscCall(ISGetIndices(zero, (const PetscInt **)&idx));
322:     PetscCall(VecSetValues(nat, n, idx, vals, INSERT_VALUES));
323:     PetscCall(ISRestoreIndices(zero, (const PetscInt **)&idx));
324:     PetscCall(PetscFree(vals));
325:     PetscCall(VecAssemblyBegin(nat));
326:     PetscCall(VecAssemblyEnd(nat));
327:     PetscCall(DMCreateGlobalVector(da, &glob));
328:     PetscCall(DMDANaturalToGlobalBegin(da, nat, INSERT_VALUES, glob));
329:     PetscCall(DMDANaturalToGlobalEnd(da, nat, INSERT_VALUES, glob));
330:     PetscCall(VecDestroy(&nat));
331:     PetscCall(ISDestroy(&zero));
332:     PetscCall(VecGetLocalSize(glob, &n));
333:     PetscCall(PetscMalloc1(n, &idx));
334:     PetscCall(VecGetOwnershipRange(glob, &st, NULL));
335:     PetscCall(VecGetArray(glob, &vals));
336:     for (i = 0, j = 0; i < n; i++)
337:       if (PetscRealPart(vals[i]) == 1.0) idx[j++] = i + st;
338:     PetscCall(VecRestoreArray(glob, &vals));
339:     PetscCall(VecDestroy(&glob));
340:     PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, j, idx, PETSC_OWN_POINTER, &zero));
341:     PetscCall(MatZeroRowsColumnsIS(A, zero, 1.0, NULL, NULL));
342:   } else {
343:     switch (user.pde) {
344:     case PDE_POISSON:
345:       PetscCall(MatNullSpaceCreate(PETSC_COMM_WORLD, PETSC_TRUE, 0, NULL, &nullsp));
346:       break;
347:     case PDE_ELASTICITY:
348:       PetscCall(MatNullSpaceCreateRigidBody(xcoor, &nullsp));
349:       break;
350:     }
351:     /* with periodic BC and Elasticity, just the displacements are in the nullspace
352:        this is no harm since we eliminate all the components of the rhs */
353:     PetscCall(MatSetNullSpace(A, nullsp));
354:   }

356:   if (user.test) {
357:     Mat AA;

359:     PetscCall(MatConvert(A, MATAIJ, MAT_INITIAL_MATRIX, &AA));
360:     PetscCall(MatViewFromOptions(AA, NULL, "-assembled_view"));
361:     PetscCall(MatDestroy(&AA));
362:   }

364:   /* Attach near null space for elasticity */
365:   if (user.pde == PDE_ELASTICITY) {
366:     MatNullSpace nearnullsp;

368:     PetscCall(MatNullSpaceCreateRigidBody(xcoor, &nearnullsp));
369:     PetscCall(MatSetNearNullSpace(A, nearnullsp));
370:     PetscCall(MatNullSpaceDestroy(&nearnullsp));
371:   }

373:   /* we may want to use MG for the local solvers: attach local nearnullspace to the local matrices */
374:   PetscCall(DMGetCoordinatesLocal(da, &xcoorl));
375:   PetscCall(PetscObjectTypeCompare((PetscObject)A, MATIS, &ismatis));
376:   if (ismatis) {
377:     MatNullSpace lnullsp = NULL;
378:     Mat          lA;

380:     PetscCall(MatISGetLocalMat(A, &lA));
381:     if (user.pde == PDE_ELASTICITY) {
382:       Vec                    lc;
383:       ISLocalToGlobalMapping l2l;
384:       IS                     is;
385:       const PetscScalar     *a;
386:       const PetscInt        *idxs;
387:       PetscInt               n, bs;

389:       /* when using a DMDA, the local matrices have an additional local-to-local map
390:          that maps from the DA local ordering to the ordering induced by the elements */
391:       PetscCall(MatCreateVecs(lA, &lc, NULL));
392:       PetscCall(MatGetLocalToGlobalMapping(lA, &l2l, NULL));
393:       PetscCall(VecSetLocalToGlobalMapping(lc, l2l));
394:       PetscCall(VecSetOption(lc, VEC_IGNORE_NEGATIVE_INDICES, PETSC_TRUE));
395:       PetscCall(VecGetLocalSize(xcoorl, &n));
396:       PetscCall(VecGetBlockSize(xcoorl, &bs));
397:       PetscCall(ISCreateStride(PETSC_COMM_SELF, n / bs, 0, 1, &is));
398:       PetscCall(ISGetIndices(is, &idxs));
399:       PetscCall(VecGetArrayRead(xcoorl, &a));
400:       PetscCall(VecSetValuesBlockedLocal(lc, n / bs, idxs, a, INSERT_VALUES));
401:       PetscCall(VecAssemblyBegin(lc));
402:       PetscCall(VecAssemblyEnd(lc));
403:       PetscCall(VecRestoreArrayRead(xcoorl, &a));
404:       PetscCall(ISRestoreIndices(is, &idxs));
405:       PetscCall(ISDestroy(&is));
406:       PetscCall(MatNullSpaceCreateRigidBody(lc, &lnullsp));
407:       PetscCall(VecDestroy(&lc));
408:     } else if (user.pde == PDE_POISSON) {
409:       PetscCall(MatNullSpaceCreate(PETSC_COMM_SELF, PETSC_TRUE, 0, NULL, &lnullsp));
410:     }
411:     PetscCall(MatSetNearNullSpace(lA, lnullsp));
412:     PetscCall(MatNullSpaceDestroy(&lnullsp));
413:     PetscCall(MatISRestoreLocalMat(A, &lA));
414:   }

416:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
417:   PetscCall(KSPSetOperators(ksp, A, A));
418:   PetscCall(KSPSetType(ksp, KSPCG));
419:   PetscCall(KSPGetPC(ksp, &pc));
420:   if (user.use_composite_pc) {
421:     PC  pcksp, pcjacobi;
422:     KSP ksprich;
423:     PetscCall(PCSetType(pc, PCCOMPOSITE));
424:     PetscCall(PCCompositeSetType(pc, PC_COMPOSITE_MULTIPLICATIVE));
425:     PetscCall(PCCompositeAddPCType(pc, PCBDDC));
426:     PetscCall(PCCompositeAddPCType(pc, PCKSP));
427:     PetscCall(PCCompositeGetPC(pc, 1, &pcksp));
428:     PetscCall(PCKSPGetKSP(pcksp, &ksprich));
429:     PetscCall(KSPSetType(ksprich, KSPRICHARDSON));
430:     PetscCall(KSPSetTolerances(ksprich, PETSC_DEFAULT, PETSC_DEFAULT, PETSC_DEFAULT, 1));
431:     PetscCall(KSPSetNormType(ksprich, KSP_NORM_NONE));
432:     PetscCall(KSPSetConvergenceTest(ksprich, KSPConvergedSkip, NULL, NULL));
433:     PetscCall(KSPGetPC(ksprich, &pcjacobi));
434:     PetscCall(PCSetType(pcjacobi, PCJACOBI));
435:   } else {
436:     PetscCall(PCSetType(pc, PCBDDC));
437:   }
438:   /* PetscCall(PCBDDCSetDirichletBoundaries(pc,zero)); */
439:   PetscCall(KSPSetFromOptions(ksp));
440:   PetscCall(PetscLogStagePush(stages[0]));
441:   PetscCall(KSPSetUp(ksp));
442:   PetscCall(PetscLogStagePop());

444:   PetscCall(MatCreateVecs(A, &x, &b));
445:   if (user.random_initial_guess) {
446:     /* Solving A x = 0 with random initial guess allows Arnoldi to run for more iterations, thereby yielding a more
447:      * complete Hessenberg matrix and more accurate eigenvalues. */
448:     PetscCall(VecZeroEntries(b));
449:     PetscCall(VecSetRandom(x, NULL));
450:     if (user.random_real) PetscCall(VecRealPart(x));
451:     if (nullsp) PetscCall(MatNullSpaceRemove(nullsp, x));
452:     PetscCall(KSPSetInitialGuessNonzero(ksp, PETSC_TRUE));
453:     PetscCall(KSPSetComputeEigenvalues(ksp, PETSC_TRUE));
454:     PetscCall(KSPGMRESSetRestart(ksp, 100));
455:   } else {
456:     PetscCall(VecSetRandom(b, NULL));
457:     if (user.random_real) PetscCall(VecRealPart(x));
458:     if (nullsp) PetscCall(MatNullSpaceRemove(nullsp, b));
459:   }
460:   PetscCall(PetscLogStagePush(stages[1]));
461:   PetscCall(KSPSolve(ksp, b, x));
462:   PetscCall(PetscLogStagePop());

464:   /* cleanup */
465:   PetscCall(VecDestroy(&x));
466:   PetscCall(VecDestroy(&b));
467:   PetscCall(ISDestroy(&zero));
468:   PetscCall(PetscFree(e_glo));
469:   PetscCall(MatNullSpaceDestroy(&nullsp));
470:   PetscCall(KSPDestroy(&ksp));
471:   PetscCall(MatDestroy(&A));
472:   PetscCall(DMDestroy(&da));
473:   PetscCall(PetscFinalize());
474:   return 0;
475: }

477: /*TEST

479:  test:
480:    nsize: 8
481:    filter: grep -v "variant HERMITIAN"
482:    suffix: bddc_1
483:    args: -pde_type Poisson -dim 3 -dirichlet 0 -ksp_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged
484:  test:
485:    nsize: 8
486:    filter: grep -v "variant HERMITIAN"
487:    suffix: bddc_2
488:    args: -pde_type Poisson -dim 3 -dirichlet 0 -ksp_view -use_global -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged
489:  test:
490:    nsize: 8
491:    filter: grep -v "variant HERMITIAN"
492:    suffix: bddc_elast
493:    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
494:  test:
495:    nsize: 8
496:    filter: grep -v "variant HERMITIAN"
497:    suffix: bddc_elast_3lev
498:    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
499:  testset:
500:    nsize: 8
501:    requires: hpddm slepc defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
502:    # on some architectures, this test will converge in 19 or 21 iterations
503:    filter: grep -v "variant HERMITIAN" | grep -v " tolerance"  | sed -e "s/CONVERGED_RTOL iterations [1-2][91]\{0,1\}$/CONVERGED_RTOL iterations 20/g"
504:    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
505:    test:
506:      args: -pc_bddc_coarse_pc_hpddm_coarse_mat_type baij -options_left no
507:      suffix: bddc_elast_3lev_hpddm_baij
508:    test:
509:      requires: !complex
510:      suffix: bddc_elast_3lev_hpddm
511:  test:
512:    nsize: 8
513:    requires: !single
514:    filter: grep -v "variant HERMITIAN"
515:    suffix: bddc_elast_4lev
516:    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
517:  test:
518:    nsize: 8
519:    filter: grep -v "variant HERMITIAN"
520:    suffix: bddc_elast_deluxe_layers
521:    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
522:  test:
523:    nsize: 8
524:    filter: grep -v "variant HERMITIAN" | sed -e "s/iterations 1[0-9]/iterations 10/g"
525:    suffix: bddc_elast_dir_approx
526:    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
527:  test:
528:    nsize: 8
529:    filter: grep -v "variant HERMITIAN" | sed -e "s/iterations 1[0-9]/iterations 10/g"
530:    suffix: bddc_elast_neu_approx
531:    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
532:  test:
533:    nsize: 8
534:    filter: grep -v "variant HERMITIAN" | sed -e "s/iterations 1[0-9]/iterations 10/g"
535:    suffix: bddc_elast_both_approx
536:    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
537:  test:
538:    nsize: 8
539:    filter: grep -v "variant HERMITIAN"
540:    suffix: fetidp_1
541:    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
542:  test:
543:    nsize: 8
544:    filter: grep -v "variant HERMITIAN"
545:    suffix: fetidp_2
546:    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
547:  test:
548:    nsize: 8
549:    filter: grep -v "variant HERMITIAN"
550:    suffix: fetidp_elast
551:    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
552:  testset:
553:    nsize: 8
554:    requires: hpddm slepc defined(PETSC_HAVE_DYNAMIC_LIBRARIES) defined(PETSC_USE_SHARED_LIBRARIES)
555:    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 -pc_hpddm_levels_1_st_pc_factor_shift_type INBLOCKS
556:    test:
557:      args: -matis_localmat_type {{aij baij sbaij}shared output} -pc_hpddm_coarse_mat_type {{baij sbaij}shared output}
558:      suffix: hpddm
559:      output_file: output/ex71_hpddm.out
560:    test:
561:      args: -matis_localmat_type sbaij -pc_hpddm_coarse_mat_type sbaij -pc_hpddm_levels_1_st_share_sub_ksp -pc_hpddm_levels_1_eps_type lapack -pc_hpddm_levels_1_eps_smallest_magnitude -pc_hpddm_levels_1_st_type shift
562:      suffix: hpddm_lapack
563:      output_file: output/ex71_hpddm.out
564:  testset:
565:    nsize: 9
566:    args: -test_assembly -assembled_view -pc_bddc_coarse_redundant_pc_type svd -ksp_error_if_not_converged
567:    test:
568:      args: -dim 1 -cells 12 -pde_type Poisson
569:      suffix: dmda_matis_poiss_1d_loc
570:      output_file: output/ex71_dmda_matis_poiss_1d.out
571:    test:
572:      args: -dim 1 -cells 12 -pde_type Poisson -use_global
573:      suffix: dmda_matis_poiss_1d_glob
574:      output_file: output/ex71_dmda_matis_poiss_1d.out
575:    test:
576:      args: -dim 1 -cells 12 -pde_type Elasticity
577:      suffix: dmda_matis_elast_1d_loc
578:      output_file: output/ex71_dmda_matis_elast_1d.out
579:    test:
580:      args: -dim 1 -cells 12 -pde_type Elasticity -use_global
581:      suffix: dmda_matis_elast_1d_glob
582:      output_file: output/ex71_dmda_matis_elast_1d.out
583:    test:
584:      args: -dim 2 -cells 5,7 -pde_type Poisson
585:      suffix: dmda_matis_poiss_2d_loc
586:      output_file: output/ex71_dmda_matis_poiss_2d.out
587:    test:
588:      args: -dim 2 -cells 5,7 -pde_type Poisson -use_global
589:      suffix: dmda_matis_poiss_2d_glob
590:      output_file: output/ex71_dmda_matis_poiss_2d.out
591:    test:
592:      args: -dim 2 -cells 5,7 -pde_type Elasticity
593:      suffix: dmda_matis_elast_2d_loc
594:      output_file: output/ex71_dmda_matis_elast_2d.out
595:    test:
596:      args: -dim 2 -cells 5,7 -pde_type Elasticity -use_global
597:      suffix: dmda_matis_elast_2d_glob
598:      output_file: output/ex71_dmda_matis_elast_2d.out
599:    test:
600:      args: -dim 3 -cells 3,3,3 -pde_type Poisson
601:      suffix: dmda_matis_poiss_3d_loc
602:      output_file: output/ex71_dmda_matis_poiss_3d.out
603:    test:
604:      args: -dim 3 -cells 3,3,3 -pde_type Poisson -use_global
605:      suffix: dmda_matis_poiss_3d_glob
606:      output_file: output/ex71_dmda_matis_poiss_3d.out
607:    test:
608:      args: -dim 3 -cells 3,3,3 -pde_type Elasticity
609:      suffix: dmda_matis_elast_3d_loc
610:      output_file: output/ex71_dmda_matis_elast_3d.out
611:    test:
612:      args: -dim 3 -cells 3,3,3 -pde_type Elasticity -use_global
613:      suffix: dmda_matis_elast_3d_glob
614:      output_file: output/ex71_dmda_matis_elast_3d.out
615:  test:
616:    nsize: 8
617:    filter: grep -v "variant HERMITIAN"
618:    suffix: bddc_elast_deluxe_layers_adapt
619:    requires: mumps !complex
620:    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
621:  # gitlab runners have a quite old MKL (2016) which interacts badly with AMD machines (not Intel-based ones!)
622:  # this is the reason behind the filtering rule
623:  test:
624:    nsize: 8
625:    suffix: bddc_elast_deluxe_layers_adapt_mkl_pardiso
626:    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"
627:    requires: mkl_pardiso !complex
628:    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
629:  test:
630:    nsize: 8
631:    filter: grep -v "variant HERMITIAN"
632:    suffix: bddc_cusparse
633:    # no kokkos since it seems kokkos's resource demand is too much with 8 ranks and the test will fail on cuda related initialization.
634:    requires: cuda !kokkos
635:    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
636:  test:
637:    nsize: 8
638:    filter: grep -v "variant HERMITIAN"
639:    suffix: bddc_elast_deluxe_layers_adapt_cuda
640:    requires: !complex mumps cuda defined(PETSC_HAVE_CUSOLVERDNDPOTRI)
641:    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}}
642:  test:
643:    nsize: 8
644:    filter: grep -v "variant HERMITIAN" | grep -v "I-node routines" | sed -e "s/seqaijcusparse/seqaij/g"
645:    suffix: bddc_elast_deluxe_layers_adapt_cuda_approx
646:    requires: !complex mumps cuda defined(PETSC_HAVE_CUSOLVERDNDPOTRI)
647:    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
648:  test:
649:    nsize: 8
650:    suffix: bddc_elast_deluxe_layers_adapt_mkl_pardiso_cuda
651:    requires: !complex mkl_pardiso cuda defined(PETSC_HAVE_CUSOLVERDNDPOTRI)
652:    filter: sed -e "s/CONVERGED_RTOL iterations 6/CONVERGED_RTOL iterations 5/g"
653:    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}}

655:  testset:
656:    nsize: 2
657:    output_file: output/ex71_aij_dmda_preall.out
658:    filter: sed -e "s/CONVERGED_RTOL iterations 7/CONVERGED_RTOL iterations 6/g"
659:    args: -pde_type Poisson -dim 1 -cells 6 -pc_type none -ksp_converged_reason
660:    test:
661:      suffix: aijviennacl_dmda_preall
662:      requires: viennacl
663:      args: -dm_mat_type aijviennacl -dm_preallocate_only {{0 1}} -dirichlet {{0 1}}
664:    # -dm_preallocate_only 0 is broken
665:    test:
666:      suffix: aijcusparse_dmda_preall
667:      requires: cuda
668:      args: -dm_mat_type aijcusparse -dm_preallocate_only -dirichlet {{0 1}}
669:    test:
670:      suffix: aij_dmda_preall
671:      args: -dm_mat_type aij -dm_preallocate_only {{0 1}} -dirichlet {{0 1}}
672:  testset:
673:    nsize: 4
674:    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
675:    test:
676:      args: -sub_0_pc_bddc_interface_ext_type lump
677:      suffix: composite_bddc_lumped
678:    test:
679:      requires: !single
680:      args: -sub_0_pc_bddc_interface_ext_type dirichlet
681:      suffix: composite_bddc_dirichlet

683:  testset:
684:    nsize: 8
685:    filter: grep -v "variant HERMITIAN"
686:    args: -cells 7,9,8 -dim 3 -ksp_view -ksp_error_if_not_converged -pc_type mg -pc_mg_levels 2 -pc_mg_galerkin -pc_mg_adapt_interp_coarse_space gdsw -mg_levels_pc_type asm -mg_levels_sub_pc_type icc -mg_coarse_redundant_pc_type cholesky
687:    test:
688:      suffix: gdsw_poisson
689:      args: -pde_type Poisson
690:    test:
691:      requires: mumps !complex
692:      suffix: gdsw_poisson_adaptive
693:      args: -pde_type Poisson -mg_levels_gdsw_tolerance 0.01 -ksp_monitor_singular_value -mg_levels_gdsw_userdefined {{0 1}separate output} -mg_levels_gdsw_pseudo_pc_type qr
694:    test:
695:      suffix: gdsw_elast
696:      args: -pde_type Elasticity
697:    test:
698:      requires: hpddm
699:      suffix: gdsw_elast_hpddm
700:      args: -pde_type Elasticity -mg_levels_gdsw_ksp_type hpddm -mg_levels_gdsw_ksp_hpddm_type cg
701:    test:
702:      requires: mumps !complex
703:      suffix: gdsw_elast_adaptive
704:      args: -pde_type Elasticity -mg_levels_gdsw_tolerance 0.01 -ksp_monitor_singular_value -mg_levels_gdsw_userdefined {{0 1}separate output}

706: TEST*/