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    multi_element;
111:   PetscBool    dirbc;
112:   PetscBool    per[3];
113:   PetscBool    test;
114:   PetscScalar *elemMat;
115:   PetscBool    use_composite_pc;
116:   PetscBool    random_initial_guess;
117:   PetscBool    random_real;
118: } AppCtx;

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

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

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

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

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

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

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

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

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

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

273:   if (user.multi_element) {
274:     ISLocalToGlobalMapping mapn;
275:     PetscInt              *e_glo = NULL;

277:     PetscCall(PetscMalloc1(nel * nen, &e_glo));
278:     PetscCall(ISLocalToGlobalMappingApplyBlock(map, nen * nel, e_loc, e_glo));
279:     PetscCall(ISLocalToGlobalMappingCreate(PetscObjectComm((PetscObject)map), user.dof, nen * nel, e_glo, PETSC_OWN_POINTER, &mapn));
280:     PetscCall(MatISSetAllowRepeated(A, PETSC_TRUE));
281:     PetscCall(MatSetLocalToGlobalMapping(A, mapn, mapn));
282:     PetscCall(ISLocalToGlobalMappingViewFromOptions(mapn, NULL, "-multi_view"));
283:     PetscCall(MatSetDM(A, NULL));
284:     PetscCall(ISLocalToGlobalMappingDestroy(&mapn));
285:   }

287:   /* we reorder the indices since the element matrices are given in lexicographic order,
288:      whereas the elements indices returned by DMDAGetElements follow the usual FEM ordering
289:      i.e., element matrices     DMDA ordering
290:                2---3              3---2
291:               /   /              /   /
292:              0---1              0---1
293:   */
294:   for (i = 0; i < nel; ++i) {
295:     PetscInt ord[8] = {0, 1, 3, 2, 4, 5, 7, 6};
296:     PetscInt j, idxs[8];

298:     PetscCheck(nen <= 8, PETSC_COMM_WORLD, PETSC_ERR_SUP, "Not coded");
299:     if (!user.useglobal) {
300:       if (user.multi_element) {
301:         for (j = 0; j < nen; j++) idxs[j] = i * nen + ord[j];
302:       } else {
303:         for (j = 0; j < nen; j++) idxs[j] = e_loc[i * nen + ord[j]];
304:       }
305:       PetscCall(MatSetValuesBlockedLocal(A, nen, idxs, nen, idxs, user.elemMat, ADD_VALUES));
306:     } else {
307:       for (j = 0; j < nen; j++) idxs[j] = e_glo[i * nen + ord[j]];
308:       PetscCall(MatSetValuesBlocked(A, nen, idxs, nen, idxs, user.elemMat, ADD_VALUES));
309:     }
310:   }
311:   PetscCall(DMDARestoreElements(da, &nel, &nen, &e_loc));
312:   PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
313:   PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
314:   PetscCall(MatViewFromOptions(A, NULL, "-A_mat_view"));
315:   PetscCall(MatSetOption(A, MAT_SPD, PETSC_TRUE));
316:   PetscCall(MatSetOption(A, MAT_SPD_ETERNAL, PETSC_TRUE));

