November 2011 - petsc-dev - mailman3.cels.anl.gov

Re: [petsc-dev] Seeking Advice on Petsc Preconditioners to Try
by Matthew Knepley 30 Nov '11

30 Nov '11

On Wed, Nov 30, 2011 at 12:41 AM, Dave Nystrom <dnystrom1(a)comcast.net>wrote: > I have a linear system in a code that I have interfaced to petsc that is > taking about 80 percent of the run time per timestep. This linear system > is > a symmetric block banded matrix where the blocks are 2x2. The matrix looks > as follows: > > 1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 0 > 1X X Y Y Y > 2X X X Y Y Y > 3 X X X Y Y Y > 4 X X X Y Y Y > 5 X X X Y Y Y > 6 X X X Y Y Y > 7 X X X Y Y Y > 8 X X X Y Y Y > 9 X X X Y Y Y > 0 X X X Y Y Y > 1 X X X Y Y Y > 2 X X X Y Y Y > 3Z X X X Y Y Y > 4Z Z X X X Y Y Y > 5Z Z Z X X X Y Y Y > 6 Z Z Z X X X Y Y Y > 7 Z Z Z X X X Y Y Y > 8 Z Z Z X X X Y Y Y > 9 Z Z Z X X X Y Y Y > 0 Z Z Z X X X Y Y Y > > So in my diagram above, X, Y and Z are 2x2 blocks. The symmetry of the > matrix requires that X_ij = transpose(X_ji) and Y_ij = transpose(Z_ji). So > far, I have just input this matrix to petsc without indicating that it was > block banded with 2x2 blocks. I have also not told petsc that the matrix > is > symmetric. And I have allowed petsc to decide the best way to store the > matrix. > > I can solve this linear system over the course of a run using -ksp_type > preonly -pc_type lu. But that will not scale very well to larger problems > that I want to solve. I can also solve this system over the course of a > run > using -ksp_type cg -pc_type jacobi -vec_type cusp -mat_type aijcusp. > However, over the course of a run, the iteration count ranges from 771 to > 47300. I have also tried sacusp, ainvcusp, sacusppoly, ilu(k) and icc(k) > with k=0. The sacusppoly preconditioner fails because of a thrust error > related to an invalid device pointer, if I am remembering correctly. I > reported this problem to petsc-maint a while back and have also reported it > for the cusp bugtracker. But it does not appear that anyone has really > looked into the bug. For the other preconditioners of sacusp, ilu(k) and > icc(k), they do not result in convergence to a solution and the runs fail. > All preconditioners are custom. Have you done a literature search for PCs known to work for this problem? Can yu say anything about the spectrum of the operator? conditioning? what is the principal symbol (if its a PDE)? The pattern is not enough to recommend a PC. Matt > I'm wondering if there are suggestions of other preconditioners in petsc > that > I should try. The only third party package that I have tried is the > txpetscgpu package. I have not tried hypre or any of the multigrid > preconditioners yet. I'm not sure how difficult it is to try those > packages. Anyway, so far I have not found a preconditioner available in > petsc that provides a robust solution to this problem and would be > interested > in any suggestions that anyone might have of things to try. > > I'd be happy to provide additional info and am planning on packaging up a > couple of examples of the matrix and rhs for people I am interacting with > at > Tech-X and EMPhotonics. So I'd be happy to provide the matrix examples for > this forum as well if anyone wants a copy. > > Thanks, > > Dave > > -- > Dave Nystrom > > phone: 505-661-9943 (home office) > 505-662-6893 (home) > skype: dave.nystrom76 > email: dnystrom1(a)comcast.net > smail: 219 Loma del Escolar > Los Alamos, NM 87544 > -- 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

1 0

Re: [petsc-dev] Performance of Petsc Cholesky vs LU
by Dave Nystrom 30 Nov '11

30 Nov '11

Jed Brown writes: > On Tue, Nov 29, 2011 at 23:53, Dave Nystrom <dnystrom1(a)comcast.net> wrote: > > > I have a resistive mhd code that I have recently interfaced to petsc which > > has 7 linear solves that are all symmetric. I recently tried using > > -pc_type > > cholesky -ksp_type preonly for a run and found that it was taking about 6 > > times as long per linear solve as when I was using -pc_type lu -ksp_type > > preonly. > > > > Try -pc_factor_mat_ordering_type nd > > > > I was wondering if that was reasonable behavior. I would not have > > thought that using a cholesky direct solve would take longer than an LU > > direct solve in petsc for the serial case and was hoping it would be > > faster. > > Does this behavior seem reasonable? > > > > Try this: > > -pc_type cholesky -pc_factor_mat_ordering_type nd Thanks. I'll give this a try and report back on the results. > Barry, why is natural ordering still the default for Cholesky? It is so > slow that it is worthless.

1 0

