[petsc-users] why a certain option cannot be used
Dave May
dave.mayhem23 at gmail.com
Wed Jul 16 17:40:20 CDT 2014
You need to specify the fieldsplit preconditioning type. The choices are
-pc_fieldsplit_type <MULTIPLICATIVE> (choose one of) ADDITIVE
MULTIPLICATIVE SYMMETRIC_MULTIPLICATIVE SPECIAL SCHUR
Maybe what you want to run is the following
./ex70 -ksp_view -pc_type fieldsplit -fieldsplit_0_pc_type jacobi
-fieldsplit_1_pc_type jacobi -pc_fieldsplit_type SCHUR
-pc_fieldsplit_schur_precondition user
The option " -fieldsplit_1_user_pc" seems meaningless to me.
I don't see anything in ex70.c which would cause this option to have any
effect. The lines
ierr = PetscOptionsHasName(NULL, "-user_pc", &s.userPC);CHKERRQ(ierr);
ierr = PetscOptionsHasName(NULL, "-user_ksp",
&s.userKSP);CHKERRQ(ierr);
imply that the option -user_pc and -user_ksp are meaningful and modify how
the schur complement preconditioner will be configured. The first option
triggers this setup to occur;
if (s->userPC) {
ierr = PCFieldSplitSchurPrecondition(pc, PC_FIELDSPLIT_SCHUR_PRE_USER,
s->myS);CHKERRQ(ierr);
}
the second option will force this setup to occur
if (s->userKSP) {
ierr = PCSetUp(pc);CHKERRQ(ierr);
ierr = PCFieldSplitGetSubKSP(pc, &n, &subksp);CHKERRQ(ierr);
ierr = KSPSetOperators(subksp[1], s->myS, s->myS,
SAME_PRECONDITIONER);CHKERRQ(ierr);
ierr = PetscFree(subksp);CHKERRQ(ierr);
}
On 17 July 2014 00:24, Sun, Hui <hus003 at ucsd.edu> wrote:
> Thank you Dave. It does give me the other options. However I want to be
> able to use the user defined fieldsplit_1_pc_type. In example ex70.c, there
> is a user defined Schur complement preconditioner. I want to use that. Then
> how should I set the values for the options?
>
> I tried to see which values could be set by typing
>
> ./ex70 -pc_type fieldsplit -help 2>/dev/null | grep -A5 fieldsplit_1_pc
>
> And I get the following output:
>
> -fieldsplit_1_pc_type <bjacobi>: Preconditioner (one of) none jacobi
> pbjacobi bjacobi sor lu shell mg
>
> eisenstat ilu icc cholesky asm gasm ksp composite redundant nn mat
> fieldsplit galerkin exotic hmpi asa cp lsc redistribute svd gamg tfs
> (PCSetType)
>
> -fieldsplit_1_pc_use_amat: <FALSE> use Amat (instead of Pmat) to define
> preconditioner in nested inner solves (PCSetUseAmat)
>
> Block Jacobi options
>
> -fieldsplit_1_pc_bjacobi_blocks <-1>: Total number of blocks
> (PCBJacobiSetTotalBlocks)
>
> Krylov Method (KSP) options
> -------------------------------------------------
>
> -fieldsplit_1_ksp_type <preonly>: Krylov method (one of) cg groppcg
> pipecg cgne nash stcg gltr richardson
>
> chebyshev gmres tcqmr bcgs ibcgs fbcgs fbcgsr bcgsl cgs tfqmr cr
> pipecr lsqr preonly qcg bicg fgmres minres symmlq lgmres lcd gcr pgmres
> specest dgmres (KSPSetType)
>
> -fieldsplit_1_ksp_max_it <10000>: Maximum number of iterations
> (KSPSetTolerances)
> -fieldsplit_1_ksp_rtol <1e-05>: Relative decrease in residual norm
> (KSPSetTolerances)
>
> I cannot find anything here about the user defined pc. Do you have any
> idea what might be the choice?
>
>
>
> ------------------------------
> *From:* Dave May [dave.mayhem23 at gmail.com]
> *Sent:* Wednesday, July 16, 2014 3:11 PM
>
> *To:* Sun, Hui
> *Cc:* Matthew Knepley; petsc-users at mcs.anl.gov
> *Subject:* Re: [petsc-users] why a certain option cannot be used
>
> To see all the fieldsplit options, run with
> ./ex70 -pc_type fieldsplit -help
>
> -help will only display the options relevant to the current configuration
> of your KSP and PC.
>
>
> On 17 July 2014 00:10, Dave May <dave.mayhem23 at gmail.com> wrote:
>
>> Run with
>> -help
>>
>>
>> On 17 July 2014 00:07, Sun, Hui <hus003 at ucsd.edu> wrote:
>>
>>> Thank you Dave, it does solve the problem. Another question is: since
>>> some options like
>>>
>>> -fieldsplit_0_pc_type
>>>
>>> does not show up if I type
>>>
>>> ./ex70 -help
>>>
>>> how then can I have a complete list of all the options?
>>>
>>> Best,
>>> Hui
>>>
>>>
>>> ------------------------------
>>> *From:* Dave May [dave.mayhem23 at gmail.com]
>>> *Sent:* Wednesday, July 16, 2014 2:57 PM
>>> *To:* Sun, Hui
>>> *Cc:* Matthew Knepley; petsc-users at mcs.anl.gov
>>>
>>> *Subject:* Re: [petsc-users] why a certain option cannot be used
>>>
>>> You need to specify on the command line
>>> -pc_type fieldsplit -fieldsplit_0_pc_type xxx -fieldsplit_1_pc_type yyy
>>>
>>> Where xxx and yyy are the desired preconditioners for the (1,1) block
>>> and (2,2) block. Note that the default PC (ilu(0)) applied to both splits
>>> will fail for this example.
>>>
>>> The comments in the code are misleading and incorrect. Without the
>>> option -pc_type fieldsplit, all fieldsplit options have no effect.
>>>
>>> When debugging solvers, always run with - ksp_view (as Matt
>>> recommends)
>>>
>>>
>>>
>>> On Wednesday, 16 July 2014, Sun, Hui <hus003 at ucsd.edu> wrote:
>>>
>>>> No I don't think I'm using fieldsplitting although I want to use it.
>>>> The output of
>>>>
>>>> ./ex70 -nx 32 -ny 48 -fieldsplit_1_pc_type jacobi -ksp_view
>>>>
>>>> is the following:
>>>>
>>>> KSP Object: 1 MPI processes
>>>>
>>>> type: gmres
>>>>
>>>> GMRES: restart=30, using Classical (unmodified) Gram-Schmidt
>>>> Orthogonalization with no iterative refinement
>>>>
>>>> GMRES: happy breakdown tolerance 1e-30
>>>>
>>>> maximum iterations=10000, initial guess is zero
>>>>
>>>> tolerances: relative=1e-05, absolute=1e-50, divergence=10000
>>>>
>>>> left preconditioning
>>>>
>>>> using PRECONDITIONED norm type for convergence test
>>>>
>>>> PC Object: 1 MPI processes
>>>>
>>>> type: none
>>>>
>>>> linear system matrix = precond matrix:
>>>>
>>>> Matrix Object: 1 MPI processes
>>>>
>>>> type: nest
>>>>
>>>> rows=4608, cols=4608
>>>>
>>>> Matrix object:
>>>>
>>>> type=nest, rows=2, cols=2
>>>>
>>>> MatNest structure:
>>>>
>>>> (0,0) : prefix="a00_", type=mpiaij, rows=3072, cols=3072
>>>>
>>>> (0,1) : prefix="a01_", type=mpiaij, rows=3072, cols=1536
>>>>
>>>> (1,0) : prefix="a10_", type=mpiaij, rows=1536, cols=3072
>>>>
>>>> (1,1) : prefix="a11_", type=mpiaij, rows=1536, cols=1536
>>>>
>>>> residual u = 2.6315e-05
>>>>
>>>> residual p = 0.000229872
>>>>
>>>> residual [u,p] = 0.000231373
>>>>
>>>> discretization error u = 0.00633503
>>>>
>>>> discretization error p = 0.121534
>>>>
>>>> discretization error [u,p] = 0.121699
>>>>
>>>> WARNING! There are options you set that were not used!
>>>>
>>>> WARNING! could be spelling mistake, etc!
>>>>
>>>> Option left: name:-fieldsplit_1_pc_type value: jacobi
>>>>
>>>>
>>>> ------------------------------
>>>> *From:* Sun, Hui
>>>> *Sent:* Wednesday, July 16, 2014 2:40 PM
>>>> *To:* Matthew Knepley
>>>> *Cc:* petsc-users at mcs.anl.gov <http://UrlBlockedError.aspx>
>>>> *Subject:* RE: [petsc-users] why a certain option cannot be used
>>>>
>>>> Thank you Matt. I've corrected it to
>>>>
>>>> ./ex70 -nx 32 -ny 48 -fieldsplit_1_pc_type jacobi
>>>>
>>>>
>>>> However the output is still:
>>>>
>>>>
>>>> residual u = 2.6315e-05
>>>>
>>>> residual p = 0.000229872
>>>>
>>>> residual [u,p] = 0.000231373
>>>>
>>>> discretization error u = 0.00633503
>>>>
>>>> discretization error p = 0.121534
>>>>
>>>> discretization error [u,p] = 0.121699
>>>>
>>>> WARNING! There are options you set that were not used!
>>>>
>>>> WARNING! could be spelling mistake, etc!
>>>>
>>>> Option left: name:-fieldsplit_1_pc_type value: jacobi
>>>>
>>>>
>>>> If I type
>>>>
>>>>
>>>> ./ex70 -help | grep fieldsplit_1_
>>>>
>>>>
>>>> the output is that nothing coming out.
>>>>
>>>>
>>>>
>>>>
>>>>
>>>> ------------------------------
>>>> *From:* Matthew Knepley [knepley at gmail.com
>>>> <http://UrlBlockedError.aspx>]
>>>> *Sent:* Wednesday, July 16, 2014 2:32 PM
>>>> *To:* Sun, Hui
>>>> *Cc:* petsc-users at mcs.anl.gov <http://UrlBlockedError.aspx>
>>>> *Subject:* Re: [petsc-users] why a certain option cannot be used
>>>>
>>>> On Wed, Jul 16, 2014 at 4:25 PM, Matthew Knepley <knepley at gmail.com
>>>> <http://UrlBlockedError.aspx>> wrote:
>>>>
>>>>> On Wed, Jul 16, 2014 at 4:22 PM, Sun, Hui <hus003 at ucsd.edu
>>>>> <http://UrlBlockedError.aspx>> wrote:
>>>>>
>>>>>> Thank you Matt for answering, however even with that I don't think
>>>>>> it works. For example, the command
>>>>>>
>>>>>> ./ex70 -nx 32 -ny 48 -fieldsplit_1_type jacobi
>>>>>>
>>>>>
>>>>> You have mistyped. Look carefully at my last message:
>>>>>
>>>>> -fieldsplit_1_pc_type jacobi
>>>>>
>>>>
>>>> I also want to note that you can see all the available options with
>>>> -help. It does produce
>>>> a lot of output, but you can segregate it by prefix ("fieldsplit_1_").
>>>> You can also see the
>>>> prefix of each solver component using -ksp_view.
>>>>
>>>> Matt
>>>>
>>>>
>>>>> Matt
>>>>>
>>>>>
>>>>>> gives me the following output
>>>>>>
>>>>>>
>>>>>> residual u = 2.6315e-05
>>>>>>
>>>>>> residual p = 0.000229872
>>>>>>
>>>>>> residual [u,p] = 0.000231373
>>>>>>
>>>>>> discretization error u = 0.00633503
>>>>>>
>>>>>> discretization error p = 0.121534
>>>>>>
>>>>>> discretization error [u,p] = 0.121699
>>>>>>
>>>>>> WARNING! There are options you set that were not used!
>>>>>>
>>>>>> WARNING! could be spelling mistake, etc!
>>>>>>
>>>>>> Option left: name:-fieldsplit_1_type value: jacobi
>>>>>>
>>>>>>
>>>>>> Furthermore, if I look into the possible options with keyword
>>>>>> fieldsplit involved by the command:
>>>>>>
>>>>>>
>>>>>> ./ex70 -help | grep -A5 fieldsplit
>>>>>>
>>>>>>
>>>>>> I get this output:
>>>>>>
>>>>>>
>>>>>> eisenstat ilu icc cholesky asm gasm ksp composite redundant nn mat
>>>>>> fieldsplit galerkin exotic hmpi asa cp lsc redistribute svd gamg tfs
>>>>>> (PCSetType)
>>>>>>
>>>>>> -pc_use_amat: <FALSE> use Amat (instead of Pmat) to define
>>>>>> preconditioner in nested inner solves (PCSetUseAmat)
>>>>>>
>>>>>> Krylov Method (KSP) options
>>>>>> -------------------------------------------------
>>>>>>
>>>>>> -ksp_type <gmres>: Krylov method (one of) cg groppcg pipecg cgne
>>>>>> nash stcg gltr richardson
>>>>>>
>>>>>> chebyshev gmres tcqmr bcgs ibcgs fbcgs fbcgsr bcgsl cgs tfqmr
>>>>>> cr pipecr lsqr preonly qcg bicg fgmres minres symmlq lgmres lcd gcr pgmres
>>>>>> specest dgmres (KSPSetType)
>>>>>>
>>>>>> -ksp_max_it <10000>: Maximum number of iterations (KSPSetTolerances)
>>>>>>
>>>>>>
>>>>>> So I don't see any option that's similar to -fieldsplit_1_type.
>>>>>>
>>>>>>
>>>>>> Hui
>>>>>>
>>>>>>
>>>>>>
>>>>>> ------------------------------
>>>>>> *From:* Matthew Knepley [knepley at gmail.com
>>>>>> <http://UrlBlockedError.aspx>]
>>>>>> *Sent:* Wednesday, July 16, 2014 2:16 PM
>>>>>> *To:* Sun, Hui
>>>>>> *Cc:* petsc-users at mcs.anl.gov <http://UrlBlockedError.aspx>
>>>>>> *Subject:* Re: [petsc-users] why a certain option cannot be used
>>>>>>
>>>>>> On Wed, Jul 16, 2014 at 4:09 PM, Sun, Hui <hus003 at ucsd.edu
>>>>>> <http://UrlBlockedError.aspx>> wrote:
>>>>>>
>>>>>>> I want to solve Stokes equation using user defined schur
>>>>>>> complement preconditioner. So I'm reading and testing
>>>>>>> snes/examples/tutorial/ex55.c, and ex70.c. In these examples, there are
>>>>>>> comments about the usage, for example:
>>>>>>>
>>>>>>> mpiexec -n 2 ./stokes -nx 32 -ny 48 -fieldsplit_1_user_pc.
>>>>>>>
>>>>>>> However the option -fieldsplit_1_user_pc is not recognized by the
>>>>>>> executable. The output of the above command is:
>>>>>>>
>>>>>>
>>>>>> This ("user_pc") is just a mnemonic. What it means is that you can
>>>>>> use any PC in this slot. What he should have written is
>>>>>>
>>>>>> -fieldsplit_1_pc_type <user pc>
>>>>>>
>>>>>> Thanks,
>>>>>>
>>>>>> Matt
>>>>>>
>>>>>>
>>>>>>> residual u = 2.6315e-05
>>>>>>>
>>>>>>> residual p = 0.000229872
>>>>>>>
>>>>>>> residual [u,p] = 0.000231373
>>>>>>>
>>>>>>> discretization error u = 0.00633503
>>>>>>>
>>>>>>> discretization error p = 0.121534
>>>>>>>
>>>>>>> discretization error [u,p] = 0.121699
>>>>>>>
>>>>>>> WARNING! There are options you set that were not used!
>>>>>>>
>>>>>>> WARNING! could be spelling mistake, etc!
>>>>>>>
>>>>>>> Option left: name:-fieldsplit_1_user_pc (no value)
>>>>>>>
>>>>>>>
>>>>>>> I'm wondering what's going on. Is that because I need some
>>>>>>> specific options during my configuring of the package? By the way, the
>>>>>>> version of PETSc I'm using is 3.4.4.
>>>>>>>
>>>>>>>
>>>>>>> Best,
>>>>>>>
>>>>>>> Hui
>>>>>>>
>>>>>>
>>>>>>
>>>>>>
>>>>>> --
>>>>>> 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
>>>>>
>>>>
>>>>
>>>>
>>>> --
>>>> 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
>>>>
>>>
>>
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140717/1cb642ba/attachment-0001.html>
More information about the petsc-users
mailing list