<div><div><div class="gmail_quote"><div dir="auto">On Sun, 6 May 2018 at 10:40, Robert Speck <<a href="mailto:r.speck@fz-juelich.de" target="_blank">r.speck@fz-juelich.de</a>> wrote:<br></div></div></div><div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Hi!<br>
<br>
I would like to do a matrix-vector multiplication (besides using linear<br>
solvers and so on) with petsc4py. I took the matrix from this example<br></blockquote></div></div></div><div><div><div class="gmail_quote"><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
(<a href="https://bitbucket.org/petsc/petsc4py/src/master/demo/kspsolve/petsc-mat.py" rel="noreferrer" target="_blank">https://bitbucket.org/petsc/petsc4py/src/master/demo/kspsolve/petsc-mat.py</a>)</blockquote><div dir="auto"><br></div><div dir="auto">This example only produces a matrix. And from the code the matrix produced is identical in serial or parallel. </div></div></div></div><div><div><div class="gmail_quote"><div dir="auto"><br></div><div dir="auto"><br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><br>
and applied it to a PETSc Vector. All works well in serial, but in<br>
parallel (in particular if ordering becomes relevant) the resulting<br>
vector looks very different. </blockquote><div dir="auto"><br></div></div></div></div><div><div><div class="gmail_quote"><div dir="auto">Given this, the way you defined the x vector in y = A x must be different when run on 1 versus N mpi ranks. </div></div></div></div><div><div><div class="gmail_quote"><div dir="auto"><br></div><div dir="auto"><br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">Using the shell matrix of this example<br>
(<a href="https://bitbucket.org/petsc/petsc4py/src/master/demo/poisson2d/poisson2d.py" rel="noreferrer" target="_blank">https://bitbucket.org/petsc/petsc4py/src/master/demo/poisson2d/poisson2d.py</a>)<br>
helps, but then I cannot use matrix-based preconditioners for KSP<br>
directly (right?). I also tried using DMDA for creating vectors and<br>
matrix and for taking care of their ordering (which seems to be my<br>
problem here), but that did not help either.<br>
<br>
So, my question is this: How do I do easy parallel matrix-vector<br>
multiplication with petsc4py in a way that allows me to use parallel<br>
linear solvers etc. later on? I want to deal with spatial decomposition<br>
as little as possible. </blockquote><div dir="auto"><br></div><div dir="auto">What's the context - are you solving a PDE?</div><div dir="auto"><br></div><div dir="auto">Assuming you are using your own grid object (e.g. as you might have if solving a PDE), and assuming you are not solving a 1D problem, you actually have to "deal" with the spatial decomposition otherwise performance could be quite terrible - even for something simple like a 5 point Laplacian on a structured grid in 2D</div><div dir="auto"><br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">What data structures should I use? DMDA or<br>
PETSc.Vec() and PETSc.Mat() or something else?</blockquote><div dir="auto"><br></div></div></div></div><div><div><div class="gmail_quote"><div dir="auto">The mat vec product is not causing you a problem. Your issue appears to be that you do not have a way to label entries in a vector in a consistent manner.</div><div dir="auto"><br></div><div dir="auto">What's the objective? Are you solving a PDE? If yes, structured grid? If yes again, use the DMDA. It takes care of all the local-to-global and global-to-local mapping you need.</div><div dir="auto"><br></div><div dir="auto">Thanks,</div><div dir="auto"> Dave</div></div></div></div><div><div><div class="gmail_quote"><div dir="auto"><br></div><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex"><br>
<br>
Thanks!<br>
-Robert-<br>
<br>
--<br>
Dr. Robert Speck<br>
Juelich Supercomputing Centre<br>
Institute for Advanced Simulation<br>
Forschungszentrum Juelich GmbH<br>
52425 Juelich, Germany<br>
<br>
Tel: +49 2461 61 1644<br>
Fax: +49 2461 61 6656<br>
<br>
Email: <a href="mailto:r.speck@fz-juelich.de" target="_blank">r.speck@fz-juelich.de</a><br>
Website: <a href="http://www.fz-juelich.de/ias/jsc/speck_r" rel="noreferrer" target="_blank">http://www.fz-juelich.de/ias/jsc/speck_r</a><br>
PinT: <a href="http://www.fz-juelich.de/ias/jsc/pint" rel="noreferrer" target="_blank">http://www.fz-juelich.de/ias/jsc/pint</a><br>
<br>
<br>
<br>
------------------------------------------------------------------------------------------------<br>
------------------------------------------------------------------------------------------------<br>
Forschungszentrum Juelich GmbH<br>
52425 Juelich<br>
Sitz der Gesellschaft: Juelich<br>
Eingetragen im Handelsregister des Amtsgerichts Dueren Nr. HR B 3498<br>
Vorsitzender des Aufsichtsrats: MinDir Dr. Karl Eugen Huthmacher<br>
Geschaeftsfuehrung: Prof. Dr.-Ing. Wolfgang Marquardt (Vorsitzender),<br>
Karsten Beneke (stellv. Vorsitzender), Prof. Dr.-Ing. Harald Bolt,<br>
Prof. Dr. Sebastian M. Schmidt<br>
------------------------------------------------------------------------------------------------<br>
------------------------------------------------------------------------------------------------<br>
<br>
</blockquote></div></div></div>