Actual source code: ex29.c
petsc-3.7.4 2016-10-02
2: static char help[] = "Tests VecSetValues() and VecSetValuesBlocked() on MPI vectors.\n\
3: Where atleast a couple of mallocs will occur in the stash code.\n\n";
5: #include <petscvec.h>
9: int main(int argc,char **argv)
10: {
12: PetscMPIInt size;
13: PetscInt i,j,r,n = 50,repeat = 1,bs;
14: PetscScalar val,*vals,zero=0.0;
15: PetscBool subset = PETSC_FALSE,flg;
16: Vec x,y;
18: PetscInitialize(&argc,&argv,(char*)0,help);
19: MPI_Comm_size(PETSC_COMM_WORLD,&size);
20: bs = size;
22: PetscOptionsGetInt(NULL,NULL,"-repeat",&repeat,NULL);
23: PetscOptionsGetBool(NULL,NULL,"-subset",&subset,NULL);
24: PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);
25: VecCreate(PETSC_COMM_WORLD,&x);
26: VecSetSizes(x,PETSC_DECIDE,n*bs);
27: VecSetBlockSize(x,bs);
28: VecSetFromOptions(x);
29: VecDuplicate(x,&y);
31: if (subset) {VecSetOption(x,VEC_SUBSET_OFF_PROC_ENTRIES,PETSC_TRUE);}
33: for (r=0; r<repeat; r++) {
34: /* Assemble the full vector on the first and last iteration, otherwise don't set any values */
35: for (i=0; i<n*bs*(!r || !(repeat-1-r)); i++) {
36: val = i*1.0;
37: VecSetValues(x,1,&i,&val,INSERT_VALUES);
38: }
39: VecAssemblyBegin(x);
40: VecAssemblyEnd(x);
41: if (!r) {VecCopy(x,y);} /* Save result of first assembly */
42: }
44: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
45: VecEqual(x,y,&flg);
46: if (!flg) {PetscPrintf(PETSC_COMM_WORLD,"Vectors from repeat assembly do not match.");}
48: /* Create a new vector because the old stash is a subset. */
49: VecDestroy(&x);
50: VecDuplicate(y,&x);
51: if (subset) {VecSetOption(x,VEC_SUBSET_OFF_PROC_ENTRIES,PETSC_TRUE);}
53: /* Now do the blocksetvalues */
54: VecSet(x,zero);
55: PetscMalloc1(bs,&vals);
56: for (r=0; r<repeat; r++) {
57: /* Assemble the full vector on the first and last iteration, otherwise don't set any values */
58: for (i=0; i<n*(!r || !(repeat-1-r)); i++) {
59: for (j=0; j<bs; j++) vals[j] = (i*bs+j)*1.0;
60: VecSetValuesBlocked(x,1,&i,vals,INSERT_VALUES);
61: }
62: VecAssemblyBegin(x);
63: VecAssemblyEnd(x);
64: if (!r) {VecCopy(x,y);} /* Save result of first assembly */
65: }
67: VecView(x,PETSC_VIEWER_STDOUT_WORLD);
68: VecEqual(x,y,&flg);
69: if (!flg) {PetscPrintf(PETSC_COMM_WORLD,"Vectors from repeat block assembly do not match.");}
71: VecDestroy(&x);
72: VecDestroy(&y);
73: PetscFree(vals);
74: PetscFinalize();
75: return 0;
76: }