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

Re: [petsc-dev] Fwd: Changes in PETSc
by Matthew Knepley 21 Aug '12

21 Aug '12

On Tue, Aug 21, 2012 at 11:24 AM, Barry Smith <bsmith(a)mcs.anl.gov> wrote: > > On Aug 21, 2012, at 10:53 AM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > > > On Tue, Aug 21, 2012 at 10:39 AM, Matthew Knepley < > knepley(a)ci.uchicago.edu> wrote: > > Wolfgang is not wrong here. I was against the name change :) > > > > I complained about this in May when it was pushed. > > Then who the fuck changed this. It sure as hell wasn't me. I believe Hong. Matt > > Barry > > > > > > > > Matt > > > > ---------- Forwarded message ---------- > > From: Wolfgang Bangerth <bangerth(a)math.tamu.edu> > > Date: Mon, Aug 20, 2012 at 10:30 PM > > Subject: Changes in PETSc > > To: Matthew Knepley <knepley(a)ci.uchicago.edu>, "Toby D. Young" < > tyoung(a)ippt.gov.pl> > > > > > > > > Hi Matt, > > I ran across something today that I believe perfectly illustrates the > kind of thing that drives everyone completely mad about PETSc's habit of > changing everything in every release. I was staring for a long time at a > patch by one of our contributors labeled "fix catastrophic typo" that reads > like this: > > ........................... > > - int ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYSHEV)); > > + // set the type of solver. > > + int ierr; > > + > > +#ifdef DEAL_II_PETSC_VERSION_LT(3,3,0) > > + ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYCHEV)); > > > > They can get rid of the const_cast. > > > > Also, instead of the magic version nonsense, they can write > > > > #ifndef KSPCHEBYSHEV > > # define KSPCHEBYSHEV KSPCHEBYCHEV > > #endif > > > > and then always use the new version in their code. > > > > Also, a package like deal.II would be better off snooping the list of > valid parameters using PetscFListGet() so that they would *automatically* > support *all* KSPs and PCs in PETSc. > > > > +#else > > + ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYSHEV)); > > +#endif > > ........................... > > It took me minutes to find the difference: before 3.3, the name was > spelled chebyChev, after that chebyShev. That's really a completely > unnecessary change -- but hugely disruptive because not compatible, > incredible difficult to see, and it also has no real value at all: whether > you want to spell the name the French or the American way really makes no > difference at all. Worst of all, it would have been absolutely trivial to > make the change backward compatible by simply doing > > #define KSPCHEBYSHEV KSPCHEBYCHEV > > but that didn't happen. > > > > So why do it? Sorry for venting, but it just doesn't make sense to me... > > Best > > W. > > > > -- > > ------------------------------------------------------------------------ > > Wolfgang Bangerth email: bangerth(a)math.tamu.edu > > www: http://www.math.tamu.edu/~bangerth/ > > > > > > > > > > -- > > 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 > > > > -- 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] Fwd: Changes in PETSc
by Barry Smith 21 Aug '12

21 Aug '12

On Aug 21, 2012, at 10:53 AM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > On Tue, Aug 21, 2012 at 10:39 AM, Matthew Knepley <knepley(a)ci.uchicago.edu> wrote: > Wolfgang is not wrong here. I was against the name change :) > > I complained about this in May when it was pushed. Then who the fuck changed this. It sure as hell wasn't me. Barry > > > Matt > > ---------- Forwarded message ---------- > From: Wolfgang Bangerth <bangerth(a)math.tamu.edu> > Date: Mon, Aug 20, 2012 at 10:30 PM > Subject: Changes in PETSc > To: Matthew Knepley <knepley(a)ci.uchicago.edu>, "Toby D. Young" <tyoung(a)ippt.gov.pl> > > > > Hi Matt, > I ran across something today that I believe perfectly illustrates the kind of thing that drives everyone completely mad about PETSc's habit of changing everything in every release. I was staring for a long time at a patch by one of our contributors labeled "fix catastrophic typo" that reads like this: > ........................... > - int ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYSHEV)); > + // set the type of solver. > + int ierr; > + > +#ifdef DEAL_II_PETSC_VERSION_LT(3,3,0) > + ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYCHEV)); > > They can get rid of the const_cast. > > Also, instead of the magic version nonsense, they can write > > #ifndef KSPCHEBYSHEV > # define KSPCHEBYSHEV KSPCHEBYCHEV > #endif > > and then always use the new version in their code. > > Also, a package like deal.II would be better off snooping the list of valid parameters using PetscFListGet() so that they would *automatically* support *all* KSPs and PCs in PETSc. > > +#else > + ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYSHEV)); > +#endif > ........................... > It took me minutes to find the difference: before 3.3, the name was spelled chebyChev, after that chebyShev. That's really a completely unnecessary change -- but hugely disruptive because not compatible, incredible difficult to see, and it also has no real value at all: whether you want to spell the name the French or the American way really makes no difference at all. Worst of all, it would have been absolutely trivial to make the change backward compatible by simply doing > #define KSPCHEBYSHEV KSPCHEBYCHEV > but that didn't happen. > > So why do it? Sorry for venting, but it just doesn't make sense to me... > Best > W. > > -- > ------------------------------------------------------------------------ > Wolfgang Bangerth email: bangerth(a)math.tamu.edu > www: http://www.math.tamu.edu/~bangerth/ > > > > > -- > 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] FieldSplit fixes in 3.3
by Jed Brown 21 Aug '12

21 Aug '12

