Re: [petsc-users] petsc4py: parallel matrix-vector multiplication
OK, thanks, I see. This is basically what happens in the poisson2d.py example, too, right? I tried it with the shell matrix (?) used in the poisson2d example and it works right away, but then I fail to see how to make use of the preconditioners for KSP (see my original message).. Thanks again! -Robert- On 06.05.18 16:52, Jed Brown wrote: > Robert Speck writes: > >> Thanks for your reply and help. Yes, this is going to be a PDE solver >> for structured grids. The first goal would be IDC (or Crank-Nicholson) >> for the heat equation, which would require both solving a linear system >> and application of the matrix. >> >> The code I wrote for testing parallel matrix-vector multiplication can >> be found here: >> https://github.com/Parallel-in-Time/pySDC/blob/petsc/pySDC/playgrounds/PETSc/playground_matmult.py >> >> Both vectors and matrix come from the DMDA, but I guess filling them is >> done in a wrong way? Or do I need to convert global/natural vectors to >> local ones somewhere? > > Global and Natural are not the same (see user's manual for details). > The matrix acts on a Global vector. See > petsc4py/demo/bratu3d/bratu3d.py for examples of efficiently setting > values (and computing residuals) using Global vectors. It should be > simpler/cleaner code than you currently have. > >> >> Best >> -Robert- >> >> >> >> On 06.05.18 14:44, Dave May wrote: >>> On Sun, 6 May 2018 at 10:40, Robert Speck >> <mailto:r.sp...@fz-juelich.de>> wrote: >>> >>> Hi! >>> >>> I would like to do a matrix-vector multiplication (besides using linear >>> solvers and so on) with petsc4py. I took the matrix from this example >>> >>> >>> (https://bitbucket.org/petsc/petsc4py/src/master/demo/kspsolve/petsc-mat.py) >>> >>> >>> This example only produces a matrix. And from the code the matrix >>> produced is identical in serial or parallel. >>> >>> >>> >>> and applied it to a PETSc Vector. All works well in serial, but in >>> parallel (in particular if ordering becomes relevant) the resulting >>> vector looks very different. >>> >>> >>> Given this, the way you defined the x vector in y = A x must be >>> different when run on 1 versus N mpi ranks. >>> >>> >>> Using the shell matrix of this example >>> >>> (https://bitbucket.org/petsc/petsc4py/src/master/demo/poisson2d/poisson2d.py) >>> helps, but then I cannot use matrix-based preconditioners for KSP >>> directly (right?). I also tried using DMDA for creating vectors and >>> matrix and for taking care of their ordering (which seems to be my >>> problem here), but that did not help either. >>> >>> So, my question is this: How do I do easy parallel matrix-vector >>> multiplication with petsc4py in a way that allows me to use parallel >>> linear solvers etc. later on? I want to deal with spatial decomposition >>> as little as possible. >>> >>> >>> What's the context - are you solving a PDE? >>> >>> 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 >>> >>> What data structures should I use? DMDA or >>> PETSc.Vec() and PETSc.Mat() or something else? >>> >>> >>> 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. >>> >>> 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. >>> >>> Thanks, >>> Dave >>> >>> >>> >>> Thanks! >>> -Robert- >>> >>> -- >>> Dr. Robert Speck >>> Juelich Supercomputing Centre >>> Institute for Advanced Simulation >>> Forschungszentrum Juelich GmbH >>> 52425 Juelich, Germany >>> >>> Tel: +49 2461 61 1644 >>> Fax: +49 2461 61
Re: [petsc-users] petsc4py: parallel matrix-vector multiplication
Thanks for your reply and help. Yes, this is going to be a PDE solver for structured grids. The first goal would be IDC (or Crank-Nicholson) for the heat equation, which would require both solving a linear system and application of the matrix. The code I wrote for testing parallel matrix-vector multiplication can be found here: https://github.com/Parallel-in-Time/pySDC/blob/petsc/pySDC/playgrounds/PETSc/playground_matmult.py Both vectors and matrix come from the DMDA, but I guess filling them is done in a wrong way? Or do I need to convert global/natural vectors to local ones somewhere? Best -Robert- On 06.05.18 14:44, Dave May wrote: > On Sun, 6 May 2018 at 10:40, Robert Speck <mailto:r.sp...@fz-juelich.de>> wrote: > > Hi! > > I would like to do a matrix-vector multiplication (besides using linear > solvers and so on) with petsc4py. I took the matrix from this example > > > (https://bitbucket.org/petsc/petsc4py/src/master/demo/kspsolve/petsc-mat.py) > > > This example only produces a matrix. And from the code the matrix > produced is identical in serial or parallel. > > > > and applied it to a PETSc Vector. All works well in serial, but in > parallel (in particular if ordering becomes relevant) the resulting > vector looks very different. > > > Given this, the way you defined the x vector in y = A x must be > different when run on 1 versus N mpi ranks. > > > Using the shell matrix of this example > > (https://bitbucket.org/petsc/petsc4py/src/master/demo/poisson2d/poisson2d.py) > helps, but then I cannot use matrix-based preconditioners for KSP > directly (right?). I also tried using DMDA for creating vectors and > matrix and for taking care of their ordering (which seems to be my > problem here), but that did not help either. > > So, my question is this: How do I do easy parallel matrix-vector > multiplication with petsc4py in a way that allows me to use parallel > linear solvers etc. later on? I want to deal with spatial decomposition > as little as possible. > > > What's the context - are you solving a PDE? > > 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 > > What data structures should I use? DMDA or > PETSc.Vec() and PETSc.Mat() or something else? > > > 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. > > 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. > > Thanks, > Dave > > > > Thanks! > -Robert- > > -- > Dr. Robert Speck > Juelich Supercomputing Centre > Institute for Advanced Simulation > Forschungszentrum Juelich GmbH > 52425 Juelich, Germany > > Tel: +49 2461 61 1644 > Fax: +49 2461 61 6656 > > Email: r.sp...@fz-juelich.de <mailto:r.sp...@fz-juelich.de> > Website: http://www.fz-juelich.de/ias/jsc/speck_r > PinT: http://www.fz-juelich.de/ias/jsc/pint > > > > > > > > Forschungszentrum Juelich GmbH > 52425 Juelich > Sitz der Gesellschaft: Juelich > Eingetragen im Handelsregister des Amtsgerichts Dueren Nr. HR B 3498 > Vorsitzender des Aufsichtsrats: MinDir Dr. Karl Eugen Huthmacher > Geschaeftsfuehrung: Prof. Dr.-Ing. Wolfgang Marquardt (Vorsitzender), > Karsten Beneke (stellv. Vorsitzender), Prof. Dr.-Ing. Harald Bolt, > Prof. Dr. Sebastian M. Schmidt > > > > > -- Dr. Robert Speck Juelich Supercomputing Centre Institute for Advanced Simulation Forschungszentrum Juelich GmbH 52425 Juelich, Germany Tel: +49 2461 61 1644 Fax: +49 2461 61 6656 Email: r.sp...@fz-juelich.de Website: http://www.fz-juelich.de/ias/jsc/speck_r PinT:http://www.fz-juelich.de/ias/jsc/pint
[petsc-users] petsc4py: parallel matrix-vector multiplication
Hi! I would like to do a matrix-vector multiplication (besides using linear solvers and so on) with petsc4py. I took the matrix from this example (https://bitbucket.org/petsc/petsc4py/src/master/demo/kspsolve/petsc-mat.py) and applied it to a PETSc Vector. All works well in serial, but in parallel (in particular if ordering becomes relevant) the resulting vector looks very different. Using the shell matrix of this example (https://bitbucket.org/petsc/petsc4py/src/master/demo/poisson2d/poisson2d.py) helps, but then I cannot use matrix-based preconditioners for KSP directly (right?). I also tried using DMDA for creating vectors and matrix and for taking care of their ordering (which seems to be my problem here), but that did not help either. So, my question is this: How do I do easy parallel matrix-vector multiplication with petsc4py in a way that allows me to use parallel linear solvers etc. later on? I want to deal with spatial decomposition as little as possible. What data structures should I use? DMDA or PETSc.Vec() and PETSc.Mat() or something else? Thanks! -Robert- -- Dr. Robert Speck Juelich Supercomputing Centre Institute for Advanced Simulation Forschungszentrum Juelich GmbH 52425 Juelich, Germany Tel: +49 2461 61 1644 Fax: +49 2461 61 6656 Email: r.sp...@fz-juelich.de Website: http://www.fz-juelich.de/ias/jsc/speck_r PinT:http://www.fz-juelich.de/ias/jsc/pint Forschungszentrum Juelich GmbH 52425 Juelich Sitz der Gesellschaft: Juelich Eingetragen im Handelsregister des Amtsgerichts Dueren Nr. HR B 3498 Vorsitzender des Aufsichtsrats: MinDir Dr. Karl Eugen Huthmacher Geschaeftsfuehrung: Prof. Dr.-Ing. Wolfgang Marquardt (Vorsitzender), Karsten Beneke (stellv. Vorsitzender), Prof. Dr.-Ing. Harald Bolt, Prof. Dr. Sebastian M. Schmidt
Re: [petsc-users] petsc4py: reuse setup for multiple solver calls?
Thank you for your answer! Please see below for comments/questions. On 05.04.18 12:53, Matthew Knepley wrote: > On Thu, Apr 5, 2018 at 6:39 AM, Robert Speck <mailto:r.sp...@fz-juelich.de>> wrote: > > Hi! > > I would like to use petsc4py for my own Python library. Installation > went well, first tests (serial and parallel) look good. > > Here is what I want to do: I have my own time-stepping playground and I > want petsc4py to be one of my backbones for the data types and (linear > or non-linear, serial or parallel) solvers. I don't want to use PETSc's > time-steppers, at least not for now. So, I saw in some examples, in > particular the ones shipped with petsc4py, that the standard way of > running one of PETSc's solvers is a bunch of setup routines, then > setting the right-hand side and solve. > > Now, I don't want to rerun the whole setup part each time I call the > solver. I know that I can change the right-hand side without having to > do so, but what if I change a parameter of my operator like, say, the > time-step size or some material parameter? > > Take the simplest case: Say I have my own implicit Euler written in > Python. I know the right-hand side F of my ODE, so in each step I want > to solve "I - dt*F". But the time-step changes every now and then, so I > cannot pre-assemble everything once and for all (or I don't want to). > What do I need to rerun before I can use the solver again, what can I > reuse? Could I just assemble F and combine it with the identity and the > parameter dt right before I call the solver? How would that look like? > > I'm pretty new to PETSc and to petsc4py, so please forgive any stupidity > or ignorance in these questions. I'm happy to take any advice, links to > examples or previous questions. Thanks! > > > For linear solves which stay the same size, you just have to call > SetOperators > again with the new operator. OK, this sounds straightforward. Thanks! > > For nonlinear solves which stay the same size, you do nothing. "nothing" in terms of "nothing you can do" or "nothing you have to do"? > > If the system size changes, it generally better to create the object. What does this mean? Thanks again! -Robert- Forschungszentrum Juelich GmbH 52425 Juelich Sitz der Gesellschaft: Juelich Eingetragen im Handelsregister des Amtsgerichts Dueren Nr. HR B 3498 Vorsitzender des Aufsichtsrats: MinDir Dr. Karl Eugen Huthmacher Geschaeftsfuehrung: Prof. Dr.-Ing. Wolfgang Marquardt (Vorsitzender), Karsten Beneke (stellv. Vorsitzender), Prof. Dr.-Ing. Harald Bolt, Prof. Dr. Sebastian M. Schmidt
[petsc-users] petsc4py: reuse setup for multiple solver calls?
Hi! I would like to use petsc4py for my own Python library. Installation went well, first tests (serial and parallel) look good. Here is what I want to do: I have my own time-stepping playground and I want petsc4py to be one of my backbones for the data types and (linear or non-linear, serial or parallel) solvers. I don't want to use PETSc's time-steppers, at least not for now. So, I saw in some examples, in particular the ones shipped with petsc4py, that the standard way of running one of PETSc's solvers is a bunch of setup routines, then setting the right-hand side and solve. Now, I don't want to rerun the whole setup part each time I call the solver. I know that I can change the right-hand side without having to do so, but what if I change a parameter of my operator like, say, the time-step size or some material parameter? Take the simplest case: Say I have my own implicit Euler written in Python. I know the right-hand side F of my ODE, so in each step I want to solve "I - dt*F". But the time-step changes every now and then, so I cannot pre-assemble everything once and for all (or I don't want to). What do I need to rerun before I can use the solver again, what can I reuse? Could I just assemble F and combine it with the identity and the parameter dt right before I call the solver? How would that look like? I'm pretty new to PETSc and to petsc4py, so please forgive any stupidity or ignorance in these questions. I'm happy to take any advice, links to examples or previous questions. Thanks! Kind regards -Robert- -- Dr. Robert Speck Juelich Supercomputing Centre Institute for Advanced Simulation Forschungszentrum Juelich GmbH 52425 Juelich, Germany Tel: +49 2461 61 1644 Fax: +49 2461 61 6656 Email: r.sp...@fz-juelich.de Website: http://www.fz-juelich.de/ias/jsc/speck_r PinT:http://www.fz-juelich.de/ias/jsc/pint Forschungszentrum Juelich GmbH 52425 Juelich Sitz der Gesellschaft: Juelich Eingetragen im Handelsregister des Amtsgerichts Dueren Nr. HR B 3498 Vorsitzender des Aufsichtsrats: MinDir Dr. Karl Eugen Huthmacher Geschaeftsfuehrung: Prof. Dr.-Ing. Wolfgang Marquardt (Vorsitzender), Karsten Beneke (stellv. Vorsitzender), Prof. Dr.-Ing. Harald Bolt, Prof. Dr. Sebastian M. Schmidt
