2: static char help[] = "Tests MatConvert(), MatLoad() for MATELEMENTAL interface.\n\n";
  3: /*
  4:  Example:
  5:    mpiexec -n <np> ./ex173 -fA <A_data> -fB <B_data> -orig_mat_type <type> -orig_mat_type <mat_type>
  6: */
  7: 
  8:  #include <petscmat.h>
  9:  #include <petscmatelemental.h>

 11: int main(int argc,char **args)
 12: {
 13:   Mat            A,Ae,B,Be;
 15:   PetscViewer    view;
 16:   char           file[2][PETSC_MAX_PATH_LEN];
 17:   PetscBool      flg,flgB,isElemental,isDense,isAij,isSbaij;
 18:   PetscScalar    one = 1.0;
 19:   PetscMPIInt    rank,size;
 20:   PetscInt       M,N;

 22:   PetscInitialize(&argc,&args,(char*)0,help);if (ierr) return ierr;
 23: #if !defined(PETSC_HAVE_ELEMENTAL)
 24:   SETERRQ(PETSC_COMM_WORLD,1,"This example requires ELEMENTAL");
 25: #endif
 26:   MPI_Comm_rank(PETSC_COMM_WORLD,&rank);
 27:   MPI_Comm_size(PETSC_COMM_WORLD,&size);

 29:   /* Load PETSc matrices */
 30:   PetscOptionsGetString(NULL,NULL,"-fA",file[0],PETSC_MAX_PATH_LEN,NULL);
 31:   PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[0],FILE_MODE_READ,&view);
 32:   MatCreate(PETSC_COMM_WORLD,&A);
 33:   MatSetOptionsPrefix(A,"orig_");
 34:   MatSetType(A,MATAIJ);
 35:   MatSetFromOptions(A);
 36:   MatLoad(A,view);
 37:   PetscViewerDestroy(&view);

 39:   PetscOptionsGetString(NULL,NULL,"-fB",file[1],PETSC_MAX_PATH_LEN,&flgB);
 40:   if (flgB) {
 41:     PetscViewerBinaryOpen(PETSC_COMM_WORLD,file[1],FILE_MODE_READ,&view);
 42:     MatCreate(PETSC_COMM_WORLD,&B);
 43:     MatSetOptionsPrefix(B,"orig_");
 44:     MatSetType(B,MATAIJ);
 45:     MatSetFromOptions(B);
 46:     MatLoad(B,view);
 47:     PetscViewerDestroy(&view);
 48:   } else {
 49:     /* Create matrix B = I */
 50:     PetscInt rstart,rend,i;
 51:     MatGetSize(A,&M,&N);
 52:     MatGetOwnershipRange(A,&rstart,&rend);

 54:     MatCreate(PETSC_COMM_WORLD,&B);
 55:     MatSetOptionsPrefix(B,"orig_");
 56:     MatSetSizes(B,PETSC_DECIDE,PETSC_DECIDE,M,N);
 57:     MatSetType(B,MATAIJ);
 58:     MatSetFromOptions(B);
 59:     MatSetUp(B);
 60:     for (i=rstart; i<rend; i++) {
 61:       MatSetValues(B,1,&i,1,&i,&one,ADD_VALUES);
 62:     }
 63:     MatAssemblyBegin(B,MAT_FINAL_ASSEMBLY);
 64:     MatAssemblyEnd(B,MAT_FINAL_ASSEMBLY);
 65:   }

 67:   PetscObjectTypeCompare((PetscObject)A,MATELEMENTAL,&isElemental);
 68:   if (isElemental) {
 69:     Ae = A;
 70:     Be = B;
 71:     isDense = isAij = isSbaij = PETSC_FALSE;
 72:   } else { /* Convert AIJ/DENSE/SBAIJ matrices into Elemental matrices */
 73:     if (size == 1) {
 74:       PetscObjectTypeCompare((PetscObject)A,MATSEQDENSE,&isDense);
 75:       PetscObjectTypeCompare((PetscObject)A,MATSEQAIJ,&isAij);
 76:       PetscObjectTypeCompare((PetscObject)A,MATSEQSBAIJ,&isSbaij);
 77:     } else {
 78:       PetscObjectTypeCompare((PetscObject)A,MATMPIDENSE,&isDense);
 79:       PetscObjectTypeCompare((PetscObject)A,MATMPIAIJ,&isAij);
 80:        PetscObjectTypeCompare((PetscObject)A,MATMPISBAIJ,&isSbaij);
 81:     }

 83:     if (!rank) {
 84:       if (isDense) {
 85:         printf(" Convert DENSE matrices A and B into Elemental matrix... \n");
 86:       } else if (isAij) {
 87:         printf(" Convert AIJ matrices A and B into Elemental matrix... \n");
 88:       } else if (isSbaij) {
 89:         printf(" Convert SBAIJ matrices A and B into Elemental matrix... \n");
 90:       } else SETERRQ(PetscObjectComm((PetscObject)A),PETSC_ERR_SUP,"Not supported yet");
 91:     }
 92:     MatConvert(A, MATELEMENTAL, MAT_INITIAL_MATRIX, &Ae);
 93:     MatConvert(B, MATELEMENTAL, MAT_INITIAL_MATRIX, &Be);

 95:     /* Test accuracy */
 96:     MatMultEqual(A,Ae,5,&flg);
 97:     if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"A != A_elemental.");
 98:     MatMultEqual(B,Be,5,&flg);
 99:     if (!flg) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_NOTSAMETYPE,"B != B_elemental.");
100:   }

102:   if (!isElemental) {
103:     MatDestroy(&Ae);
104:     MatDestroy(&Be);

106:     /* Test MAT_REUSE_MATRIX which is only supported for inplace conversion */
107:     MatConvert(A, MATELEMENTAL, MAT_INPLACE_MATRIX, &A);
108:     //MatView(A,PETSC_VIEWER_STDOUT_WORLD);
109:   }

111:   MatDestroy(&A);
112:   MatDestroy(&B);
113:   PetscFinalize();
114:   return ierr;
115: }