[petsc-users] MatSetValues vs MatSetValuesBlocked

[Jed Brown]

Nidish <nb25 at rice.edu> writes:

> It's a 1D model with displacements and rotations as DoFs at each node.
>
> I couldn't find much in the manual on MatSetBlockSize - could you 
> provide some more information on its use?
>
> I thought since I've setup the system using DMDACreate1d (I've given 
> 2dofs per node and a stencil width of 1 there), the matrix object should 
> have the nonzero elements preallocated. Here's the call to DMDACreate1d:
>
>    DMDACreate1d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, N, 2, 1, NULL, &mshdm);

Ah, that will set the block size, but then it'll expect elstiff to be an 8x8 matrix where you've only passed 4x4.

    idx[0] = 2*e; idx[1] = 2*e+1; idx[2] = 2*e+2; idx[3] = 2*e+3;

    MatSetValuesBlocked(jac, 4, (const PetscInt*)idx, 4, (const PetscInt*)idx,
    			(const PetscScalar*)elstiff, ADD_VALUES);

You don't need the casts in either case, BTW.  You probably want something like this.

    idx[0] = e; idx[1] = e + 1;

    MatSetValuesBlocked(jac, 2, idx, 2, idx, elstiff, ADD_VALUES);

Also, it might be more convenient to call MatSetValuesBlockedStencil(), especially if you move to a multi-dimensional problem at some point.