On Mon, Aug 20, 2012 at 4:02 AM, Dmitry Karpeev <[email protected]> wrote:
On Sun, Aug 19, 2012 at 1:17 PM, Jed Brown <[email protected]> wrote:
On Sun, Aug 19, 2012 at 12:48 PM, Dmitry Karpeev <[email protected]>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.