1*b758ae8cSBarry Smith static char help[] = "Tests inserting new block into SBAIJ and BAIJ matrix \n ";
2*b758ae8cSBarry Smith
3*b758ae8cSBarry Smith #include <petscdmda.h>
4*b758ae8cSBarry Smith
main(int argc,char ** argv)5*b758ae8cSBarry Smith int main(int argc, char **argv)
6*b758ae8cSBarry Smith {
7*b758ae8cSBarry Smith DM dm;
8*b758ae8cSBarry Smith Mat A;
9*b758ae8cSBarry Smith PetscInt idm = 0, idn = 8;
10*b758ae8cSBarry Smith PetscScalar v[] = {1, 2, 3, 4};
11*b758ae8cSBarry Smith
12*b758ae8cSBarry Smith PetscFunctionBeginUser;
13*b758ae8cSBarry Smith PetscCall(PetscInitialize(&argc, &argv, NULL, help));
14*b758ae8cSBarry Smith PetscCall(DMDACreate2d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, 4, 4, PETSC_DECIDE, PETSC_DECIDE, 2, 1, NULL, NULL, &dm));
15*b758ae8cSBarry Smith PetscCall(DMSetFromOptions(dm));
16*b758ae8cSBarry Smith PetscCall(DMSetUp(dm));
17*b758ae8cSBarry Smith PetscCall(DMCreateMatrix(dm, &A));
18*b758ae8cSBarry Smith PetscCall(MatSetOption(A, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_FALSE));
19*b758ae8cSBarry Smith PetscCall(MatSetValuesBlocked(A, 1, &idm, 1, &idn, v, INSERT_VALUES));
20*b758ae8cSBarry Smith PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
21*b758ae8cSBarry Smith PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
22*b758ae8cSBarry Smith PetscCall(MatDestroy(&A));
23*b758ae8cSBarry Smith PetscCall(DMDestroy(&dm));
24*b758ae8cSBarry Smith PetscCall(PetscFinalize());
25*b758ae8cSBarry Smith return 0;
26*b758ae8cSBarry Smith }
27*b758ae8cSBarry Smith
28*b758ae8cSBarry Smith /*TEST
29*b758ae8cSBarry Smith
30*b758ae8cSBarry Smith test:
31*b758ae8cSBarry Smith args: -dm_mat_type {{aij baij sbaij}separate output} -mat_view
32*b758ae8cSBarry Smith
33*b758ae8cSBarry Smith test:
34*b758ae8cSBarry Smith suffix: 2
35*b758ae8cSBarry Smith nsize: 2
36*b758ae8cSBarry Smith args: -dm_mat_type {{aij baij sbaij}separate output} -mat_view
37*b758ae8cSBarry Smith
38*b758ae8cSBarry Smith test:
39*b758ae8cSBarry Smith suffix: 3
40*b758ae8cSBarry Smith nsize: 3
41*b758ae8cSBarry Smith args: -dm_mat_type {{aij baij sbaij}separate output} -mat_view
42*b758ae8cSBarry Smith
43*b758ae8cSBarry Smith TEST*/
44