Hong, as far as I can tell, none of these changes should be merged to petsc-dev right now. If you go through with the merge, use hg diff -r SHA1_OF_UPSTREAM_PETSC_DEV to make sure that none of this 3.3 stuff gets merged. Otherwise, one of us can do the merge. On Tue, Aug 21, 2012 at 11:08 AM, Hong Zhang <hzhang(a)mcs.anl.gov> wrote: > Dmitry : > I pushed a bugfix of mumps interface to 3.3. When merging it to petsc-dev, > I get > > merging src/ksp/pc/impls/fieldsplit/fieldsplit.c failed! > > Hong > > >> >> On Mon, Aug 20, 2012 at 8:36 AM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: >> >>> On Mon, Aug 20, 2012 at 4:02 AM, Dmitry Karpeev <karpeev(a)mcs.anl.gov>wrote: >>> >>>> >>>> >>>> On Sun, Aug 19, 2012 at 1:17 PM, Jed Brown <jedbrown(a)mcs.anl.gov>wrote: >>>> >>>>> On Sun, Aug 19, 2012 at 12:48 PM, Dmitry Karpeev <karpeev(a)mcs.anl.gov>wrote: >>>>> >>>>>> The main reason for this being a patch in 3.3 is that recursive >>>>>> FieldSplit is broken there (which is why I couldn't enable it in >>>>>> libMesh/Moose correctly). >>>>> >>>>> >>>>> The current interface doesn't explicitly support all the ways to get >>>>> information into the split, but several people have worked around the >>>>> limitations to get the necessary information in there. Having you and Matt >>>>> changing things to rely on mutually incompatible side-effects does not >>>>> help. >>>>> >>>> I don't think the present fix is about the API or working through side >>>> effects (nullspace stuff excluded and partially backed out here: >>>> http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/1723d4624521) >>>> >>> >>> What do you call this external DM futzing? >>> >> It's no different than holding onto the ISs set in >> PCFieldSplitSetDefaults() or directly via PCFieldSplitSetIS() until they >> can be used in PCSetUp_FieldSplit(). I don't see how this can be avoided >> here or in the future, unless we do something to fundamentally change the >> set up process. >> >>> >>> >>>> The current interface assumes that the A00 solver and the Schur inner >>>> solver are to be set up identically >>>> >>> >>> If they are being set up identically, they should literally be the same >>> object. This is often the biggest setup cost in the whole problem (e.g. >>> AMG). >>> >> I agree, but that's available only in petsc-dev, not in petsc-3.3. The >> duplicate setup has been occurring all along, now it is actually consistent >> with the DM being forwarded to the inner solver. >> >>> >>> >>>> (something we relaxed in petsc-dev). In order to use fieldsplit on >>>> these solvers the DM for the corresponding split must be forwarded to both >>>> of these A00 solvers. It was being set only on the "outer" A00 solver. >>>> Partly this is an artifact of our "split" setup process: some of it occurs >>>> in PCFieldSplitSetDefaults(), some later in PCSetUp_FieldSplit(). Some >>>> objects obtained in PCFieldSplitSetDefaults() need to be cached until >>>> PCSetUp_FieldSplit(), and the splits' DMs are among those. I now put in >>>> proper reference counting for them here: >>>> http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/188af9799779 >>>> >>>> >>>>> Maybe we can reach some agreement on the proper way to do things. >>>>> >>>> That's fine with me. The reason I wanted to fix this particular >>>> problem now is that some dependent packages (e.g., Moose) only rely on >>>> release versions of petsc and would not be able to use this functionality >>>> properly for quite some time. >>>> >>> >>> Why do packages do this? Because petsc-dev is sometimes unstable. >>> >>> Now what happens when code that is pushed to the release changes >>> behavior and uses uninitialized memory? The release becomes unstable. >>> >>> Even worse, our nightly testing is targeted at petsc-dev because the >>> expectation is that only trivial and well-tested bug-fix patches are pushed >>> to release. If we are changing that model by introducing significant >>> changes, we MUST have nightly tests for the release, we have to actually >>> look at the results, and I would be strongly in favor of switching from >>> 3.3p3 numbering to 3.3.1. >>> >>> In all cases, these subminor or patch releases MUST be binary and >>> source-level backward compatible. >>> >>> >>>> Let's fix that "down the line". Setting a DM should never force it to >>>>> be used or cause an error due to unsupported operation in a case where not >>>>> having a DM is also acceptable. >>>>> >>>> Currently the inner and outer A00 solvers are essentially identified >>>> (duplicated), so they should be set up identically, including the DM. >>>> Especially when that DM defines a recursive split. This is in part >>>> because there is no way to configure those two differently programmatically >>>> or from the command line. So if the outer A00 is using >>>> -fieldsplit_0_pc_type_fieldsplit, the inner A00 will as well. >>>> However, the DM on the inner A00 will be rather different (or absent) >>>> and produce splits incompatible with the rest of the options (e.g., a >>>> single default split, while -fieldsplit_0_pc_fieldsplit_type schur). >>>> >>>> Incidentally, this raises this question down the road: in petsc-dev the >>>> inner and outer A00 can be configure separately using different prefixes, >>>> but how should the inner A00 DM be configured? >>>> >>> >>> The KSPs should only be different if the solver is different. Having >>> them separate, but using the same DM is one reason I was not wild about >>> having physics in the DM, but I don't see a good way to plumb in the extra >>> information. In any case, we need to come up with a solution before pushing >>> code. >>> >>> >>>> Okay, I'll back out that part of the second patch, since A11's >>>>>> nullspace is assumed to be meant for S. I'm not sure what's incorrect >>>>>> about the comment, though: if a vector is in the A11's kernel, but not in >>>>>> S's, I call it a "false positive". This terminology may be wanting, but >>>>>> what's incorrect about it? >>>>>> >>>>> >>>>> "force an inconsistent rhs" >>>>> >>>> The KSP solve will fail only if the search space X is contracted to the >>>> point where the rhs b is not in the range AX, hence, >>>> inconsistent. >>>> >>> >>> Okay, I just would have said that it changes the operator. Also, the >>> Krylov methods are not necessarily correct if the matrix does not have the >>> stated null space (since it's expected that the matrix has the null space >>> and only the preconditioner needs to be filtered). >>> >>> >>>> I backed out this hunk, however. >>>> >>>>> >>>>> S may well be nonsingular. >>>>> >>>>> >>>>>> >>>>>> >>>>>>> I think that's confusing and we need a specific API for it, but if >>>>>>> people are using it that way, we shouldn't change it in 3.3. >>>>>>> >>>>>>> I'd prefer to backout the first patch, split it into separate >>>>>>> pieces, and review each before pushing. I think it's a potentially big >>>>>>> behavior change for 3.3. >>>>>>> >>>>>> Which pieces of the first patch do you think are too big for 3.3? >>>>>> The splits' and Schur KSPs have to be set up (maybe not all in every >>>>>> situation) in order for recursive splitting to work. Should we declare >>>>>> that capability unavailable for 3.3? I can eliminate the unneeded >>>>>> KSPSetUp() calls. >>>>>> >>>>> >>>>> 1. I hate those hacky dm reference-non-references. >>>>> >>>>> 2. I don't want unused KSPs to be set up. (We should eventually fix >>>>> the model so they don't exist...) >>>>> >>>>> 3. Two of the three IncrementTabLevels that you introduced are not >>>>> necessary because the KSP's tab level was incremented before getting the >>>>> PC. >>>>> >>>> Where? >>>> >>> >>> ierr = >>> KSPCreate(((PetscObject)pc)->comm,&jac->kspschur);CHKERRQ(ierr); >>> ierr = >>> PetscLogObjectParent((PetscObject)pc,(PetscObject)jac->kspschur);CHKERRQ(ierr); >>> ierr = >>> PetscObjectIncrementTabLevel((PetscObject)jac->kspschur,(PetscObject)pc,1);CHKERRQ(ierr); >>> + { >>> + PC pcschur; >>> + ierr = KSPGetPC(jac->kspschur, &pcschur); >>> CHKERRQ(ierr); >>> + ierr = >>> PetscObjectIncrementTabLevel((PetscObject)pcschur,(PetscObject)pc,1);CHKERRQ(ierr); >>> + } >>> >>> >>> ierr = >>> KSPCreate(((PetscObject)pc)->comm,&ilink->ksp);CHKERRQ(ierr); >>> ierr = >>> PetscObjectIncrementTabLevel((PetscObject)ilink->ksp,(PetscObject)pc,1);CHKERRQ(ierr); >>> + { >>> + PC ilinkpc; >>> + ierr = KSPGetPC(ilink->ksp, &ilinkpc); CHKERRQ(ierr); >>> + ierr = >>> PetscObjectIncrementTabLevel((PetscObject)ilinkpc,(PetscObject)pc,1);CHKERRQ(ierr); >>> + } >>> >>> >>> This is the pattern used everywhere else in PETSc so that sub-pcs always >>> have the correct tab level. >>> >> >> >

1 0

Re: [petsc-dev] FieldSplit fixes in 3.3
by Hong Zhang 21 Aug '12

21 Aug '12

Dmitry : I pushed a bugfix of mumps interface to 3.3. When merging it to petsc-dev, I get merging src/ksp/pc/impls/fieldsplit/fieldsplit.c failed! Hong > > On Mon, Aug 20, 2012 at 8:36 AM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > >> On Mon, Aug 20, 2012 at 4:02 AM, Dmitry Karpeev <karpeev(a)mcs.anl.gov>wrote: >> >>> >>> >>> On Sun, Aug 19, 2012 at 1:17 PM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: >>> >>>> On Sun, Aug 19, 2012 at 12:48 PM, Dmitry Karpeev <karpeev(a)mcs.anl.gov>wrote: >>>> >>>>> The main reason for this being a patch in 3.3 is that recursive >>>>> FieldSplit is broken there (which is why I couldn't enable it in >>>>> libMesh/Moose correctly). >>>> >>>> >>>> The current interface doesn't explicitly support all the ways to get >>>> information into the split, but several people have worked around the >>>> limitations to get the necessary information in there. Having you and Matt >>>> changing things to rely on mutually incompatible side-effects does not >>>> help. >>>> >>> I don't think the present fix is about the API or working through side >>> effects (nullspace stuff excluded and partially backed out here: >>> http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/1723d4624521) >>> >> >> What do you call this external DM futzing? >> > It's no different than holding onto the ISs set in > PCFieldSplitSetDefaults() or directly via PCFieldSplitSetIS() until they > can be used in PCSetUp_FieldSplit(). I don't see how this can be avoided > here or in the future, unless we do something to fundamentally change the > set up process. > >> >> >>> The current interface assumes that the A00 solver and the Schur inner >>> solver are to be set up identically >>> >> >> If they are being set up identically, they should literally be the same >> object. This is often the biggest setup cost in the whole problem (e.g. >> AMG). >> > I agree, but that's available only in petsc-dev, not in petsc-3.3. The > duplicate setup has been occurring all along, now it is actually consistent > with the DM being forwarded to the inner solver. > >> >> >>> (something we relaxed in petsc-dev). In order to use fieldsplit on >>> these solvers the DM for the corresponding split must be forwarded to both >>> of these A00 solvers. It was being set only on the "outer" A00 solver. >>> Partly this is an artifact of our "split" setup process: some of it occurs >>> in PCFieldSplitSetDefaults(), some later in PCSetUp_FieldSplit(). Some >>> objects obtained in PCFieldSplitSetDefaults() need to be cached until >>> PCSetUp_FieldSplit(), and the splits' DMs are among those. I now put in >>> proper reference counting for them here: >>> http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/188af9799779 >>> >>> >>>> Maybe we can reach some agreement on the proper way to do things. >>>> >>> That's fine with me. The reason I wanted to fix this particular problem >>> now is that some dependent packages (e.g., Moose) only rely on release >>> versions of petsc and would not be able to use this functionality properly >>> for quite some time. >>> >> >> Why do packages do this? Because petsc-dev is sometimes unstable. >> >> Now what happens when code that is pushed to the release changes behavior >> and uses uninitialized memory? The release becomes unstable. >> >> Even worse, our nightly testing is targeted at petsc-dev because the >> expectation is that only trivial and well-tested bug-fix patches are pushed >> to release. If we are changing that model by introducing significant >> changes, we MUST have nightly tests for the release, we have to actually >> look at the results, and I would be strongly in favor of switching from >> 3.3p3 numbering to 3.3.1. >> >> In all cases, these subminor or patch releases MUST be binary and >> source-level backward compatible. >> >> >>> Let's fix that "down the line". Setting a DM should never force it to be >>>> used or cause an error due to unsupported operation in a case where not >>>> having a DM is also acceptable. >>>> >>> Currently the inner and outer A00 solvers are essentially identified >>> (duplicated), so they should be set up identically, including the DM. >>> Especially when that DM defines a recursive split. This is in part >>> because there is no way to configure those two differently programmatically >>> or from the command line. So if the outer A00 is using >>> -fieldsplit_0_pc_type_fieldsplit, the inner A00 will as well. >>> However, the DM on the inner A00 will be rather different (or absent) >>> and produce splits incompatible with the rest of the options (e.g., a >>> single default split, while -fieldsplit_0_pc_fieldsplit_type schur). >>> >>> Incidentally, this raises this question down the road: in petsc-dev the >>> inner and outer A00 can be configure separately using different prefixes, >>> but how should the inner A00 DM be configured? >>> >> >> The KSPs should only be different if the solver is different. Having them >> separate, but using the same DM is one reason I was not wild about having >> physics in the DM, but I don't see a good way to plumb in the extra >> information. In any case, we need to come up with a solution before pushing >> code. >> >> >>> Okay, I'll back out that part of the second patch, since A11's >>>>> nullspace is assumed to be meant for S. I'm not sure what's incorrect >>>>> about the comment, though: if a vector is in the A11's kernel, but not in >>>>> S's, I call it a "false positive". This terminology may be wanting, but >>>>> what's incorrect about it? >>>>> >>>> >>>> "force an inconsistent rhs" >>>> >>> The KSP solve will fail only if the search space X is contracted to the >>> point where the rhs b is not in the range AX, hence, >>> inconsistent. >>> >> >> Okay, I just would have said that it changes the operator. Also, the >> Krylov methods are not necessarily correct if the matrix does not have the >> stated null space (since it's expected that the matrix has the null space >> and only the preconditioner needs to be filtered). >> >> >>> I backed out this hunk, however. >>> >>>> >>>> S may well be nonsingular. >>>> >>>> >>>>> >>>>> >>>>>> I think that's confusing and we need a specific API for it, but if >>>>>> people are using it that way, we shouldn't change it in 3.3. >>>>>> >>>>>> I'd prefer to backout the first patch, split it into separate >>>>>> pieces, and review each before pushing. I think it's a potentially big >>>>>> behavior change for 3.3. >>>>>> >>>>> Which pieces of the first patch do you think are too big for 3.3? The >>>>> splits' and Schur KSPs have to be set up (maybe not all in every situation) >>>>> in order for recursive splitting to work. Should we declare that >>>>> capability unavailable for 3.3? I can eliminate the unneeded KSPSetUp() >>>>> calls. >>>>> >>>> >>>> 1. I hate those hacky dm reference-non-references. >>>> >>>> 2. I don't want unused KSPs to be set up. (We should eventually fix the >>>> model so they don't exist...) >>>> >>>> 3. Two of the three IncrementTabLevels that you introduced are not >>>> necessary because the KSP's tab level was incremented before getting the >>>> PC. >>>> >>> Where? >>> >> >> ierr = >> KSPCreate(((PetscObject)pc)->comm,&jac->kspschur);CHKERRQ(ierr); >> ierr = >> PetscLogObjectParent((PetscObject)pc,(PetscObject)jac->kspschur);CHKERRQ(ierr); >> ierr = >> PetscObjectIncrementTabLevel((PetscObject)jac->kspschur,(PetscObject)pc,1);CHKERRQ(ierr); >> + { >> + PC pcschur; >> + ierr = KSPGetPC(jac->kspschur, &pcschur); >> CHKERRQ(ierr); >> + ierr = >> PetscObjectIncrementTabLevel((PetscObject)pcschur,(PetscObject)pc,1);CHKERRQ(ierr); >> + } >> >> >> ierr = >> KSPCreate(((PetscObject)pc)->comm,&ilink->ksp);CHKERRQ(ierr); >> ierr = >> PetscObjectIncrementTabLevel((PetscObject)ilink->ksp,(PetscObject)pc,1);CHKERRQ(ierr); >> + { >> + PC ilinkpc; >> + ierr = KSPGetPC(ilink->ksp, &ilinkpc); CHKERRQ(ierr); >> + ierr = >> PetscObjectIncrementTabLevel((PetscObject)ilinkpc,(PetscObject)pc,1);CHKERRQ(ierr); >> + } >> >> >> This is the pattern used everywhere else in PETSc so that sub-pcs always >> have the correct tab level. >> > >

1 0

Re: [petsc-dev] Fwd: Changes in PETSc
by Jed Brown 21 Aug '12

21 Aug '12

On Tue, Aug 21, 2012 at 10:39 AM, Matthew Knepley <knepley(a)ci.uchicago.edu>wrote: > Wolfgang is not wrong here. I was against the name change :) I complained about this in May when it was pushed. > > Matt > > ---------- Forwarded message ---------- > From: Wolfgang Bangerth <bangerth(a)math.tamu.edu> > Date: Mon, Aug 20, 2012 at 10:30 PM > Subject: Changes in PETSc > To: Matthew Knepley <knepley(a)ci.uchicago.edu>, "Toby D. Young" < > tyoung(a)ippt.gov.pl> > > > > Hi Matt, > I ran across something today that I believe perfectly illustrates the kind > of thing that drives everyone completely mad about PETSc's habit of > changing everything in every release. I was staring for a long time at a > patch by one of our contributors labeled "fix catastrophic typo" that reads > like this: > ........................... > - int ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYSHEV)); > + // set the type of solver. > + int ierr; > + > +#ifdef DEAL_II_PETSC_VERSION_LT(3,3,**0) > + ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYCHEV)); > They can get rid of the const_cast. Also, instead of the magic version nonsense, they can write #ifndef KSPCHEBYSHEV # define KSPCHEBYSHEV KSPCHEBYCHEV #endif and then always use the new version in their code. Also, a package like deal.II would be better off snooping the list of valid parameters using PetscFListGet() so that they would *automatically* support *all* KSPs and PCs in PETSc. > +#else > + ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYSHEV)); > +#endif > ........................... > It took me minutes to find the difference: before 3.3, the name was > spelled chebyChev, after that chebyShev. That's really a completely > unnecessary change -- but hugely disruptive because not compatible, > incredible difficult to see, and it also has no real value at all: whether > you want to spell the name the French or the American way really makes no > difference at all. Worst of all, it would have been absolutely trivial to > make the change backward compatible by simply doing > #define KSPCHEBYSHEV KSPCHEBYCHEV > but that didn't happen. > > So why do it? Sorry for venting, but it just doesn't make sense to me... > Best > W. > > -- > ------------------------------**------------------------------** > ------------ > Wolfgang Bangerth email: bangerth(a)math.tamu.edu > www: http://www.math.tamu.edu/~**bangerth/<http://www.math.tamu.edu/~bangerth/> > > > > > -- > 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

Fwd: Changes in PETSc
by Matthew Knepley 21 Aug '12

21 Aug '12

Wolfgang is not wrong here. I was against the name change :) Matt ---------- Forwarded message ---------- From: Wolfgang Bangerth <bangerth(a)math.tamu.edu> Date: Mon, Aug 20, 2012 at 10:30 PM Subject: Changes in PETSc To: Matthew Knepley <knepley(a)ci.uchicago.edu>, "Toby D. Young" < tyoung(a)ippt.gov.pl> Hi Matt, I ran across something today that I believe perfectly illustrates the kind of thing that drives everyone completely mad about PETSc's habit of changing everything in every release. I was staring for a long time at a patch by one of our contributors labeled "fix catastrophic typo" that reads like this: ........................... - int ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYSHEV)); + // set the type of solver. + int ierr; + +#ifdef DEAL_II_PETSC_VERSION_LT(3,3,**0) + ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYCHEV)); +#else + ierr = KSPSetType (ksp, const_cast<char *>(KSPCHEBYSHEV)); +#endif ........................... It took me minutes to find the difference: before 3.3, the name was spelled chebyChev, after that chebyShev. That's really a completely unnecessary change -- but hugely disruptive because not compatible, incredible difficult to see, and it also has no real value at all: whether you want to spell the name the French or the American way really makes no difference at all. Worst of all, it would have been absolutely trivial to make the change backward compatible by simply doing #define KSPCHEBYSHEV KSPCHEBYCHEV but that didn't happen. So why do it? Sorry for venting, but it just doesn't make sense to me... Best W. -- ------------------------------**------------------------------**------------ Wolfgang Bangerth email: bangerth(a)math.tamu.edu www: http://www.math.tamu.edu/~**bangerth/<http://www.math.tamu.edu/~bangerth/> -- 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] FieldSplit fixes in 3.3
by Dmitry Karpeev 21 Aug '12

