1*26cec326SBarry Smith static char help[] = "Test of setting values in a matrix without preallocation\n\n"; 2*26cec326SBarry Smith 3*26cec326SBarry Smith #include <petscmat.h> 4*26cec326SBarry Smith 5*26cec326SBarry Smith PetscErrorCode ex1_nonsquare_bs1(void) 6*26cec326SBarry Smith { 7*26cec326SBarry Smith Mat A; 8*26cec326SBarry Smith PetscInt M, N, m, n, bs = 1; 9*26cec326SBarry Smith char type[16]; 10*26cec326SBarry Smith PetscBool flg; 11*26cec326SBarry Smith 12*26cec326SBarry Smith /* 13*26cec326SBarry Smith Create the matrix 14*26cec326SBarry Smith */ 15*26cec326SBarry Smith PetscFunctionBegin; 16*26cec326SBarry Smith M = 10; 17*26cec326SBarry Smith N = 12; 18*26cec326SBarry Smith PetscCall(MatCreate(PETSC_COMM_WORLD, &A)); 19*26cec326SBarry Smith PetscCall(PetscOptionsGetString(NULL, NULL, "-type", type, sizeof(type), &flg)); 20*26cec326SBarry Smith if (flg) PetscCall(MatSetType(A, type)); 21*26cec326SBarry Smith else PetscCall(MatSetType(A, MATAIJ)); 22*26cec326SBarry Smith PetscCall(MatSetSizes(A, PETSC_DECIDE, PETSC_DECIDE, M, N)); 23*26cec326SBarry Smith PetscCall(PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL)); 24*26cec326SBarry Smith PetscCall(MatSetBlockSize(A, bs)); 25*26cec326SBarry Smith PetscCall(MatSetFromOptions(A)); 26*26cec326SBarry Smith 27*26cec326SBarry Smith /* 28*26cec326SBarry Smith Get the sizes of the matrix 29*26cec326SBarry Smith */ 30*26cec326SBarry Smith PetscCall(MatGetLocalSize(A, &m, &n)); 31*26cec326SBarry Smith 32*26cec326SBarry Smith /* 33*26cec326SBarry Smith Insert non-zero pattern (e.g. perform a sweep over the grid). 34*26cec326SBarry Smith You can use MatSetValues(), MatSetValuesBlocked() or MatSetValue(). 35*26cec326SBarry Smith */ 36*26cec326SBarry Smith { 37*26cec326SBarry Smith PetscInt ii, jj; 38*26cec326SBarry Smith PetscScalar vv = 22.0; 39*26cec326SBarry Smith 40*26cec326SBarry Smith ii = 3; 41*26cec326SBarry Smith jj = 3; 42*26cec326SBarry Smith PetscCall(MatSetValue(A, ii, jj, vv, INSERT_VALUES)); 43*26cec326SBarry Smith 44*26cec326SBarry Smith ii = 7; 45*26cec326SBarry Smith jj = 4; 46*26cec326SBarry Smith PetscCall(MatSetValue(A, ii, jj, vv, INSERT_VALUES)); 47*26cec326SBarry Smith PetscCall(MatSetValue(A, jj, ii, vv, INSERT_VALUES)); 48*26cec326SBarry Smith 49*26cec326SBarry Smith ii = 9; 50*26cec326SBarry Smith jj = 7; 51*26cec326SBarry Smith PetscCall(MatSetValue(A, ii, jj, vv, INSERT_VALUES)); 52*26cec326SBarry Smith PetscCall(MatSetValue(A, jj, ii, vv, INSERT_VALUES)); 53*26cec326SBarry Smith } 54*26cec326SBarry Smith PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 55*26cec326SBarry Smith PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 56*26cec326SBarry Smith 57*26cec326SBarry Smith PetscCall(PetscViewerPushFormat(PETSC_VIEWER_STDOUT_WORLD, PETSC_VIEWER_ASCII_COMMON)); 58*26cec326SBarry Smith PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 59*26cec326SBarry Smith 60*26cec326SBarry Smith /* 61*26cec326SBarry Smith Insert same location non-zero values into A. 62*26cec326SBarry Smith */ 63*26cec326SBarry Smith { 64*26cec326SBarry Smith PetscInt ii, jj; 65*26cec326SBarry Smith PetscScalar vv; 66*26cec326SBarry Smith 67*26cec326SBarry Smith ii = 3; 68*26cec326SBarry Smith jj = 3; 69*26cec326SBarry Smith vv = 0.3; 70*26cec326SBarry Smith PetscCall(MatSetValue(A, ii, jj, vv, INSERT_VALUES)); 71*26cec326SBarry Smith 72*26cec326SBarry Smith ii = 7; 73*26cec326SBarry Smith jj = 4; 74*26cec326SBarry Smith vv = 3.3; 75*26cec326SBarry Smith PetscCall(MatSetValue(A, ii, jj, vv, INSERT_VALUES)); 76*26cec326SBarry Smith PetscCall(MatSetValue(A, jj, ii, vv, INSERT_VALUES)); 77*26cec326SBarry Smith 78*26cec326SBarry Smith ii = 9; 79*26cec326SBarry Smith jj = 7; 80*26cec326SBarry Smith vv = 4.3; 81*26cec326SBarry Smith PetscCall(MatSetValue(A, ii, jj, vv, INSERT_VALUES)); 82*26cec326SBarry Smith PetscCall(MatSetValue(A, jj, ii, vv, INSERT_VALUES)); 83*26cec326SBarry Smith } 84*26cec326SBarry Smith PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY)); 85*26cec326SBarry Smith PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY)); 86*26cec326SBarry Smith 87*26cec326SBarry Smith PetscCall(MatView(A, PETSC_VIEWER_STDOUT_WORLD)); 88*26cec326SBarry Smith 89*26cec326SBarry Smith PetscCall(MatDestroy(&A)); 90*26cec326SBarry Smith PetscFunctionReturn(PETSC_SUCCESS); 91*26cec326SBarry Smith } 92*26cec326SBarry Smith 93*26cec326SBarry Smith int main(int argc, char **args) 94*26cec326SBarry Smith { 95*26cec326SBarry Smith PetscFunctionBeginUser; 96*26cec326SBarry Smith PetscCall(PetscInitialize(&argc, &args, (char *)0, help)); 97*26cec326SBarry Smith PetscCall(ex1_nonsquare_bs1()); 98*26cec326SBarry Smith PetscCall(PetscFinalize()); 99*26cec326SBarry Smith return 0; 100*26cec326SBarry Smith } 101*26cec326SBarry Smith 102*26cec326SBarry Smith /*TEST 103*26cec326SBarry Smith 104*26cec326SBarry Smith testset: 105*26cec326SBarry Smith args: -bs {{1 2}} -type {{aij baij sbaij}} 106*26cec326SBarry Smith filter: grep -v "type:" 107*26cec326SBarry Smith test: 108*26cec326SBarry Smith test: 109*26cec326SBarry Smith suffix: 2 110*26cec326SBarry Smith nsize: 2 111*26cec326SBarry Smith 112*26cec326SBarry Smith TEST*/ 113