May 2012 - petsc-dev - mailman3.cels.anl.gov

Re: [petsc-dev] anyone feel like fixing a bug?
by Barry Smith 23 May '12

23 May '12

On May 23, 2012, at 9:26 AM, Jed Brown wrote: > On Wed, May 23, 2012 at 9:27 AM, Mark F. Adams <mark.adams(a)columbia.edu> wrote: > I have gone through and done this manual inheritance, but 1) I may have missed some and 2) there were some that I could not figure out how to do. I'm back online today and can work on this. I did leave some comments in the form of SetBlockSize statements that are commented out in places were I had trouble. > > The problem is fixed. Thanks Jed. Just goes to show with the right model the solution is always easy :-) Barry

1 0

Re: [petsc-dev] Multigrid is confusing
by Jed Brown 23 May '12

23 May '12

If ex19 is what matters... $ base='mpiexec.hydra -n 4 ./ex19 -da_refine 3 -lidvelocity 100 -grashof 1e4 -pc_type mg -ksp_type fgmres -pc_mg_type full -ksp_converged_reason' $ $base lid velocity = 100, prandtl # = 1, grashof # = 10000 Linear solve converged due to CONVERGED_RTOL iterations 23 Linear solve converged due to CONVERGED_RTOL iterations 12 Linear solve converged due to CONVERGED_RTOL iterations 12 Linear solve converged due to CONVERGED_RTOL iterations 17 Linear solve converged due to CONVERGED_RTOL iterations 15 Linear solve converged due to CONVERGED_RTOL iterations 13 Linear solve converged due to CONVERGED_RTOL iterations 19 Linear solve converged due to CONVERGED_RTOL iterations 21 Linear solve converged due to CONVERGED_RTOL iterations 18 Number of SNES iterations = 9 $ $base -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_estimate_eigenvalues 0.1,1.1 -mg_levels_ksp_max_it 1 -mg_levels_pc_type pbjacobi lid velocity = 100, prandtl # = 1, grashof # = 10000 Linear solve converged due to CONVERGED_RTOL iterations 25 Linear solve converged due to CONVERGED_RTOL iterations 27 Linear solve converged due to CONVERGED_RTOL iterations 27 Linear solve converged due to CONVERGED_RTOL iterations 30 Linear solve converged due to CONVERGED_RTOL iterations 28 Linear solve converged due to CONVERGED_RTOL iterations 42 Linear solve converged due to CONVERGED_RTOL iterations 49 Linear solve converged due to CONVERGED_RTOL iterations 52 Linear solve converged due to CONVERGED_RTOL iterations 48 Number of SNES iterations = 9 $ $base -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_estimate_eigenvalues 0.1,1.1 -mg_levels_ksp_max_it 1 -mg_levels_pc_type sor lid velocity = 100, prandtl # = 1, grashof # = 10000 Linear solve converged due to CONVERGED_RTOL iterations 22 Linear solve converged due to CONVERGED_RTOL iterations 16 Linear solve converged due to CONVERGED_RTOL iterations 17 Linear solve converged due to CONVERGED_RTOL iterations 19 Linear solve converged due to CONVERGED_RTOL iterations 19 Linear solve converged due to CONVERGED_RTOL iterations 17 Linear solve converged due to CONVERGED_RTOL iterations 20 Linear solve converged due to CONVERGED_RTOL iterations 18 Linear solve converged due to CONVERGED_RTOL iterations 17 Number of SNES iterations = 9 $ $base -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_estimate_eigenvalues 0.1,1.1 -mg_levels_ksp_max_it 2 -mg_levels_pc_type pbjacobi lid velocity = 100, prandtl # = 1, grashof # = 10000 Linear solve converged due to CONVERGED_RTOL iterations 15 Linear solve converged due to CONVERGED_RTOL iterations 13 Linear solve converged due to CONVERGED_RTOL iterations 16 Linear solve converged due to CONVERGED_RTOL iterations 15 Linear solve converged due to CONVERGED_RTOL iterations 16 Linear solve converged due to CONVERGED_RTOL iterations 16 Linear solve converged due to CONVERGED_RTOL iterations 18 Linear solve converged due to CONVERGED_RTOL iterations 17 Linear solve converged due to CONVERGED_RTOL iterations 16 Number of SNES iterations = 9 $ $base -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_estimate_eigenvalues 0.1,1.1 -mg_levels_ksp_max_it 2 -mg_levels_pc_type sor lid velocity = 100, prandtl # = 1, grashof # = 10000 Linear solve converged due to CONVERGED_RTOL iterations 10 Linear solve converged due to CONVERGED_RTOL iterations 7 Linear solve converged due to CONVERGED_RTOL iterations 8 Linear solve converged due to CONVERGED_RTOL iterations 8 Linear solve converged due to CONVERGED_RTOL iterations 9 Linear solve converged due to CONVERGED_RTOL iterations 9 Linear solve converged due to CONVERGED_RTOL iterations 8 Linear solve converged due to CONVERGED_RTOL iterations 8 Linear solve converged due to CONVERGED_RTOL iterations 7 Number of SNES iterations = 9 $ $base -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_estimate_eigenvalues 0.1,1.1 -mg_levels_ksp_max_it 3 -mg_levels_pc_type pbjacobi lid velocity = 100, prandtl # = 1, grashof # = 10000 Linear solve converged due to CONVERGED_RTOL iterations 12 Linear solve converged due to CONVERGED_RTOL iterations 9 Linear solve converged due to CONVERGED_RTOL iterations 11 Linear solve converged due to CONVERGED_RTOL iterations 11 Linear solve converged due to CONVERGED_RTOL iterations 11 Linear solve converged due to CONVERGED_RTOL iterations 12 Linear solve converged due to CONVERGED_RTOL iterations 12 Linear solve converged due to CONVERGED_RTOL iterations 11 Linear solve converged due to CONVERGED_RTOL iterations 11 Number of SNES iterations = 9 $ $base -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_estimate_eigenvalues 0.1,1.1 -mg_levels_ksp_max_it 3 -mg_levels_pc_type sor lid velocity = 100, prandtl # = 1, grashof # = 10000 Linear solve converged due to CONVERGED_RTOL iterations 9 Linear solve converged due to CONVERGED_RTOL iterations 6 Linear solve converged due to CONVERGED_RTOL iterations 6 Linear solve converged due to CONVERGED_RTOL iterations 7 Linear solve converged due to CONVERGED_RTOL iterations 7 Linear solve converged due to CONVERGED_RTOL iterations 6 Linear solve converged due to CONVERGED_RTOL iterations 7 Linear solve converged due to CONVERGED_RTOL iterations 6 Linear solve converged due to CONVERGED_RTOL iterations 6 Number of SNES iterations = 9