Seeking Advice on Petsc Preconditioners to Try
by Dave Nystrom 30 Nov '11

30 Nov '11

I have a linear system in a code that I have interfaced to petsc that is taking about 80 percent of the run time per timestep. This linear system is a symmetric block banded matrix where the blocks are 2x2. The matrix looks as follows: 1 2 3 4 5 6 7 8 9 0 1 2 3 4 5 6 7 8 9 0 1X X Y Y Y 2X X X Y Y Y 3 X X X Y Y Y 4 X X X Y Y Y 5 X X X Y Y Y 6 X X X Y Y Y 7 X X X Y Y Y 8 X X X Y Y Y 9 X X X Y Y Y 0 X X X Y Y Y 1 X X X Y Y Y 2 X X X Y Y Y 3Z X X X Y Y Y 4Z Z X X X Y Y Y 5Z Z Z X X X Y Y Y 6 Z Z Z X X X Y Y Y 7 Z Z Z X X X Y Y Y 8 Z Z Z X X X Y Y Y 9 Z Z Z X X X Y Y Y 0 Z Z Z X X X Y Y Y So in my diagram above, X, Y and Z are 2x2 blocks. The symmetry of the matrix requires that X_ij = transpose(X_ji) and Y_ij = transpose(Z_ji). So far, I have just input this matrix to petsc without indicating that it was block banded with 2x2 blocks. I have also not told petsc that the matrix is symmetric. And I have allowed petsc to decide the best way to store the matrix. I can solve this linear system over the course of a run using -ksp_type preonly -pc_type lu. But that will not scale very well to larger problems that I want to solve. I can also solve this system over the course of a run using -ksp_type cg -pc_type jacobi -vec_type cusp -mat_type aijcusp. However, over the course of a run, the iteration count ranges from 771 to 47300. I have also tried sacusp, ainvcusp, sacusppoly, ilu(k) and icc(k) with k=0. The sacusppoly preconditioner fails because of a thrust error related to an invalid device pointer, if I am remembering correctly. I reported this problem to petsc-maint a while back and have also reported it for the cusp bugtracker. But it does not appear that anyone has really looked into the bug. For the other preconditioners of sacusp, ilu(k) and icc(k), they do not result in convergence to a solution and the runs fail. I'm wondering if there are suggestions of other preconditioners in petsc that I should try. The only third party package that I have tried is the txpetscgpu package. I have not tried hypre or any of the multigrid preconditioners yet. I'm not sure how difficult it is to try those packages. Anyway, so far I have not found a preconditioner available in petsc that provides a robust solution to this problem and would be interested in any suggestions that anyone might have of things to try. I'd be happy to provide additional info and am planning on packaging up a couple of examples of the matrix and rhs for people I am interacting with at Tech-X and EMPhotonics. So I'd be happy to provide the matrix examples for this forum as well if anyone wants a copy. Thanks, Dave -- Dave Nystrom phone: 505-661-9943 (home office) 505-662-6893 (home) skype: dave.nystrom76 email: dnystrom1(a)comcast.net smail: 219 Loma del Escolar Los Alamos, NM 87544

1 0

Re: [petsc-dev] Performance of Petsc Cholesky vs LU
by Jed Brown 30 Nov '11

30 Nov '11

On Tue, Nov 29, 2011 at 23:53, Dave Nystrom <dnystrom1(a)comcast.net> wrote: > I have a resistive mhd code that I have recently interfaced to petsc which > has 7 linear solves that are all symmetric. I recently tried using > -pc_type > cholesky -ksp_type preonly for a run and found that it was taking about 6 > times as long per linear solve as when I was using -pc_type lu -ksp_type > preonly. > Try -pc_factor_mat_ordering_type nd > I was wondering if that was reasonable behavior. I would not have > thought that using a cholesky direct solve would take longer than an LU > direct solve in petsc for the serial case and was hoping it would be > faster. > Does this behavior seem reasonable? > Try this: -pc_type cholesky -pc_factor_mat_ordering_type nd Barry, why is natural ordering still the default for Cholesky? It is so slow that it is worthless.

1 0

Performance of Petsc Cholesky vs LU
by Dave Nystrom 29 Nov '11

29 Nov '11

I have a resistive mhd code that I have recently interfaced to petsc which has 7 linear solves that are all symmetric. I recently tried using -pc_type cholesky -ksp_type preonly for a run and found that it was taking about 6 times as long per linear solve as when I was using -pc_type lu -ksp_type preonly. I was wondering if that was reasonable behavior. I would not have thought that using a cholesky direct solve would take longer than an LU direct solve in petsc for the serial case and was hoping it would be faster. Does this behavior seem reasonable? Thanks, Dave -- Dave Nystrom phone: 505-661-9943 (home office) 505-662-6893 (home) skype: dave.nystrom76 email: dnystrom1(a)comcast.net smail: 219 Loma del Escolar Los Alamos, NM 87544

