pyccel / psydac

Python 3 library for isogeometric analysis
https://pyccel.github.io/psydac/
MIT License
47 stars 17 forks source link

Inefficient conversion from PETSC to array #389

Open e-moral-sanchez opened 6 months ago

e-moral-sanchez commented 6 months ago

At the moment when we convert a distributed StencilVector or StencilMatrix or BlockLinearOperator to PETSc.Vec or PETSc.Mat format with the topetsc() method, the corresponding (dense) arrays cannot be compared directly. This is because PETSc distributes the vectors or matrices differently than Psydac.

In a distributed Psydac StencilVector or StencilMatrix a process may own non-contiguous rows (for instance, if the domain is 2D and it is distributed by 4 processes). At the moment, a process with the PETSC distributed Vec or Mat object will always own contiguous rows. This is possible because when calling the assembly of the PETSc Vec or Mat object, there is an exchange of global communication.

The reason why PETSc does this is because it does not have any information about how the domain is decomposed. Consequently, it simply decomposes the vector or matrix by rows uniformly.

The major consequence of this issue is that the function petsc_to_array becomes unnecessarily expensive, since we need to gather the array from all the processes and redistribute it according to Psydac's domain partition.

A potential solution might be to create a PETSc.DM object that contains the domain decomposition information from Psydac.

Note: In order to get the global indices of the owned rows and columns of the PETSc.Mat object, one can use the method getOwnershipIS()

yguclu commented 5 months ago

Duplicate of #247?

e-moral-sanchez commented 5 months ago

I agree.