1 0

Re: [petsc-dev] anyone feel like fixing a bug?
by Mark F. Adams 23 May '12

23 May '12

On May 23, 2012, at 9:36 AM, Barry Smith wrote: > > On May 22, 2012, at 9:32 PM, Jed Brown wrote: > >> Well, here's the problem: >> >> MatSetBlockSize() used to call mat->ops->setblocksize() which the MPI matrices forwarded down to their seq submatrices. >> >> But then this patch removed that behavior so MatSetBlockSize() is processed immediately and not forwarded down. >> >> http://petsc.cs.iit.edu/petsc/petsc-dev/rev/80ec81eb9db9 >> >> Since the Seq inner matrices have different PetscLayouts, they don't inherit the block sizes set on their parent. >> >> Barry, why did you remove this MatOp? How do you want the MPI matrices to forward block sizes now? > > Blocksize is now suppose to be known "early" now. That is, you cannot change it after MatSetUp() is called. It is the job of each matrix implementation to use the already set block size in the appropriate way. Generally this means before the "Seq" matrices are created or finalized. One could possible arrange moving the blocksize information down during the MatSetUp() of that implementation or the MatXXXXSetPreallocation() or when the inner Mat objects are created. I have gone through and done this manual inheritance, but 1) I may have missed some and 2) there were some that I could not figure out how to do. I'm back online today and can work on this. I did leave some comments in the form of SetBlockSize statements that are commented out in places were I had trouble. Mark > > If you tell me the EXACT change of events of the "inner" block size not being set properly I can take a look at it. > > Barry > >> >> On Sat, May 19, 2012 at 10:57 AM, Mark F. Adams <mark.adams(a)columbia.edu> wrote: >> I'm going to be tied up for a few days, should be able to get to this Wed., if not sooner, but if anyone wants to take a crack at it, here it is. This is run with 8 procs. >> >> Mark >> >> ==51037== Command: ./ex56 -ne 19 -alpha 1.e-3 -ksp_type cg -ksp_monitor -ksp_max_it 15 -pc_type gamg -pc_gamg_type agg -pc_gamg_agg_nsmooths 1 -pc_gamg_verbose 2 -ksp_converged_reason -options_left -theta -30.0 -epsilon .01 -blob_center 0.,0. -pc_gamg_threshold .0 -pc_gamg_coarse_eq_limit 100 -use_coordinates -out_verbose 2 -aggmg_smooths 1 -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_estimate_eigenvalues 0,0.2,0,1.1 -mg_levels_pc_type pbjacobi -mg_levels_ksp_max_it 1 >> >> ==51039== Invalid read of size 8 >> ==51039== at 0x91E38B: PCApply_PBJacobi_3 (pbjacobi.c:81) >> ==51039== by 0xA0FE23: PCApply (precon.c:384) >> ==51039== by 0xADA9B2: KSPInitialResidual (itres.c:64) >> ==51039== by 0xA7EA39: KSPSolve_GMRES (gmres.c:230) >> ==51039== by 0xAB3CA2: KSPSolve (itfunc.c:446) >> ==51039== by 0xAA6FF8: KSPSolve_Chebyshev (cheby.c:269) >> ==51039== by 0xAB3CA2: KSPSolve (itfunc.c:446) >> ==51039== by 0x9ECD2F: PCMGMCycle_Private (mg.c:20) >> ==51039== by 0x9F0BEF: PCApply_MG (mg.c:321) >> ==51039== by 0xA0FE23: PCApply (precon.c:384) >> ==51039== by 0xA43EBF: KSPSolve_CG (cg.c:139) >> ==51039== by 0xAB3CA2: KSPSolve (itfunc.c:446) >> ==51039== by 0x4805: main (in ./ex56) >> ==51039== Address 0x27f4c90 is 0 bytes after a block of size 24,032 alloc'd >> ==51039== at 0x1B2AF16: malloc (vg_replace_malloc.c:236) >> ==51039== by 0x8A8C7: PetscMallocAlign (mal.c:37) >> ==51039== by 0x61E602: MatInvertBlockDiagonal_SeqAIJ (aij.c:2816) >> ==51039== by 0x2F0AC3: MatInvertBlockDiagonal (matrix.c:9280) >> ==51039== by 0x6D1494: MatInvertBlockDiagonal_MPIAIJ (mpiaij.c:2981) >> ==51039== by 0x2F0AC3: MatInvertBlockDiagonal (matrix.c:9280) >> ==51039== by 0x920A85: PCSetUp_PBJacobi (pbjacobi.c:213) >> ==51039== by 0xA1783B: PCSetUp (precon.c:832) >> ==51039== by 0xAB2125: KSPSetUp (itfunc.c:278) >> ==51039== by 0x9F536C: PCSetUp_MG (mg.c:678) >> ==51039== by 0x98C26C: PCSetUp_GAMG (gamg.c:984) >> ==51039== by 0xA1783B: PCSetUp (precon.c:832) >> ==51039== by 0xAB2125: KSPSetUp (itfunc.c:278) >> ==51039== by 0x46AD: main (in ./ex56) >> ==51039== >> >> > >