21 Aug '12

On Mon, Aug 20, 2012 at 8:36 AM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > On Mon, Aug 20, 2012 at 4:02 AM, Dmitry Karpeev <karpeev(a)mcs.anl.gov>wrote: > >> >> >> On Sun, Aug 19, 2012 at 1:17 PM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: >> >>> On Sun, Aug 19, 2012 at 12:48 PM, Dmitry Karpeev <karpeev(a)mcs.anl.gov>wrote: >>> >>>> The main reason for this being a patch in 3.3 is that recursive >>>> FieldSplit is broken there (which is why I couldn't enable it in >>>> libMesh/Moose correctly). >>> >>> >>> The current interface doesn't explicitly support all the ways to get >>> information into the split, but several people have worked around the >>> limitations to get the necessary information in there. Having you and Matt >>> changing things to rely on mutually incompatible side-effects does not >>> help. >>> >> I don't think the present fix is about the API or working through side >> effects (nullspace stuff excluded and partially backed out here: >> http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/1723d4624521) >> > > What do you call this external DM futzing? > It's no different than holding onto the ISs set in PCFieldSplitSetDefaults() or directly via PCFieldSplitSetIS() until they can be used in PCSetUp_FieldSplit(). I don't see how this can be avoided here or in the future, unless we do something to fundamentally change the set up process. > > >> The current interface assumes that the A00 solver and the Schur inner >> solver are to be set up identically >> > > If they are being set up identically, they should literally be the same > object. This is often the biggest setup cost in the whole problem (e.g. > AMG). > I agree, but that's available only in petsc-dev, not in petsc-3.3. The duplicate setup has been occurring all along, now it is actually consistent with the DM being forwarded to the inner solver. > > >> (something we relaxed in petsc-dev). In order to use fieldsplit on these >> solvers the DM for the corresponding split must be forwarded to both of >> these A00 solvers. It was being set only on the "outer" A00 solver. >> Partly this is an artifact of our "split" setup process: some of it occurs >> in PCFieldSplitSetDefaults(), some later in PCSetUp_FieldSplit(). Some >> objects obtained in PCFieldSplitSetDefaults() need to be cached until >> PCSetUp_FieldSplit(), and the splits' DMs are among those. I now put in >> proper reference counting for them here: >> http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/188af9799779 >> >> >>> Maybe we can reach some agreement on the proper way to do things. >>> >> That's fine with me. The reason I wanted to fix this particular problem >> now is that some dependent packages (e.g., Moose) only rely on release >> versions of petsc and would not be able to use this functionality properly >> for quite some time. >> > > Why do packages do this? Because petsc-dev is sometimes unstable. > > Now what happens when code that is pushed to the release changes behavior > and uses uninitialized memory? The release becomes unstable. > > Even worse, our nightly testing is targeted at petsc-dev because the > expectation is that only trivial and well-tested bug-fix patches are pushed > to release. If we are changing that model by introducing significant > changes, we MUST have nightly tests for the release, we have to actually > look at the results, and I would be strongly in favor of switching from > 3.3p3 numbering to 3.3.1. > > In all cases, these subminor or patch releases MUST be binary and > source-level backward compatible. > > >> Let's fix that "down the line". Setting a DM should never force it to be >>> used or cause an error due to unsupported operation in a case where not >>> having a DM is also acceptable. >>> >> Currently the inner and outer A00 solvers are essentially identified >> (duplicated), so they should be set up identically, including the DM. >> Especially when that DM defines a recursive split. This is in part >> because there is no way to configure those two differently programmatically >> or from the command line. So if the outer A00 is using >> -fieldsplit_0_pc_type_fieldsplit, the inner A00 will as well. >> However, the DM on the inner A00 will be rather different (or absent) and >> produce splits incompatible with the rest of the options (e.g., a single >> default split, while -fieldsplit_0_pc_fieldsplit_type schur). >> >> Incidentally, this raises this question down the road: in petsc-dev the >> inner and outer A00 can be configure separately using different prefixes, >> but how should the inner A00 DM be configured? >> > > The KSPs should only be different if the solver is different. Having them > separate, but using the same DM is one reason I was not wild about having > physics in the DM, but I don't see a good way to plumb in the extra > information. In any case, we need to come up with a solution before pushing > code. > > >> Okay, I'll back out that part of the second patch, since A11's >>>> nullspace is assumed to be meant for S. I'm not sure what's incorrect >>>> about the comment, though: if a vector is in the A11's kernel, but not in >>>> S's, I call it a "false positive". This terminology may be wanting, but >>>> what's incorrect about it? >>>> >>> >>> "force an inconsistent rhs" >>> >> The KSP solve will fail only if the search space X is contracted to the >> point where the rhs b is not in the range AX, hence, >> inconsistent. >> > > Okay, I just would have said that it changes the operator. Also, the > Krylov methods are not necessarily correct if the matrix does not have the > stated null space (since it's expected that the matrix has the null space > and only the preconditioner needs to be filtered). > > >> I backed out this hunk, however. >> >>> >>> S may well be nonsingular. >>> >>> >>>> >>>> >>>>> I think that's confusing and we need a specific API for it, but if >>>>> people are using it that way, we shouldn't change it in 3.3. >>>>> >>>>> I'd prefer to backout the first patch, split it into separate pieces, >>>>> and review each before pushing. I think it's a potentially big behavior >>>>> change for 3.3. >>>>> >>>> Which pieces of the first patch do you think are too big for 3.3? The >>>> splits' and Schur KSPs have to be set up (maybe not all in every situation) >>>> in order for recursive splitting to work. Should we declare that >>>> capability unavailable for 3.3? I can eliminate the unneeded KSPSetUp() >>>> calls. >>>> >>> >>> 1. I hate those hacky dm reference-non-references. >>> >>> 2. I don't want unused KSPs to be set up. (We should eventually fix the >>> model so they don't exist...) >>> >>> 3. Two of the three IncrementTabLevels that you introduced are not >>> necessary because the KSP's tab level was incremented before getting the >>> PC. >>> >> Where? >> > > ierr = > KSPCreate(((PetscObject)pc)->comm,&jac->kspschur);CHKERRQ(ierr); > ierr = > PetscLogObjectParent((PetscObject)pc,(PetscObject)jac->kspschur);CHKERRQ(ierr); > ierr = > PetscObjectIncrementTabLevel((PetscObject)jac->kspschur,(PetscObject)pc,1);CHKERRQ(ierr); > + { > + PC pcschur; > + ierr = KSPGetPC(jac->kspschur, &pcschur); CHKERRQ(ierr); > + ierr = > PetscObjectIncrementTabLevel((PetscObject)pcschur,(PetscObject)pc,1);CHKERRQ(ierr); > + } > > > ierr = > KSPCreate(((PetscObject)pc)->comm,&ilink->ksp);CHKERRQ(ierr); > ierr = > PetscObjectIncrementTabLevel((PetscObject)ilink->ksp,(PetscObject)pc,1);CHKERRQ(ierr); > + { > + PC ilinkpc; > + ierr = KSPGetPC(ilink->ksp, &ilinkpc); CHKERRQ(ierr); > + ierr = > PetscObjectIncrementTabLevel((PetscObject)ilinkpc,(PetscObject)pc,1);CHKERRQ(ierr); > + } > > > This is the pattern used everywhere else in PETSc so that sub-pcs always > have the correct tab level. >

1 0

Re: [petsc-dev] [petsc-users] GAMG setup scalability
by Hong Zhang 20 Aug '12

20 Aug '12

John: Thanks for reporting your experience with petsc mat-mat operations, which are consistent with current status of these operations. Petsc mat-mat operations were written many years ago for users who wanted these operations for convenience on small size problems, not intended for huge size and large number of cores. As Mark started using these operations on 10k+ cores for gamg, I spent several weeks reworking on some of these operations. See attached slides on current status of these operations. In summary: 1. C=A*B were rewritten with two versions: a) fast but not memory scalable: default b) slow but scalable: use option '-matmatmult_scalable' I have a rough idea to merge these two version to a single one that is fast and scalable, but need time to work on it :-( 2. SPARSE inner-product is more expensive than sparse outer-product. As you observed, C=A*B^t is worse than C=A*B. Reading your data, I'm surprised that your C=A*B^t is so much slower than PtAP. How large is your matrices? Can you run your np=1 case with option '-matmattrnsmult_color' and send us the output of -log_summary? 3. Current PtAP is a simple cleanup of old version, thus it needs update - has been on my ToDo for almost a year :-( According to the analysis (see slides), PtAP should only take 2x of time(A*B), but it is taking 8x A*B on 14k cores (see slide 8). Improving this operations is on the top of our ToDo list. Hong On Fri, Aug 17, 2012 at 11:27 AM, Mark F. Adams <mark.adams(a)columbia.edu> > wrote: > > John: could you perhaps test MatMatMult and MatMatTransposeMult so we > could > > get some concrete numbers? > > For MatMatTransposeMult I get this error: > > [0]PETSC ERROR: --------------------- Error Message > ------------------------------------ > [0]PETSC ERROR: No support for this operation for this object type! > [0]PETSC ERROR: MatMatTransposeMult not supported for A of type mpiaij! > [0]PETSC ERROR: > ------------------------------------------------------------------------ > [0]PETSC ERROR: Petsc Development HG revision: > f9c6cac2d69c724a2258d4e0dc2f51b0825aa875 HG Date: Thu Aug 16 08:37:18 > 2012 -0700 > > For MatMatMult, my problem is not symmetric. It doesn't look like > MatTranspose() followed by MatMatMult is a viable option. > > John >

1 0

Re: [petsc-dev] FieldSplit fixes in 3.3
by Jed Brown 20 Aug '12

20 Aug '12

On Mon, Aug 20, 2012 at 4:02 AM, Dmitry Karpeev <karpeev(a)mcs.anl.gov> wrote: > > > On Sun, Aug 19, 2012 at 1:17 PM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > >> On Sun, Aug 19, 2012 at 12:48 PM, Dmitry Karpeev <karpeev(a)mcs.anl.gov>wrote: >> >>> The main reason for this being a patch in 3.3 is that recursive >>> FieldSplit is broken there (which is why I couldn't enable it in >>> libMesh/Moose correctly). >> >> >> The current interface doesn't explicitly support all the ways to get >> information into the split, but several people have worked around the >> limitations to get the necessary information in there. Having you and Matt >> changing things to rely on mutually incompatible side-effects does not >> help. >> > I don't think the present fix is about the API or working through side > effects (nullspace stuff excluded and partially backed out here: > http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/1723d4624521) > What do you call this external DM futzing? > The current interface assumes that the A00 solver and the Schur inner > solver are to be set up identically > If they are being set up identically, they should literally be the same object. This is often the biggest setup cost in the whole problem (e.g. AMG). > (something we relaxed in petsc-dev). In order to use fieldsplit on these > solvers the DM for the corresponding split must be forwarded to both of > these A00 solvers. It was being set only on the "outer" A00 solver. > Partly this is an artifact of our "split" setup process: some of it occurs > in PCFieldSplitSetDefaults(), some later in PCSetUp_FieldSplit(). Some > objects obtained in PCFieldSplitSetDefaults() need to be cached until > PCSetUp_FieldSplit(), and the splits' DMs are among those. I now put in > proper reference counting for them here: > http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/188af9799779 > > >> Maybe we can reach some agreement on the proper way to do things. >> > That's fine with me. The reason I wanted to fix this particular problem > now is that some dependent packages (e.g., Moose) only rely on release > versions of petsc and would not be able to use this functionality properly > for quite some time. > Why do packages do this? Because petsc-dev is sometimes unstable. Now what happens when code that is pushed to the release changes behavior and uses uninitialized memory? The release becomes unstable. Even worse, our nightly testing is targeted at petsc-dev because the expectation is that only trivial and well-tested bug-fix patches are pushed to release. If we are changing that model by introducing significant changes, we MUST have nightly tests for the release, we have to actually look at the results, and I would be strongly in favor of switching from 3.3p3 numbering to 3.3.1. In all cases, these subminor or patch releases MUST be binary and source-level backward compatible. > Let's fix that "down the line". Setting a DM should never force it to be >> used or cause an error due to unsupported operation in a case where not >> having a DM is also acceptable. >> > Currently the inner and outer A00 solvers are essentially identified > (duplicated), so they should be set up identically, including the DM. > Especially when that DM defines a recursive split. This is in part > because there is no way to configure those two differently programmatically > or from the command line. So if the outer A00 is using > -fieldsplit_0_pc_type_fieldsplit, the inner A00 will as well. > However, the DM on the inner A00 will be rather different (or absent) and > produce splits incompatible with the rest of the options (e.g., a single > default split, while -fieldsplit_0_pc_fieldsplit_type schur). > > Incidentally, this raises this question down the road: in petsc-dev the > inner and outer A00 can be configure separately using different prefixes, > but how should the inner A00 DM be configured? > The KSPs should only be different if the solver is different. Having them separate, but using the same DM is one reason I was not wild about having physics in the DM, but I don't see a good way to plumb in the extra information. In any case, we need to come up with a solution before pushing code. > Okay, I'll back out that part of the second patch, since A11's nullspace >>> is assumed to be meant for S. I'm not sure what's incorrect about the >>> comment, though: if a vector is in the A11's kernel, but not in S's, I call >>> it a "false positive". This terminology may be wanting, but what's >>> incorrect about it? >>> >> >> "force an inconsistent rhs" >> > The KSP solve will fail only if the search space X is contracted to the > point where the rhs b is not in the range AX, hence, > inconsistent. > Okay, I just would have said that it changes the operator. Also, the Krylov methods are not necessarily correct if the matrix does not have the stated null space (since it's expected that the matrix has the null space and only the preconditioner needs to be filtered). > I backed out this hunk, however. > >> >> S may well be nonsingular. >> >> >>> >>> >>>> I think that's confusing and we need a specific API for it, but if >>>> people are using it that way, we shouldn't change it in 3.3. >>>> >>>> I'd prefer to backout the first patch, split it into separate pieces, >>>> and review each before pushing. I think it's a potentially big behavior >>>> change for 3.3. >>>> >>> Which pieces of the first patch do you think are too big for 3.3? The >>> splits' and Schur KSPs have to be set up (maybe not all in every situation) >>> in order for recursive splitting to work. Should we declare that >>> capability unavailable for 3.3? I can eliminate the unneeded KSPSetUp() >>> calls. >>> >> >> 1. I hate those hacky dm reference-non-references. >> >> 2. I don't want unused KSPs to be set up. (We should eventually fix the >> model so they don't exist...) >> >> 3. Two of the three IncrementTabLevels that you introduced are not >> necessary because the KSP's tab level was incremented before getting the >> PC. >> > Where? > ierr = KSPCreate(((PetscObject)pc)->comm,&jac->kspschur);CHKERRQ(ierr); ierr = PetscLogObjectParent((PetscObject)pc,(PetscObject)jac->kspschur);CHKERRQ(ierr); ierr = PetscObjectIncrementTabLevel((PetscObject)jac->kspschur,(PetscObject)pc,1);CHKERRQ(ierr); + { + PC pcschur; + ierr = KSPGetPC(jac->kspschur, &pcschur); CHKERRQ(ierr); + ierr = PetscObjectIncrementTabLevel((PetscObject)pcschur,(PetscObject)pc,1);CHKERRQ(ierr); + } ierr = KSPCreate(((PetscObject)pc)->comm,&ilink->ksp);CHKERRQ(ierr); ierr = PetscObjectIncrementTabLevel((PetscObject)ilink->ksp,(PetscObject)pc,1);CHKERRQ(ierr); + { + PC ilinkpc; + ierr = KSPGetPC(ilink->ksp, &ilinkpc); CHKERRQ(ierr); + ierr = PetscObjectIncrementTabLevel((PetscObject)ilinkpc,(PetscObject)pc,1);CHKERRQ(ierr); + } This is the pattern used everywhere else in PETSc so that sub-pcs always have the correct tab level.

1 0

Re: [petsc-dev] FieldSplit fixes in 3.3
by Dmitry Karpeev 20 Aug '12

20 Aug '12

On Sun, Aug 19, 2012 at 1:17 PM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: > On Sun, Aug 19, 2012 at 12:48 PM, Dmitry Karpeev <karpeev(a)mcs.anl.gov>wrote: > >> The main reason for this being a patch in 3.3 is that recursive >> FieldSplit is broken there (which is why I couldn't enable it in >> libMesh/Moose correctly). > > > The current interface doesn't explicitly support all the ways to get > information into the split, but several people have worked around the > limitations to get the necessary information in there. Having you and Matt > changing things to rely on mutually incompatible side-effects does not > help. > I don't think the present fix is about the API or working through side effects (nullspace stuff excluded and partially backed out here: http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/1723d4624521) The current interface assumes that the A00 solver and the Schur inner solver are to be set up identically (something we relaxed in petsc-dev). In order to use fieldsplit on these solvers the DM for the corresponding split must be forwarded to both of these A00 solvers. It was being set only on the "outer" A00 solver. Partly this is an artifact of our "split" setup process: some of it occurs in PCFieldSplitSetDefaults(), some later in PCSetUp_FieldSplit(). Some objects obtained in PCFieldSplitSetDefaults() need to be cached until PCSetUp_FieldSplit(), and the splits' DMs are among those. I now put in proper reference counting for them here: http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/188af9799779 > Maybe we can reach some agreement on the proper way to do things. > That's fine with me. The reason I wanted to fix this particular problem now is that some dependent packages (e.g., Moose) only rely on release versions of petsc and would not be able to use this functionality properly for quite some time. As far as I see it, there are other problems with fieldsplit that can (and should) wait for petsc-dev to be resolved (e.g., the proper handling PCReset()). > > More comments below. >> >> On Sun, Aug 19, 2012 at 9:32 AM, Jed Brown <jedbrown(a)mcs.anl.gov> wrote: >> >>> On Sun, Aug 19, 2012 at 2:31 AM, Dmitry Karpeev <karpeev(a)mcs.anl.gov>wrote: >>> >>>> I just pushed a couple of fixes to PCFIELDSPLIT in petsc-3.3 ( >>>> http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/943dde820f7f, >>> >>> >>> Yikes: >>> >>> --- a/src/ksp/pc/impls/fieldsplit/fieldsplit.c >>> +++ b/src/ksp/pc/impls/fieldsplit/fieldsplit.c >>> @@ -14,6 +14,7 @@ >>> PetscInt *fields,*fields_col; >>> VecScatter sctx; >>> IS is,is_col; >>> + DM dm; >>> PC_FieldSplitLink next,previous; >>> }; >>> >>> @@ -263,6 +264,13 @@ >>> if(dms) { >>> ierr = PetscInfo(pc,"Setting up physics based fieldsplit >>> preconditioner using the embedded DM\n");CHKERRQ(ierr); >>> for (ilink=jac->head,i=0; ilink; ilink=ilink->next,i++) { >>> + /* >>> + HACK: >>> + keep a handle to DM here without increfing it: may need the >>> DM to set up the Schur solvers; >>> + can't rely on KSPGetDM() since it will create a DMShell if >>> none was set; >>> + can't use ilink->ksp->dm without creating a dependence on >>> dmimpl.h. >>> + */ >>> + ilink->dm = dms[i]; >>> >>> I think your comment was intended to say "kspimpl.h", but you are >>> creating a tangled web of barbed wire here. You are supposed to call >>> KSPGetDM() and switch based on its *capability*, not *presence*. PETSc >>> should never do something different based on existence of a DM. >>> >> I need to pass the DM from the A00 split to the Schur solvers, but not if >> it's null. If it was PETSC_NULL in the A00 split, >> KSPGetDM() will give me a DMShell, which will generate an error down the >> line. >> > > Let's fix that "down the line". Setting a DM should never force it to be > used or cause an error due to unsupported operation in a case where not > having a DM is also acceptable. > Currently the inner and outer A00 solvers are essentially identified (duplicated), so they should be set up identically, including the DM. Especially when that DM defines a recursive split. This is in part because there is no way to configure those two differently programmatically or from the command line. So if the outer A00 is using -fieldsplit_0_pc_type_fieldsplit, the inner A00 will as well. However, the DM on the inner A00 will be rather different (or absent) and produce splits incompatible with the rest of the options (e.g., a single default split, while -fieldsplit_0_pc_fieldsplit_type schur). Incidentally, this raises this question down the road: in petsc-dev the inner and outer A00 can be configure separately using different prefixes, but how should the inner A00 DM be configured? > > >> I can't reliably get the DM from the split's KSP, >> so I have to stash it somewhere. Otherwise recursive splits don't work. >> If there is a better solution, I'm all for it. >> >>> >>> - /* set tabbing and options prefix of KSP inside the MatSchur */ >>> + /* set tabbing, options prefix and DM of KSP inside the MatSchur >>> */ >>> ierr = MatSchurComplementGetKSP(jac->schur,&ksp);CHKERRQ(ierr); >>> ierr = >>> PetscObjectIncrementTabLevel((PetscObject)ksp,(PetscObject)pc,2);CHKERRQ(ierr); >>> + { >>> + PC pcinner; >>> + ierr = KSPGetPC(ksp, &pcinner); CHKERRQ(ierr); >>> + ierr = >>> PetscObjectIncrementTabLevel((PetscObject)pcinner,(PetscObject)pc,2);CHKERRQ(ierr); >>> + } >>> >>> There is KSPIncrementTabLevel() for a reason. >>> >> In petsc-dev, but not in petsc-3.3. >> >>> >>> - ierr = PetscSNPrintf(schurprefix,sizeof >>> schurprefix,"%sfieldsplit_%s_",((PetscObject)pc)->prefix?((PetscObject)pc)->prefix:"",ilink->splitname);CHKERRQ(ierr); >>> - ierr = >>> KSPSetOptionsPrefix(jac->kspschur,schurprefix);CHKERRQ(ierr); >>> /* really want setfromoptions called in >>> PCSetFromOptions_FieldSplit(), but it is not ready yet */ >>> /* need to call this every time, since the jac->kspschur is >>> freshly created, otherwise its options never get set */ >>> - ierr = KSPSetFromOptions(jac->kspschur);CHKERRQ(ierr); >>> + ierr = PetscSNPrintf(schurprefix,sizeof >>> schurprefix,"%sfieldsplit_%s_",((PetscObject)pc)->prefix?((PetscObject)pc)->prefix:"",ilink->splitname);CHKERRQ(ierr); >>> + ierr = >>> KSPSetOptionsPrefix(jac->kspschur,schurprefix);CHKERRQ(ierr); >>> + ierr = KSPSetFromOptions(jac->kspschur); CHKERRQ(ierr); >>> + ierr = KSPSetUp(jac->kspschur); CHKERRQ(ierr); >>> >>> Why this KSPSetUp()? >>> >> I guess this KSPSetUp() can go. >> >>> You have changed the formatting here, but not explained why. (Also, >>> can't you follow PETSc coding guidelines regarding placement of CHKERRQ?) >>> >>> - } else { >>> - /* set up the individual PCs */ >>> - i = 0; >>> - ilink = jac->head; >>> - while (ilink) { >>> - ierr = >>> KSPSetOperators(ilink->ksp,jac->mat[i],jac->pmat[i],flag);CHKERRQ(ierr); >>> - /* really want setfromoptions called in >>> PCSetFromOptions_FieldSplit(), but it is not ready yet */ >>> - if (!jac->suboptionsset) {ierr = >>> KSPSetFromOptions(ilink->ksp);CHKERRQ(ierr);} >>> - ierr = KSPSetUp(ilink->ksp);CHKERRQ(ierr); >>> - i++; >>> - ilink = ilink->next; >>> - } >>> + } >>> + /* set up the individual PCs */ >>> + i = 0; >>> + ilink = jac->head; >>> + while (ilink) { >>> + ierr = KSPSetOperators(ilink->ksp,jac->mat[i],jac->pmat[i],flag); >>> CHKERRQ(ierr); >>> + /* really want setfromoptions called in >>> PCSetFromOptions_FieldSplit(), but it is not ready yet */ >>> + if (!jac->suboptionsset) {ierr = KSPSetFromOptions(ilink->ksp); >>> CHKERRQ(ierr);} >>> + ierr = KSPSetUp(ilink->ksp);CHKERRQ(ierr); >>> + i++; >>> + ilink = ilink->next; >>> >>> Why do you want this stuff done even for Schur? At least one of these >>> blocks is not used at all, so why are we setting it up? >>> >> In the previous version NONE of these were done for Schur. I suppose we >> can handle these cases more carefully. >> > > You haven't explained why that was wrong. You are certainly setting up > more than necessary here. > I backed out all of KSPSetUp() calls -- they occur inside KSPSolve() anyway, and only when necessary. KSPSetFromOptions() also occurs more sparingly now: http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/3c0d43fb5911 > > >>> >>>> http://petsc.cs.iit.edu/petsc/releases/petsc-3.3/rev/5a9ebf885615) >>>> >>> >>> /* HACK: Check for the constant null space */ >>> - ierr = MatGetNullSpace(pc->pmat, &sp);CHKERRQ(ierr); >>> + /* >>> + This replaces a previous test that checked whether the split's >>> constants are injected into the ambient problem's nullspace. >>> + That can create false negatives (declaring constants not in the >>> nullspace, when they are) as in Stokes or A = [1, 1; 1, 0]. >>> + */ >>> if (sp) { >>> >>> So now you have stale comments and conditionals based on uninitialized >>> memory (sp)? I suppose your intent is to always do this constant null space >>> check regardless of whether the larger matrix has a null space set. I agree >>> that the old code was incorrect because the restriction of the null space >>> of the parent matrix is not comparable to the null space of the submatrix. >>> >>> ierr = MatGetNullSpace(jac->pmat[1], &sp);CHKERRQ(ierr); >>> - if (sp) {ierr = MatSetNullSpace(jac->schur, sp);CHKERRQ(ierr);} >>> - /* set tabbing, options prefix and DM of KSP inside the MatSchur >>> */ >>> + /* >>> + Do we really want to attach the A11-block's nullspace to S? Note >>> that this may generate >>> + "false positives", when the A11's nullspace isn't S's: Stokes or >>> A = [1, 1; 1, 0]. >>> + Using a false nullspace may prevent KSP(S) from converging, >>> since it might force an inconsistent rhs. >>> + */ >>> + /* if (sp) {ierr = MatSetNullSpace(jac->schur, >>> sp);CHKERRQ(ierr);} */ >>> + /* set tabbing, options prefix and DM of KSP inside the >>> MatSchurComplement */ >>> >>> I agree with this change, though the comment is incorrect. The typical >>> problem is that S is nonsingular or has a small null space while A11 is has >>> huge null space (e.g. the zero matrix). The true null space projection for >>> the zero matrix prevents progress in any direction. >>> >>> >>> I think we can repair the second patch, though we really need a general >>> mechanism for providing Schur complement null spaces. Matt advocates >>> attaching it to A11 (where the user knows full well that it's not the null >>> space of A11 when they set it), making the old code "correct". >>> >> Okay, I'll back out that part of the second patch, since A11's nullspace >> is assumed to be meant for S. I'm not sure what's incorrect about the >> comment, though: if a vector is in the A11's kernel, but not in S's, I call >> it a "false positive". This terminology may be wanting, but what's >> incorrect about it? >> > > "force an inconsistent rhs" > The KSP solve will fail only if the search space X is contracted to the point where the rhs b is not in the range AX, hence, inconsistent. I backed out this hunk, however. > > S may well be nonsingular. > > >> >> >>> I think that's confusing and we need a specific API for it, but if >>> people are using it that way, we shouldn't change it in 3.3. >>> >>> I'd prefer to backout the first patch, split it into separate pieces, >>> and review each before pushing. I think it's a potentially big behavior >>> change for 3.3. >>> >> Which pieces of the first patch do you think are too big for 3.3? The >> splits' and Schur KSPs have to be set up (maybe not all in every situation) >> in order for recursive splitting to work. Should we declare that >> capability unavailable for 3.3? I can eliminate the unneeded KSPSetUp() >> calls. >> > > 1. I hate those hacky dm reference-non-references. > > 2. I don't want unused KSPs to be set up. (We should eventually fix the > model so they don't exist...) > > 3. Two of the three IncrementTabLevels that you introduced are not > necessary because the KSP's tab level was incremented before getting the > PC. > Where? > This unnecessary code is ugly and confusing. > The current code is equivalent to KSPIncrementTabLevel(), introduced in petsc-dev for this very reason. My recursive fieldsplit test cases produce what to me is very readable -ksp_monitor output with these tabbing options. It wasn't the case before. Dmitry. > > If we get rid of these three things, the patch reduces to fixing the tab > level of the PC residing inside the MatSchurComplement, which I have no > problem with. >

1 0