[petsc-users] 3-dimension to matrix

[Jed Brown]

Joon hee Choi <choi240 at purdue.edu> writes:

> Hi all,
>
> I am trying to set up a matrix with 3d coordinates(x,y,z) and
> value. If X=max(x), Y=max(y) and Z=max(z), then the matrix is as
> follows: A = [ A1 A2 ... AZ ] (each A* is X-by-Y block matrix and A
> consists of Z block matrices)

Are these decoupled diagonal blocks or a sequence of blocks in the row.
That is, is the full size of the matrix (X, Y*Z) or (X*Z, Y*Z).  From
the code below, it looks like you intend for it to be the former.  What
parallel distribution do you intend for the input vector (of dimension
Y*Z) and the output vector (of dimension X)?

> When I tried to set up the matrix using MatSetValues with each
> coordinates and value, I could get the matrix, but it took 10
> minutes. 

This is almost certainly due to not preallocating.

> So I tried to set up the matrix again using
> MatSetLocalToGlobalMapping(Mat x, ISLocalToGlobalMapping rmapping,
> ISLocalToGlobalMapping cmapping), but ISLocalToGloalMapping was minus
> because my dataset is too big (X:26M, Y:26M, Z:48M, nonzero:144M) and
> I got errors. Is there the way to set up the matrix quickly in petsc?
> If so, please let me know the way. And my code using
> MatSetLocalToGlobalMapping following:
>
>
> ISLocalToGlobalMappingCreate(PETSC_COMM_SELF, 1, &zr, PETSC_OWN_POINTER, &irow);

This has only one entry.

> ISLocalToGlobalMappingUnBlock(irow, X, &orow);

How did you set the block size of the matrix?

>
> for (z=0; z<Z; i++)
> {
>     ISLocalToGlobalMappingCreate(PETSC_COMM_SELF, 1, &z, PETSC_OWN_POINTER, &icol);

This mapping also has only one entry.  It is not intended for you to
change the local-to-global mapping on a matrix.

>     MatSetLocalToGlobalMappingBlock(A, irow, icol);
>     ISLocalToGlobalMappingUnBlock(icol, Y, &ocol);
>     MatSetLocalToGlobalMapping(A, orow, ocol);
>
>     while(!mapX[z].empty())
>     {
>        x = mapX[z].front();
>        y = mapY[z].front();

So instead of changing it, just add the offset to y here.

>        val = mapVal[z].front();
>        MatSetValuesLocal(A, 1, &x, 1, &y, &val, INSERT_VALUES);
>        mapX[z].pop_front();
>        mapY[z].pop_front();
>        mapVal[z].pop_front();
>     }
> }
>
> Thank you,
>
> Joon