1 0

Re: [petsc-dev] anyone feel like fixing a bug?
by Jed Brown 23 May '12

23 May '12

On Wed, May 23, 2012 at 9:27 AM, Mark F. Adams <mark.adams(a)columbia.edu>wrote: > I have gone through and done this manual inheritance, but 1) I may have > missed some and 2) there were some that I could not figure out how to do. > I'm back online today and can work on this. I did leave some comments in > the form of SetBlockSize statements that are commented out in places were I > had trouble. The problem is fixed.

1 0

Re: [petsc-dev] Multigrid is confusing
by Mark F. Adams 23 May '12

23 May '12

On May 22, 2012, at 11:09 PM, Barry Smith wrote: > > On May 22, 2012, at 9:58 PM, Jed Brown wrote: > >> On Tue, May 22, 2012 at 9:52 PM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: >> How about a default nonlinear Krylov accelerator? fgrmes is not much more expensive than gmres under causal circumstances. >> >> I could handle that. The extra storage concerns me more than time. >> >> I'd eventually like to make PCMG linear by default, especially considering that the convergence with GMRES(1) smoothers tends to be worse in my experience (except when coarse level operators are poorly scaled/unstable) and adds a ton of synchronization points. >> >> Oh, FGMRES by default would necessarily switch from left to right preconditioning. I'm not sure that's a good thing since it requires much tighter tolerances for ill-conditioned systems to reach a given error (especially penalty boundaries). There are advantages both ways and I keep waffling. I expect that changing the default will cause some surprises. > > > So does it make sense to use ILU(0) smoothing with Chebyshev(3) acceleration (with eigenvalues precomputed by GMRES?). I don't like ILU(0) because it does not guarantee to keep the operator positive (eg, elasticity), which is annoying if you have SPD matrix. And if you use Cheby then its obviously bad. Mark > I'm fine with switching now, I just want something reasonably robust so for example ex19 converges well with the defaults. > > > Barry > >

