Ah I get it now, MatSetBlocked has to be set node-wise. I tried this and
it works, thank you.
The other question I had was why are the arguments for MatSetValues()
and MatSetValuesBlocked() set to const PetscInt* and const PetscScalar*
instead of just PetscInt* and PetscScalar* ? I have the typecast there
so my flycheck doesn't keep throwing me warnings on emacs ;)
Thank You,
Nidish
On 8/10/20 5:16 PM, Jed Brown wrote:
Nidish <[email protected]> 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.
--
Nidish