318:   /* Boundary conditions */
319:   zero = NULL;
320:   if (user.dirbc) { /* fix one side of DMDA */
321:     Vec          nat, glob;
322:     PetscScalar *vals;
323:     PetscInt     n, *idx, j, st;

325:     n = PetscGlobalRank ? 0 : (user.cells[1] + 1) * (user.cells[2] + 1);
326:     PetscCall(ISCreateStride(PETSC_COMM_WORLD, n, 0, user.cells[0] + 1, &zero));
327:     if (user.dof > 1) { /* zero all components */
328:       const PetscInt *idx;
329:       IS              bzero;

331:       PetscCall(ISGetIndices(zero, (const PetscInt **)&idx));
332:       PetscCall(ISCreateBlock(PETSC_COMM_WORLD, user.dof, n, idx, PETSC_COPY_VALUES, &bzero));
333:       PetscCall(ISRestoreIndices(zero, (const PetscInt **)&idx));
334:       PetscCall(ISDestroy(&zero));
335:       zero = bzero;
336:     }
337:     /* map indices from natural to global */
338:     PetscCall(DMDACreateNaturalVector(da, &nat));
339:     PetscCall(ISGetLocalSize(zero, &n));
340:     PetscCall(PetscMalloc1(n, &vals));
341:     for (i = 0; i < n; i++) vals[i] = 1.0;
342:     PetscCall(ISGetIndices(zero, (const PetscInt **)&idx));
343:     PetscCall(VecSetValues(nat, n, idx, vals, INSERT_VALUES));
344:     PetscCall(ISRestoreIndices(zero, (const PetscInt **)&idx));
345:     PetscCall(PetscFree(vals));
346:     PetscCall(VecAssemblyBegin(nat));
347:     PetscCall(VecAssemblyEnd(nat));
348:     PetscCall(DMCreateGlobalVector(da, &glob));
349:     PetscCall(DMDANaturalToGlobalBegin(da, nat, INSERT_VALUES, glob));
350:     PetscCall(DMDANaturalToGlobalEnd(da, nat, INSERT_VALUES, glob));
351:     PetscCall(VecDestroy(&nat));
352:     PetscCall(ISDestroy(&zero));
353:     PetscCall(VecGetLocalSize(glob, &n));
354:     PetscCall(PetscMalloc1(n, &idx));
355:     PetscCall(VecGetOwnershipRange(glob, &st, NULL));
356:     PetscCall(VecGetArray(glob, &vals));
357:     for (i = 0, j = 0; i < n; i++)
358:       if (PetscRealPart(vals[i]) == 1.0) idx[j++] = i + st;
359:     PetscCall(VecRestoreArray(glob, &vals));
360:     PetscCall(VecDestroy(&glob));
361:     PetscCall(ISCreateGeneral(PETSC_COMM_WORLD, j, idx, PETSC_OWN_POINTER, &zero));
362:     PetscCall(MatZeroRowsColumnsIS(A, zero, 1.0, NULL, NULL));
363:   } else {
364:     switch (user.pde) {
365:     case PDE_POISSON:
366:       PetscCall(MatNullSpaceCreate(PETSC_COMM_WORLD, PETSC_TRUE, 0, NULL, &nullsp));
367:       break;
368:     case PDE_ELASTICITY:
369:       PetscCall(MatNullSpaceCreateRigidBody(xcoor, &nullsp));
370:       break;
371:     }
372:     /* with periodic BC and Elasticity, just the displacements are in the nullspace
373:        this is no harm since we eliminate all the components of the rhs */
374:     PetscCall(MatSetNullSpace(A, nullsp));
375:   }

377:   PetscCall(PetscOptionsHasName(NULL, NULL, "-assembled_view", &flg));
378:   if (flg) {
379:     Mat AA;

381:     PetscCall(MatConvert(A, MATAIJ, MAT_INITIAL_MATRIX, &AA));
382:     PetscCall(MatViewFromOptions(AA, NULL, "-assembled_view"));
383:     PetscCall(MatDestroy(&AA));
384:   }

386:   /* Attach near null space for elasticity */
387:   if (user.pde == PDE_ELASTICITY) {
388:     MatNullSpace nearnullsp;

390:     PetscCall(MatNullSpaceCreateRigidBody(xcoor, &nearnullsp));
391:     PetscCall(MatSetNearNullSpace(A, nearnullsp));
392:     PetscCall(MatNullSpaceDestroy(&nearnullsp));
393:   }

