1c4762a1bSJed Brown 2c4762a1bSJed Brown static char help[] = "Tests MatTranspose() and MatEqual() for MPIAIJ matrices.\n\n"; 3c4762a1bSJed Brown 4c4762a1bSJed Brown #include <petscmat.h> 5c4762a1bSJed Brown 6c4762a1bSJed Brown int main(int argc,char **argv) 7c4762a1bSJed Brown { 8c4762a1bSJed Brown Mat A,B; 9c4762a1bSJed Brown PetscInt m = 7,n,i,rstart,rend,cols[3]; 10c4762a1bSJed Brown PetscErrorCode ierr; 11c4762a1bSJed Brown PetscScalar v[3]; 12c4762a1bSJed Brown PetscBool equal; 13c4762a1bSJed Brown const char *eq[2]; 14c4762a1bSJed Brown 15c4762a1bSJed Brown ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 16*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetInt(NULL,NULL,"-m",&m,NULL)); 17c4762a1bSJed Brown n = m; 18c4762a1bSJed Brown 19c4762a1bSJed Brown /* ------- Assemble matrix, --------- */ 20c4762a1bSJed Brown 21*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatCreateAIJ(PETSC_COMM_WORLD,PETSC_DECIDE,PETSC_DECIDE,m,n,0,0,0,0,&A)); 22*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatSetOption(A,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_FALSE)); 23*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatGetOwnershipRange(A,&rstart,&rend)); 24c4762a1bSJed Brown if (!rstart) { 25c4762a1bSJed Brown cols[0] = 0; 26c4762a1bSJed Brown cols[1] = 1; 27c4762a1bSJed Brown v[0] = 2.0; v[1] = -1.0; 28*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatSetValues(A,1,&rstart,2,cols,v,INSERT_VALUES)); 29c4762a1bSJed Brown rstart++; 30c4762a1bSJed Brown } 31c4762a1bSJed Brown if (rend == m) { 32c4762a1bSJed Brown rend--; 33c4762a1bSJed Brown cols[0] = rend-1; 34c4762a1bSJed Brown cols[1] = rend; 35c4762a1bSJed Brown v[0] = -1.0; v[1] = 2.0; 36*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatSetValues(A,1,&rend,2,cols,v,INSERT_VALUES)); 37c4762a1bSJed Brown } 38c4762a1bSJed Brown v[0] = -1.0; v[1] = 2.0; v[2] = -1.0; 39c4762a1bSJed Brown for (i=rstart; i<rend; i++) { 40c4762a1bSJed Brown cols[0] = i-1; 41c4762a1bSJed Brown cols[1] = i; 42c4762a1bSJed Brown cols[2] = i+1; 43*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatSetValues(A,1,&i,3,cols,v,INSERT_VALUES)); 44c4762a1bSJed Brown } 45*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY)); 46*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY)); 47c4762a1bSJed Brown 48*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatTranspose(A,MAT_INITIAL_MATRIX,&B)); 49c4762a1bSJed Brown 50*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatEqual(A,B,&equal)); 51c4762a1bSJed Brown 52c4762a1bSJed Brown eq[0] = "not equal"; 53c4762a1bSJed Brown eq[1] = "equal"; 54*5f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"Matrices are %s\n",eq[equal])); 55c4762a1bSJed Brown 56*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatTranspose(A,MAT_REUSE_MATRIX,&B)); 57*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatEqual(A,B,&equal)); 58*5f80ce2aSJacob Faibussowitsch if (!equal) CHKERRQ(PetscPrintf(PETSC_COMM_WORLD,"MatTranspose with MAT_REUSE_MATRIX failed")); 59c4762a1bSJed Brown 60c4762a1bSJed Brown /* Free data structures */ 61*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatDestroy(&A)); 62*5f80ce2aSJacob Faibussowitsch CHKERRQ(MatDestroy(&B)); 63c4762a1bSJed Brown 64c4762a1bSJed Brown ierr = PetscFinalize(); 65c4762a1bSJed Brown return ierr; 66c4762a1bSJed Brown } 67c4762a1bSJed Brown 68c4762a1bSJed Brown /*TEST 69c4762a1bSJed Brown 70c4762a1bSJed Brown test: 71c4762a1bSJed Brown 72c4762a1bSJed Brown test: 73c4762a1bSJed Brown suffix: 2 74c4762a1bSJed Brown nsize: 2 75c4762a1bSJed Brown output_file: output/ex58_1.out 76c4762a1bSJed Brown 77c4762a1bSJed Brown TEST*/ 78