1 0

Re: [petsc-dev] anyone feel like fixing a bug?
by Jed Brown 23 May '12

23 May '12

On Wed, May 23, 2012 at 8:36 AM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > > On May 22, 2012, at 9:32 PM, Jed Brown wrote: > > > Well, here's the problem: > > > > MatSetBlockSize() used to call mat->ops->setblocksize() which the MPI > matrices forwarded down to their seq submatrices. > > > > But then this patch removed that behavior so MatSetBlockSize() is > processed immediately and not forwarded down. > > > > http://petsc.cs.iit.edu/petsc/petsc-dev/rev/80ec81eb9db9 > > > > Since the Seq inner matrices have different PetscLayouts, they don't > inherit the block sizes set on their parent. > > > > Barry, why did you remove this MatOp? How do you want the MPI matrices > to forward block sizes now? > > Blocksize is now suppose to be known "early" now. That is, you cannot > change it after MatSetUp() is called. It is the job of each matrix > implementation to use the already set block size in the appropriate way. > Generally this means before the "Seq" matrices are created or finalized. > One could possible arrange moving the blocksize information down during > the MatSetUp() of that implementation or the MatXXXXSetPreallocation() or > when the inner Mat objects are created. > > If you tell me the EXACT change of events of the "inner" block size not > being set properly I can take a look at it. > This is EXACTLY where the inner blocks now need to be set. ;-) http://petsc.cs.iit.edu/petsc/petsc-dev/rev/0bb703470caa

1 0

Re: [petsc-dev] cheby.c
by Mark F. Adams 23 May '12

23 May '12