395:   /* we may want to use MG for the local solvers: attach local nearnullspace to the local matrices */
396:   PetscCall(DMGetCoordinatesLocal(da, &xcoorl));
397:   PetscCall(PetscObjectTypeCompare((PetscObject)A, MATIS, &ismatis));
398:   if (ismatis) {
399:     MatNullSpace lnullsp = NULL;
400:     Mat          lA;

402:     PetscCall(MatISGetLocalMat(A, &lA));
403:     if (user.pde == PDE_ELASTICITY) {
404:       Vec                    lc;
405:       ISLocalToGlobalMapping l2l;
406:       IS                     is;
407:       const PetscScalar     *a;
408:       const PetscInt        *idxs;
409:       PetscInt               n, bs;

411:       /* when using a DMDA, the local matrices have an additional local-to-local map
412:          that maps from the DA local ordering to the ordering induced by the elements */
413:       PetscCall(MatGetLocalToGlobalMapping(lA, &l2l, NULL));
414:       if (l2l) {
415:         PetscCall(MatCreateVecs(lA, &lc, NULL));
416:         PetscCall(VecSetLocalToGlobalMapping(lc, l2l));

418:         PetscCall(VecSetOption(lc, VEC_IGNORE_NEGATIVE_INDICES, PETSC_TRUE));
419:         PetscCall(VecGetLocalSize(xcoorl, &n));
420:         PetscCall(VecGetBlockSize(xcoorl, &bs));
421:         PetscCall(ISCreateStride(PETSC_COMM_SELF, n / bs, 0, 1, &is));
422:         PetscCall(ISGetIndices(is, &idxs));
423:         PetscCall(VecGetArrayRead(xcoorl, &a));
424:         PetscCall(VecSetValuesBlockedLocal(lc, n / bs, idxs, a, INSERT_VALUES));
425:         PetscCall(VecAssemblyBegin(lc));
426:         PetscCall(VecAssemblyEnd(lc));
427:         PetscCall(VecRestoreArrayRead(xcoorl, &a));
428:         PetscCall(ISRestoreIndices(is, &idxs));
429:         PetscCall(ISDestroy(&is));
430:         PetscCall(MatNullSpaceCreateRigidBody(lc, &lnullsp));
431:         PetscCall(VecDestroy(&lc));
432:       }
433:     } else if (user.pde == PDE_POISSON) {
434:       PetscCall(MatNullSpaceCreate(PETSC_COMM_SELF, PETSC_TRUE, 0, NULL, &lnullsp));
435:     }
436:     PetscCall(MatSetNearNullSpace(lA, lnullsp));
437:     PetscCall(MatNullSpaceDestroy(&lnullsp));
438:     PetscCall(MatISRestoreLocalMat(A, &lA));
439:   }

441:   PetscCall(KSPCreate(PETSC_COMM_WORLD, &ksp));
442:   PetscCall(KSPSetOperators(ksp, A, A));
443:   PetscCall(KSPSetType(ksp, KSPCG));
444:   PetscCall(KSPGetPC(ksp, &pc));
445:   if (ismatis) {
446:     if (user.use_composite_pc) {
447:       PC  pcksp, pcjacobi;
448:       KSP ksprich;
449:       PetscCall(PCSetType(pc, PCCOMPOSITE));
450:       PetscCall(PCCompositeSetType(pc, PC_COMPOSITE_MULTIPLICATIVE));
451:       PetscCall(PCCompositeAddPCType(pc, PCBDDC));
452:       PetscCall(PCCompositeAddPCType(pc, PCKSP));
453:       PetscCall(PCCompositeGetPC(pc, 1, &pcksp));
454:       PetscCall(PCKSPGetKSP(pcksp, &ksprich));
455:       PetscCall(KSPSetType(ksprich, KSPRICHARDSON));
456:       PetscCall(KSPSetTolerances(ksprich, PETSC_DEFAULT, PETSC_DEFAULT, PETSC_DEFAULT, 1));
457:       PetscCall(KSPSetNormType(ksprich, KSP_NORM_NONE));
458:       PetscCall(KSPSetConvergenceTest(ksprich, KSPConvergedSkip, NULL, NULL));
459:       PetscCall(KSPGetPC(ksprich, &pcjacobi));
460:       PetscCall(PCSetType(pcjacobi, PCJACOBI));
461:     } else {
462:       PetscCall(PCSetType(pc, PCBDDC));
463:     }
464:   }
465:   PetscCall(KSPSetFromOptions(ksp));
466:   PetscCall(PetscLogStagePush(stages[0]));
467:   PetscCall(KSPSetUp(ksp));
468:   PetscCall(PetscLogStagePop());

