Hello,
Given this block matrix:
A = [A11,A12,A13,A14;
A21,A22,A23,A24;
A31,A32,A33,A34;
A41,A42,A43,A44];
It is trivial to precondition Ax = b with M^-1 = diag(A11^-1, A22^-1,
A33^-1, A44^-1);
My application requires a slightly fancier preconditionner which should
be M^-1 = diag(inv([A11,A12;A21,A22]),inv([A33,A34;A43,A44]));
I'm not sure what is the right tool for this.I've stopped at a 4x4 block matrix, but at scale I have a matrix with few thousands x few thousands blocks (still with the nested 2 x 2 block structure).
1) should I implement a PCSHELL myself, or use a fieldsplit preconditioner with "few thousands / 2" fields (i.e., does PCFIELDSPLIT scale relatively well with the number of fields, or do you recommend it only for "Stokes-like" problems?)? 2) I gave PCFIELDSPLIT a go, but I'm failing miserably. In the attached tarball, I'm loading matrix A on four processes. Each process owns 2 rows of A. I'm thus creating two ISes:
IS Object: 4 MPI processes type: general [0] Number of indices in set 2 [0] 0 0 [0] 1 1 [1] Number of indices in set 2 [1] 0 2 [1] 1 3 [2] Number of indices in set 0 [3] Number of indices in set 0 IS Object: 4 MPI processes type: general [0] Number of indices in set 0 [1] Number of indices in set 0 [2] Number of indices in set 2 [2] 0 4 [2] 1 5 [3] Number of indices in set 2 [3] 0 6 [3] 1 7 which looks good to me. But when I run:$ mpirun -np 4 ./a.out -f almost_full_jolivet -pc_type fieldsplit -pc_fieldsplit_type additive -fieldsplit_ksp_type preonly -fieldsplit_pc_type lu -fieldsplit_pc_factor_mat_solver_package mumps
It fails during setup: [0]PETSC ERROR: Petsc has generated inconsistent data [0]PETSC ERROR: Unhandled case, must have at least two fields, not 1 I would appreciate any feedback and/or fix for this. Thanks, and feel free to move the thread to petsc-users. Pierre
fieldsplit.tar.bz2
Description: application/bzip