1 0

Petsc Pthreads Capability
by Dave Nystrom 29 Nov '11

29 Nov '11

I have a 2d resistive mhd code which spends over 95 percent of its run time doing linear solves and I have interfaced these linear solves to petsc. I'm interested in optimizing the single node performance of this code right now and was wondering if there could be any benefit to trying out the petsc pthreads capability that is in petsc-dev. Is this work at a stage yet where I might be able to benefit from trying it with my code when running on a multi-core cpu? Thanks, Dave -- Dave Nystrom phone: 505-661-9943 (home office) 505-662-6893 (home) skype: dave.nystrom76 email: dnystrom1(a)comcast.net smail: 219 Loma del Escolar Los Alamos, NM 87544

1 0

GPU based preconditioner and PCSHELL
by Dave Nystrom 29 Nov '11

29 Nov '11

I have a 2d resistive mhd code interfaced to petsc. The code has seven different linear solves per timestep and these linear solves consume around 95 percent of the run time for a reasonably small grid of 100x301. The run time is dominated by the solution of one of the linear solves that is particularly difficult to solve. This particular linear solve actually takes about 80 percent of the total run time. All of these linear systems are symmetric. I am actually able to run a simulation where I am using jacobi preconditioning with conjugate gradient on a gpu and get a solution and this is the fastest solution which I currently get but the iteration count ranges from a minimum of 771 to a max of 47300 for the difficult linear system. We also have a native cg solver in the code which uses a full cholesky band solve of the inner set of bands and that solver has a much better iteration count although it takes a bit over 2x the run time of petsc w/ jacobi and cusp. So, I have been interested in making a custom petsc preconditioner that does this cholesky solve using the PCSHELL capability of petsc. And I'm interested in trying to do this cholesky solve on the GPU using a CULA SPARSE band solver. However, I'm wondering if the petsc PCSHELL capability only runs on the cpu and if I would need a GPU analog such as a PCGPUSHELL capability in petsc. I'm wondering if this is the case and if so, whether there is a possibility that a PCGPUSHELL capability might be added to petsc in the near term. Thanks, Dave -- Dave Nystrom phone: 505-661-9943 (home office) 505-662-6893 (home) skype: dave.nystrom76 email: dnystrom1(a)comcast.net smail: 219 Loma del Escolar Los Alamos, NM 87544

1 0

test please ignore
by Satish Balay 29 Nov '11

29 Nov '11

1 0

Re: [petsc-dev] Generality of VecScatter
by Jed Brown 29 Nov '11

29 Nov '11

On Tue, Nov 29, 2011 at 14:09, Dmitry Karpeev <karpeev(a)mcs.anl.gov> wrote: > I wouldn't consider that stackoverflow post to be authoritative, but there >> is a large body of literature on lock-free synchronization. >> > > Could you point me to some of it? > These are not bad places to start https://www.google.com/search?q=lock-free http://scholar.google.com/scholar?q=lock-free > >> Unless there is resource contention, spinning is always the lowest >> latency choice. All multi-threaded cores have some sort of pause >> instruction because it's necessary to make that design useful. For this >> purpose, we should just spin. >> >> Note that with sequential consistency, >> > > What do you mean by "sequential consistency"? If we did serialize the > threads before waking them up? > Sequential consistency has to do with the memory model. http://en.wikipedia.org/wiki/Sequential_consistency

1 0

Re: [petsc-dev] Generality of VecScatter
by Dmitry Karpeev 29 Nov '11

29 Nov '11

On Tue, Nov 29, 2011 at 1:01 PM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > On Tue, Nov 29, 2011 at 12:53, Dmitry Karpeev <karpeev(a)mcs.anl.gov> wrote: > >> I understand that any solution requires serialization: atomic reads by >> the spinning loop would serialize on that read, >> while pthread_cond_wait requires mutex serialization. sigwait serializes >> since it clears the signal atomically. >> Ultimately, all three rely on atomicity of some sort, but mutexes >> apparently have higher overhead for it? >> At least the stackoverflow page that suggests the sigwait solution >> reports a 40x improvement over the pthread_cond_wait >> solution (admittedly, I don't know the details of the sigwait stuff). >> > > I wouldn't consider that stackoverflow post to be authoritative, but there > is a large body of literature on lock-free synchronization. > Could you point me to some of it? > Unless there is resource contention, spinning is always the lowest latency > choice. All multi-threaded cores have some sort of pause instruction > because it's necessary to make that design useful. For this purpose, we > should just spin. > > Note that with sequential consistency, > What do you mean by "sequential consistency"? If we did serialize the threads before waking them up? > you can also provide information about the next work unit when you > activate a thread. For example, you can pass a function pointer. >

1 0