Cheby can be tricky, textbook algorithms (eg, Saad's book) have some curve balls. Prometheus has a Cheby implementation that I know is correct. It was worked on by several good people and implemented in several codes and we wrote a paper on it. Mark On May 23, 2012, at 9:21 AM, Barry Smith wrote: > > On May 22, 2012, at 4:17 PM, Hong Zhang wrote: > >> Barry, >> Do you recall the paper or detailed Chebyshev algorithm that you implemented the Chebyshev in petsc? > > No > >> I added 'hybrid' to it, which is not correct. Jie is helping to fix it. Having an original >> algorithm would help. > > The code has also heavily evolved since first being code with bells and whistles being added. > > Perhaps look at the templates book for the definitive algorithm :-) > > > > Barry > >> >> On Tue, May 22, 2012 at 1:21 PM, Jie Chen <jiechen(a)mcs.anl.gov> wrote: >> Then do you know who wrote the code? Hong said it was you. >> >> It was written in the 90s, before version control history. >> Sorry, my mistake. I recall that Barry said that he wrote it when Jed was in dia... ;-( >> Hong >> >> >> Hong is trying to modify the code but there is something wrong. I am helping to see what is the problem. The code itself is fine. >> >> Hard to speculate without more than "something wrong". >> >> >> Jie >> >> >> >> ----- Original Message ----- >> From: "Jed Brown" <jedbrown(a)mcs.anl.gov> >> To: "Jie Chen" <jiechen(a)mcs.anl.gov>, "For users of the development version of PETSc" <petsc-dev(a)mcs.anl.gov> >> Sent: Tuesday, May 22, 2012 1:18:07 PM >> Subject: Re: cheby.c >> >> >> On Tue, May 22, 2012 at 12:37 PM, Jie Chen < jiechen(a)mcs.anl.gov > wrote: >> >> >> Jed, I am reading your implementation of the chebyshev code. The algorithm does not seem to be the same as the version I am familiar with. Did you implement it according to some paper or book? >> >> >> >> I didn't write it. >> >> >> What do you find strange? >> >> > >

1 0

Re: [petsc-dev] anyone feel like fixing a bug?
by Barry Smith 23 May '12

23 May '12

On May 22, 2012, at 9:32 PM, Jed Brown wrote: > Well, here's the problem: > > MatSetBlockSize() used to call mat->ops->setblocksize() which the MPI matrices forwarded down to their seq submatrices. > > But then this patch removed that behavior so MatSetBlockSize() is processed immediately and not forwarded down. > > http://petsc.cs.iit.edu/petsc/petsc-dev/rev/80ec81eb9db9 > > Since the Seq inner matrices have different PetscLayouts, they don't inherit the block sizes set on their parent. > > Barry, why did you remove this MatOp? How do you want the MPI matrices to forward block sizes now? Blocksize is now suppose to be known "early" now. That is, you cannot change it after MatSetUp() is called. It is the job of each matrix implementation to use the already set block size in the appropriate way. Generally this means before the "Seq" matrices are created or finalized. One could possible arrange moving the blocksize information down during the MatSetUp() of that implementation or the MatXXXXSetPreallocation() or when the inner Mat objects are created. If you tell me the EXACT change of events of the "inner" block size not being set properly I can take a look at it. Barry > > On Sat, May 19, 2012 at 10:57 AM, Mark F. Adams <mark.adams(a)columbia.edu> wrote: > I'm going to be tied up for a few days, should be able to get to this Wed., if not sooner, but if anyone wants to take a crack at it, here it is. This is run with 8 procs. > > Mark > > ==51037== Command: ./ex56 -ne 19 -alpha 1.e-3 -ksp_type cg -ksp_monitor -ksp_max_it 15 -pc_type gamg -pc_gamg_type agg -pc_gamg_agg_nsmooths 1 -pc_gamg_verbose 2 -ksp_converged_reason -options_left -theta -30.0 -epsilon .01 -blob_center 0.,0. -pc_gamg_threshold .0 -pc_gamg_coarse_eq_limit 100 -use_coordinates -out_verbose 2 -aggmg_smooths 1 -mg_levels_ksp_type chebyshev -mg_levels_ksp_chebyshev_estimate_eigenvalues 0,0.2,0,1.1 -mg_levels_pc_type pbjacobi -mg_levels_ksp_max_it 1 > > ==51039== Invalid read of size 8 > ==51039== at 0x91E38B: PCApply_PBJacobi_3 (pbjacobi.c:81) > ==51039== by 0xA0FE23: PCApply (precon.c:384) > ==51039== by 0xADA9B2: KSPInitialResidual (itres.c:64) > ==51039== by 0xA7EA39: KSPSolve_GMRES (gmres.c:230) > ==51039== by 0xAB3CA2: KSPSolve (itfunc.c:446) > ==51039== by 0xAA6FF8: KSPSolve_Chebyshev (cheby.c:269) > ==51039== by 0xAB3CA2: KSPSolve (itfunc.c:446) > ==51039== by 0x9ECD2F: PCMGMCycle_Private (mg.c:20) > ==51039== by 0x9F0BEF: PCApply_MG (mg.c:321) > ==51039== by 0xA0FE23: PCApply (precon.c:384) > ==51039== by 0xA43EBF: KSPSolve_CG (cg.c:139) > ==51039== by 0xAB3CA2: KSPSolve (itfunc.c:446) > ==51039== by 0x4805: main (in ./ex56) > ==51039== Address 0x27f4c90 is 0 bytes after a block of size 24,032 alloc'd > ==51039== at 0x1B2AF16: malloc (vg_replace_malloc.c:236) > ==51039== by 0x8A8C7: PetscMallocAlign (mal.c:37) > ==51039== by 0x61E602: MatInvertBlockDiagonal_SeqAIJ (aij.c:2816) > ==51039== by 0x2F0AC3: MatInvertBlockDiagonal (matrix.c:9280) > ==51039== by 0x6D1494: MatInvertBlockDiagonal_MPIAIJ (mpiaij.c:2981) > ==51039== by 0x2F0AC3: MatInvertBlockDiagonal (matrix.c:9280) > ==51039== by 0x920A85: PCSetUp_PBJacobi (pbjacobi.c:213) > ==51039== by 0xA1783B: PCSetUp (precon.c:832) > ==51039== by 0xAB2125: KSPSetUp (itfunc.c:278) > ==51039== by 0x9F536C: PCSetUp_MG (mg.c:678) > ==51039== by 0x98C26C: PCSetUp_GAMG (gamg.c:984) > ==51039== by 0xA1783B: PCSetUp (precon.c:832) > ==51039== by 0xAB2125: KSPSetUp (itfunc.c:278) > ==51039== by 0x46AD: main (in ./ex56) > ==51039== > >

1 0

Re: [petsc-dev] cheby.c
by Barry Smith 23 May '12

23 May '12

On May 22, 2012, at 4:17 PM, Hong Zhang wrote: > Barry, > Do you recall the paper or detailed Chebyshev algorithm that you implemented the Chebyshev in petsc? No > I added 'hybrid' to it, which is not correct. Jie is helping to fix it. Having an original > algorithm would help. The code has also heavily evolved since first being code with bells and whistles being added. Perhaps look at the templates book for the definitive algorithm :-) Barry > > On Tue, May 22, 2012 at 1:21 PM, Jie Chen <jiechen(a)mcs.anl.gov> wrote: > Then do you know who wrote the code? Hong said it was you. > > It was written in the 90s, before version control history. > Sorry, my mistake. I recall that Barry said that he wrote it when Jed was in dia... ;-( > Hong > > > Hong is trying to modify the code but there is something wrong. I am helping to see what is the problem. The code itself is fine. > > Hard to speculate without more than "something wrong". > > > Jie > > > > ----- Original Message ----- > From: "Jed Brown" <jedbrown(a)mcs.anl.gov> > To: "Jie Chen" <jiechen(a)mcs.anl.gov>, "For users of the development version of PETSc" <petsc-dev(a)mcs.anl.gov> > Sent: Tuesday, May 22, 2012 1:18:07 PM > Subject: Re: cheby.c > > > On Tue, May 22, 2012 at 12:37 PM, Jie Chen < jiechen(a)mcs.anl.gov > wrote: > > > Jed, I am reading your implementation of the chebyshev code. The algorithm does not seem to be the same as the version I am familiar with. Did you implement it according to some paper or book? > > > > I didn't write it. > > > What do you find strange? > >

1 0

Re: [petsc-dev] A PETSc question on scicomp I can't completely answer (behavior of PETScSetxxx)
by Matthew Knepley 23 May '12

23 May '12

On Wed, May 23, 2012 at 7:23 AM, Aron Ahmadia <aron(a)ahmadia.net> wrote: > > http://scicomp.stackexchange.com/questions/2303/petscs-xxxsetxxx-methods-ow… > > My understanding is that Set for the most part acts like a > pass-by-reference, and my answer reflects that ( > http://scicomp.stackexchange.com/a/2311/9) > > If there's a better answer than that, Jed or Matt should post it or fix > mine :) > Replied. Matt > A > -- 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