#include <petsc.h>

int main(int argc, char** argv) {
    PetscInitialize(&argc, &argv, NULL, NULL);
    double one = 1;
    Mat A0, A1, A2;
    MatCreateDense(PETSC_COMM_WORLD,1,1,1,1,&one,&A0);
    MatCreateDense(PETSC_COMM_WORLD,1,1,1,1,&one,&A1);
    MatCreateDense(PETSC_COMM_WORLD,1,1,1,1,&one,&A2);
    Mat* a = new Mat[3 * 3]();
    a[0] = A0;
    a[4] = A1;
    a[8] = A2;
    Mat B;
    MatCreateNest(PETSC_COMM_WORLD,3,NULL,3,NULL,a,&B);
    KSP ksp;
    KSPCreate(PETSC_COMM_WORLD,&ksp);
    KSPSetOperators(ksp,B,B);
    KSPSetFromOptions(ksp);
    PC pc;
    KSPGetPC(ksp,&pc);
    PCSetType(pc,PCFIELDSPLIT);
    IS is0, is1, is2;
    ISCreateStride(PETSC_COMM_WORLD,1,0,1,&is0);
    ISCreateStride(PETSC_COMM_WORLD,1,1,1,&is1);
    ISCreateStride(PETSC_COMM_WORLD,1,2,1,&is2);
    PCFieldSplitSetIS(pc,NULL,is0);
    PCFieldSplitSetIS(pc,NULL,is1);
    PCFieldSplitSetIS(pc,NULL,is2);
    KSPSetUp(ksp);
    ISDestroy(&is0);
    ISDestroy(&is1);
    ISDestroy(&is2);
    KSPView(ksp, PETSC_VIEWER_STDOUT_SELF);
    KSPDestroy(&ksp);
    delete [] a;
    MatDestroy(&B);
    MatDestroy(&A0);
    MatDestroy(&A1);
    MatDestroy(&A2);
    PetscFinalize();
    return 0;
}
