[petsc-dev] MatSetValuesBlocked not inserting values

[Jed Brown]

Sure, but could you have computed them as a negative number here?

Jacob Faibussowitsch <jacob.fai at gmail.com> writes:

> Supposed to be the block row and column indices. So they are the atom global IDs. 
>
> Best regards,
>
> Jacob Faibussowitsch
> (Jacob Fai - booss - oh - vitch)
> Cell: (312) 694-3391
>
>> On Apr 8, 2020, at 4:29 PM, Jed Brown <jed at jedbrown.org> wrote:
>> 
>> What are idxm[0] and idxn[0] in this code?
>> 
>> Jacob Faibussowitsch <jacob.fai at gmail.com> writes:
>> 
>>> Hello All,
>>> 
>>> I am using MATSBAIJ to make a symmetric hamiltonian matrix, but for some reason MatSetValuesBlocked isn’t inserting anything. 
>>> 
>>> Setup:
>>> 
>>> I want a 4n x 4n symmetric matrix, where each MatSetValuesBlocked inserts a 4x4 sub-matrix of values into the global matrix, positioned by the global ID of two atoms that I am extracting. So the 4x4 interaction sub matrix between atom 1 and atom 4 will be in rows 0-4 and cols 16-20.
>>> 
>>> My non-working implementation:
>>> 
>>> PetscInt bs = 4;
>>> 
>>> ierr = MatCreateSBAIJ(comm, bs, bs*nLocalAtoms, bs*nLocalAtoms, PETSC_DETERMINE, PETSC_DETERMINE, 1, NULL, 1, NULL, H);CHKERRQ(ierr);
>>> ierr = MatSetOption(*H, MAT_SYMMETRIC, PETSC_TRUE);CHKERRQ(ierr);
>>> ierr = MatSetUp(*H);CHKERRQ(ierr);                                    
>>> .
>>> .
>>> .
>>> START LOOP {
>>> 
>>> Calculations...
>>> 
>>> const PetscScalar  valmat[16] = {Es[0], Es[1], Es[2], Es[3], Ex[0], Ex[1], Ex[2], Ex[3], Ey[0], Ey[1], Ey[2], Ey[3], Ez[0], Ez[1], Ez[2], Ez[3]};
>>> PetscInt                  idxm[1], idxn[1];
>>> PetscInt                  m = 1, n = 1;
>>> 
>>> if (globalIDs[i] < globalIDs[neighidx[0]]-1) {
>>> 	// Insert in upper triangular
>>>        idxm[0] = globalIDs[i]-1;
>>>        idxn[0] = globalIDs[neighidx[0]]-1;
>>> } else {
>>>        // swap so we insert on upper triangular                                                    
>>>        idxn[0] = globalIDs[i]-1;
>>>        idxm[0] = globalIDs[neighidx[0]]-1;
>>> }
>>> ierr = MatSetValuesBlocked(*H, m, (const PetscInt *) idxm, n, (const PetscInt *) idxn, valmat, INSERT_VALUES);CHKERRQ(ierr);
>>> } END LOOP
>>> 
>>> ierr = MatAssemblyBegin(*H, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
>>> ierr = MatAssemblyEnd(*H, MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
>>> 
>>> RETURN H
>>> 
>>> Output from -mat_view:
>>> 
>>> Mat Object: 2 MPI processes
>>>  type: mpisbaij
>>> row 0:
>>> row 1:
>>> row 2:
>>> row 3:
>>> row 4:
>>> row 5:
>>> row 6:
>>> row 7:
>>> row 8:
>>> row 9:
>>> row 10:
>>> row 11:
>>> row 12:
>>> row 13:
>>> row 14:
>>> row 15:
>>> 
>>> I also checked -info :mat: and it also reports X unneeded storage space but more importantly, 0 used. I am not sure what is going wrong here. 
>>> 
>>> Best regards,
>>> 
>>> Jacob Faibussowitsch
>>> (Jacob Fai - booss - oh - vitch)
>>> Cell: (312) 694-3391