470:   PetscCall(MatCreateVecs(A, &x, &b));
471:   if (user.random_initial_guess) {
472:     /* Solving A x = 0 with random initial guess allows Arnoldi to run for more iterations, thereby yielding a more
473:      * complete Hessenberg matrix and more accurate eigenvalues. */
474:     PetscCall(VecZeroEntries(b));
475:     PetscCall(VecSetRandom(x, NULL));
476:     if (user.random_real) PetscCall(VecRealPart(x));
477:     if (nullsp) PetscCall(MatNullSpaceRemove(nullsp, x));
478:     PetscCall(KSPSetInitialGuessNonzero(ksp, PETSC_TRUE));
479:     PetscCall(KSPSetComputeEigenvalues(ksp, PETSC_TRUE));
480:     PetscCall(KSPGMRESSetRestart(ksp, 100));
481:   } else {
482:     PetscCall(VecSetRandom(b, NULL));
483:     if (user.random_real) PetscCall(VecRealPart(x));
484:     if (nullsp) PetscCall(MatNullSpaceRemove(nullsp, b));
485:   }
486:   PetscCall(PetscLogStagePush(stages[1]));
487:   PetscCall(KSPSolve(ksp, b, x));
488:   PetscCall(PetscLogStagePop());

490:   /* cleanup */
491:   PetscCall(VecDestroy(&x));
492:   PetscCall(VecDestroy(&b));
493:   PetscCall(ISDestroy(&zero));
494:   PetscCall(PetscFree(e_glo));
495:   PetscCall(MatNullSpaceDestroy(&nullsp));
496:   PetscCall(KSPDestroy(&ksp));
497:   PetscCall(MatDestroy(&A));
498:   PetscCall(DMDestroy(&da));
499:   PetscCall(PetscFinalize());
500:   return 0;
501: }

503: /*TEST

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

682:  testset:
683:    nsize: 2
684:    output_file: output/ex71_aij_dmda_preall.out
685:    filter: sed -e "s/CONVERGED_RTOL iterations 7/CONVERGED_RTOL iterations 6/g"
686:    args: -pde_type Poisson -dim 1 -cells 6 -pc_type none -ksp_converged_reason
687:    test:
688:      suffix: aijviennacl_dmda_preall
689:      requires: viennacl
690:      args: -dm_mat_type aijviennacl -dm_preallocate_only {{0 1}} -dirichlet {{0 1}}
691:    # -dm_preallocate_only 0 is broken
692:    test:
693:      suffix: aijcusparse_dmda_preall
694:      requires: cuda
695:      args: -dm_mat_type aijcusparse -dm_preallocate_only -dirichlet {{0 1}}
696:    test:
697:      suffix: aij_dmda_preall
698:      args: -dm_mat_type aij -dm_preallocate_only {{0 1}} -dirichlet {{0 1}}
699:  testset:
700:    nsize: 4
701:    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
702:    test:
703:      args: -sub_0_pc_bddc_interface_ext_type lump
704:      suffix: composite_bddc_lumped
705:    test:
706:      requires: !single
707:      args: -sub_0_pc_bddc_interface_ext_type dirichlet
708:      suffix: composite_bddc_dirichlet

710:  testset:
711:    nsize: 8
712:    filter: grep -v "variant HERMITIAN"
713:    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
714:    test:
715:      suffix: gdsw_poisson
716:      args: -pde_type Poisson
717:    test:
718:      requires: mumps !complex
719:      suffix: gdsw_poisson_adaptive
720:      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
721:    test:
722:      suffix: gdsw_elast
723:      args: -pde_type Elasticity
724:    test:
725:      requires: hpddm
726:      suffix: gdsw_elast_hpddm
727:      args: -pde_type Elasticity -mg_levels_gdsw_ksp_type hpddm -mg_levels_gdsw_ksp_hpddm_type cg
728:    test:
729:      requires: mumps !complex
730:      suffix: gdsw_elast_adaptive
731:      args: -pde_type Elasticity -mg_levels_gdsw_tolerance 0.01 -ksp_monitor_singular_value -mg_levels_gdsw_userdefined {{0 1}separate output}

733: TEST*/