Thank your for this explanation, it makes sense. And after I updated my
code, the external CFD code runs without problems in parallel.
However, now I'm back to the problem with the creation of the
vectors/domains. By using the application ordering, I can assign the
correct points from PETSc to the corresponding points in my external
code. At least, as long as both use the same subdomain size. But
sometimes they differ, and then the KSP breaks down, because the
solution Vector it receives has a different size than what it expects.
An example:
I have an unstructured grid with 800,000 datapoints.
If I decompose this to run on 2 processors, PETSc delegates exactly
400,000 points to each process. However, the external code might assign
400,100 points to the first and 399,900 process. As a result, PETSc
expects a solution vector of size 400,000 on each process, but receives
one of 400,100 and one of 399,900, leading to a breakdown.
I suppose I could use VecScatterCreateToAll to collect all the values
from the solution vectors of my external code, and then create from
those a temporary vector that only contains the expected 400,000 values
to hand over to the KSP. But this would create a lot of communication
between the different processes and seems quite clunky.
Is there a more elegant way? Is it maybe possible to manually assign the
size of the PETSc subdomains?
Kind regards,
Michael Werner
Am 17.10.2017 um 12:31 schrieb Matthew Knepley:
On Tue, Oct 17, 2017 at 6:08 AM, Michael Werner <[email protected]
<mailto:[email protected]>> wrote:
Because usally this code is called just once. It runs one multiple
processes, but there it's still always processing the whole
domain. I can't run it on only one subdomain. As I understand it
now, when I call it from PETSc, this call is issued once per
process, so I would end up running several contesting instances of
the computation on the whole domain.
But maybe that's only because I haven't completly understood how
MPI really works in such cases...
No, it makes one call in which all processes participate. So you would
call your external CFD routine once from all processes, passing in the
MPI communicator.
Matt
Kind regards,
Michael
Am 17.10.2017 um 11:50 schrieb Matthew Knepley:
On Tue, Oct 17, 2017 at 5:46 AM, Michael Werner
<[email protected] <mailto:[email protected]>> wrote:
I'm not sure what you mean with this question?
The external CFD code, if that was what you referred to, can
be run in parallel.
Then why is it a problem that "in a parallel case, this call
obviously gets called once per process"?
Matt
Am 17.10.2017 um 11:11 schrieb Matthew Knepley:
On Tue, Oct 17, 2017 at 4:21 AM, Michael Werner
<[email protected] <mailto:[email protected]>> wrote:
That's something I'm still struggling with. In the
serial case, I can simply extract the values from the
original grid, and since the ordering of the Jacobian is
the same there is no problem. In the parallel case this
is still a more or less open question. That's why I
thought about reordering the Jacobian. As long as the
position of the individual IDs is the same for both, I
don't have to care about their absolute position.
I also wanted to thank you for your previous answer, it
seems that the application ordering might be what I'm
looking for. However, in the meantime I stumbled about
another problem, that I have to solve first. My new
problem is, that I call the external code within the
shell matrix' multiply call. But in a parallel case,
this call obviously gets called once per process. So
right now I'm trying to circumvent this, so it might
take a while before I'm able to come back to the
original problem...
I am not understanding. Is your original code parallel?
Thanks,
Matt
Kind regards,
Michael
Am 16.10.2017 um 17:25 schrieb Praveen C:
I am interested to learn more about how this works. How
are the vectors created if the ids are not contiguous
in a partition ?
Thanks
praveen
On Mon, Oct 16, 2017 at 2:02 PM, Stefano Zampini
<[email protected]
<mailto:[email protected]>> wrote:
2017-10-16 10:26 GMT+03:00 Michael Werner
<[email protected] <mailto:[email protected]>>:
Hello,
I'm having trouble with parallelizing a
matrix-free code with PETSc. In this code, I
use an external CFD code to provide the
matrix-vector product for an iterative solver
in PETSc. To increase convergence rate, I'm
using an explicitly stored Jacobian matrix to
precondition the solver. This works fine for
serial runs. However, when I try to use
multiple processes, I face the problem that
PETSc decomposes the preconditioner matrix, and
probably also the shell matrix, in a different
way than the external CFD code decomposes the grid.
The Jacobian matrix is built in a way, that its
rows and columns correspond to the global IDs
of the individual points in my CFD mesh
The CFD code decomposes the domain based on the
proximity of points to each other, so that the
resulting subgrids are coherent. However, since
its an unstructured grid, those subgrids are
not necessarily made up of points with
successive global IDs. This is a problem, since
PETSc seems to partition the matrix in coherent
slices.
I'm not sure what the best approach to this
problem might be. Is it maybe possible to
exactly tell PETSc, which rows/columns it
should assign to the individual processes?
If you are explicitly setting the values in your
Jacobians via MatSetValues(), you can create a
ISLocalToGlobalMapping
http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/IS/ISLocalToGlobalMappingCreate.html
<http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/IS/ISLocalToGlobalMappingCreate.html>
that maps the numbering you use for the Jacobians
to their counterpart in the CFD ordering, then call
MatSetLocalToGlobalMapping and then use
MatSetValuesLocal with the same arguments you are
calling MatSetValues now.
Otherwise, you can play with the application
ordering
http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/AO/index.html
<http://www.mcs.anl.gov/petsc/petsc-current/docs/manualpages/AO/index.html>
--
Stefano
--
____________________________________________________
Deutsches Zentrum für Luft- und Raumfahrt e.V. (DLR)
Institut für Aerodynamik und Strömungstechnik |Bunsenstr. 10 |
37073 Göttingen
<https://maps.google.com/?q=Bunsenstr.+10+%7C+37073+G%C3%B6ttingen&entry=gmail&source=g>
Michael Werner
Telefon 0551 709-2627 | Telefax 0551 709-2811 |[email protected]
<mailto:[email protected]>
DLR.de
--
What most experimenters take for granted before they begin
their experiments is infinitely more interesting than any
results to which their experiments lead.
-- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/
<http://www.caam.rice.edu/%7Emk51/>
--
____________________________________________________
Deutsches Zentrum für Luft- und Raumfahrt e.V. (DLR)
Institut für Aerodynamik und Strömungstechnik |Bunsenstr. 10 | 37073
Göttingen
<https://maps.google.com/?q=Bunsenstr.+10+%7C+37073+G%C3%B6ttingen&entry=gmail&source=g>
Michael Werner
Telefon 0551 709-2627 | Telefax 0551 709-2811 |[email protected]
<mailto:[email protected]>
DLR.de
--
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to
which their experiments lead.
-- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/
<http://www.caam.rice.edu/%7Emk51/>
--
____________________________________________________
Deutsches Zentrum für Luft- und Raumfahrt e.V. (DLR)
Institut für Aerodynamik und Strömungstechnik |Bunsenstr. 10 | 37073
Göttingen
<https://maps.google.com/?q=Bunsenstr.+10+%7C+37073+G%C3%B6ttingen&entry=gmail&source=g>
Michael Werner
Telefon 0551 709-2627 | Telefax 0551 709-2811 |[email protected]
<mailto:[email protected]>
DLR.de
--
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which
their experiments lead.
-- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.caam.rice.edu/%7Emk51/>
