[petsc-dev] problem when creating a new ksp that calls another one
Matthew Knepley
knepley at gmail.com
Wed Jul 15 15:05:36 CDT 2015
On Wed, Jul 15, 2015 at 3:03 PM, Barry Smith <bsmith at mcs.anl.gov> wrote:
>
> You cannot have your KSPSolve_TSIRM() change KSP options on ITS KSP and
> then call KSPSetUp() on its KSP then call KSPSolve() on its KSP. If you
> want to have an "inner" solver you need to have a different KSP for that
> inner solver.
>
Also, looking at your code, why can't you get rid of all that FGMRES code
and just use PCKSP inside?
Thanks,
Matt
> Barry
>
>
>
> > On Jul 15, 2015, at 2:50 PM, Raphaël Couturier <
> raphael.couturier at univ-fcomte.fr> wrote:
> >
> > Hello,
> >
> > Sorry for my bad example...
> >
> > So maybe it is simpler I give you my complete code called TSIRM for A
> Two-Stage Iteration with least-squares Residual Minimization algorithm to
> solve large sparse linear systems (
> https://www.researchgate.net/publication/277571149_TSIRM_A_Two-Stage_Iteration_with_least-squares_Residual_Minimization_algorithm_to_solve_large_sparse_linear_systems
> )
> > It is still under in development as a solver for PETSc (with few
> comments in the code...)
> >
> > In this code, I use FGMRES as an inner solver only for 30 iterations (I
> want to put this as a parameter after)
> > For the outer iteration I apply a minimization method on the residual
> each 12 outer iterations (this is the parameter colS in the code, this will
> also be a parameter after). The minimization is CGLS conjugate gradient
> with least square.
> > There are useless messages when running the code (but they are
> interesting to understand the behavior of the code)
> >
> > Some examples:
> >
> > With SNES ex35
> >
> > With FGMRES
> > time $PETSC_DIR/arch-linux2-c-debug/bin/mpirun -np 3 ./ex35 -snes_rtol
> 1.e-12 -snes_monitor -ksp_type fgmres -da_grid_x 660 -da_grid_y 660
> > 0 SNES Function norm 3.808595655784e+02
> > 1 SNES Function norm 3.785299355288e-03
> > 2 SNES Function norm 3.766792406275e-08
> > 3 SNES Function norm 1.428490284661e-09
> >
> > real 2m35.576s
> > user 7m44.008s
> > sys 0m2.504s
> >
> >
> > With TSIRM
> > time $PETSC_DIR/arch-linux2-c-debug/bin/mpirun -np 3 ./ex35 -snes_rtol
> 1.e-12 -snes_monitor -ksp_type tsirm -da_grid_x 660 -da_grid_y 660
> > 0 SNES Function norm 3.808595655784e+02
> > Enter tsirm
> > Size of vector 435600
> > converged 0
> > Norm of error 270.484, outer iteration 0 270.484 reason 0
> > Norm of error 236.805, outer iteration 0 236.805 reason 0
> > Norm of error 213.261, outer iteration 0 213.261 reason 0
> > ....
> > Norm of error 5.84204e-13, outer iteration 3 5.84204e-13 reason 0
> > Norm of error 4.61092e-13, outer iteration 3 4.61092e-13 reason 0
> > Norm of error 3.74294e-13, outer iteration 3 3.74294e-13 reason 2
> > -- Execution time : 22.2333 (s)
> > -- Total number of iterations : 1199
> > 3 SNES Function norm 1.432052860297e-09
> >
> > real 0m53.639s
> > user 2m39.120s
> > sys 0m1.652s
> >
> >
> >
> >
> >
> > With KSP 15
> >
> > with FGMRES
> >
> > time $PETSC_DIR/arch-linux2-c-debug/bin/mpirun -np 3 ./ex15 -m 800
> -n 800 -ksp_type fgmres -ksp_rtol 1e-6 -pc_type mg
> > Norm of error 1.79776 iterations 1102
> >
> > real 0m51.147s
> > user 2m32.788s
> > sys 0m0.508s
> >
> >
> > With TSIRM
> >
> > time $PETSC_DIR/arch-linux2-c-debug/bin/mpirun -np 3 ./ex15 -m 800
> -n 800 -ksp_type tsirm -ksp_rtol 1e-6 -pc_type mg
> > Enter tsirm
> > Size of vector 640000
> > converged 0
> > Norm of error 0.0970815, outer iteration 0 0.0970815 reason 0
> > Norm of error 0.0368809, outer iteration 0 0.0368809 reason 0
> > Norm of error 0.0223725, outer iteration 0 0.0223725 reason 0
> > Norm of error 0.0163448, outer iteration 0 0.0163448 reason 0
> > Norm of error 0.013248, outer iteration 0 0.013248 reason 0
> > Norm of error 0.0112125, outer iteration 0 0.0112125 reason 0
> > Norm of error 0.00961675, outer iteration 0 0.00961675 reason 0
> > Norm of error 0.00825048, outer iteration 0 0.00825048 reason 0
> > Norm of error 0.00707865, outer iteration 0 0.00707865 reason 0
> > Norm of error 0.00606273, outer iteration 0 0.00606273 reason 0
> > Norm of error 0.00516968, outer iteration 0 0.00516968 reason 0
> > Norm of error 0.00439978, outer iteration 0 0.00439978 reason 0
> > Norm of error 5.22567e-05, outer iteration 1 5.22567e-05 reason 2
> > -- Execution time : 18.8489 (s)
> > -- Total number of iterations : 386
> > Norm of error 0.901491 iterations 26
> >
> > real 0m19.263s
> > user 0m57.248s
> > sys 0m0.400s
> >
> >
> >
> >
> > With KSP ex45 (with two KSPs not as in this code, we obtain faster
> results TSIRM with the ASM preconditioner than FGMRES with HYPRE or with
> ASM)
> > With one KSP (this code)
> > It crashes...
> >
> > time $PETSC_DIR/arch-linux2-c-debug/bin/mpirun -np 3 ./ex45
> -da_grid_x 160 -da_grid_y 160 -da_grid_z 160 -ksp_type tsirm -ksp_rtol
> 1e-10 -pc_type asm -ksp_monitor
> > Enter tsirm
> > Size of vector 4096000
> > [1]PETSC ERROR: --------------------- Error Message
> --------------------------------------------------------------
> > [1]PETSC ERROR: Object is in wrong state
> > [1]PETSC ERROR: Vec is locked read only, argument # 1
> > [1]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html
> for trouble shooting.
> > [1]PETSC ERROR: Petsc Release Version 3.6.0, Jun, 09, 2015
> > [1]PETSC ERROR: ./ex45 on a arch-linux2-c-debug named extinction by
> couturie Wed Jul 15 21:43:51 2015
> > [1]PETSC ERROR: Configure options --download-openmpi --download-hypre
> --download-f2cblaslapack --with-fc=0
> > [1]PETSC ERROR: #1 VecGetArray() line 1646 in
> /home/couturie/work/petsc-3.6.0/src/vec/vec/interface/rvector.c
> > [1]PETSC ERROR: #2 VecGetArray3d() line 2411 in
> /home/couturie/work/petsc-3.6.0/src/vec/vec/interface/rvector.c
> > [1]PETSC ERROR: #3 DMDAVecGetArray() line 77 in
> /home/couturie/work/petsc-3.6.0/src/dm/impls/da/dagetarray.c
> > [1]PETSC ERROR: #4 ComputeRHS() line 84 in
> /home/couturie/work/petsc-3.6.0/src/ksp/ksp/examples/tutorials/ex45.c
> > [1]PETSC ERROR: #5 KSPSetUp() line 277 in
> /home/couturie/work/petsc-3.6.0/src/ksp/ksp/interface/itfunc.c
> > [1]PETSC ERROR: #6 KSPSolve_TSIRM() line 135 in
> /home/couturie/work/petsc-3.6.0/src/ksp/ksp/impls/gmres/tsirm.c
> >
> >
> >
> > Thank you in advance,
> >
> > Raphaël
> >
> >
> >> On Wed, Jul 15, 2015 at 12:04 PM, Raphaël Couturier <
> raphael.couturier at univ-fcomte.fr> wrote:
> >> Hello,
> >>
> >> I want to make a new solver inside petsc.
> >> This solver is based on a 2 stage method. So I have one outer iteration
> and one inner iteration. For the inner iteration I can use gmres or one of
> its variants.
> >> I succeeded to make my algorithm by creating a new ksp (that copies the
> linear system). So I have tried to do that without creating a new ksp. It
> works well with ksp examples that don't use dmda and with snes examples.
> >>
> >> The code you sent does not make any sense to me. How can KSPSolve(ksp,
> NULL, NULL) be meaningful?
> >>
> >> Also, can't you get exactly this effect using PCKSP?
> >>
> >> Thanks,
> >>
> >> Matt
> >>
> >> So I have tried to make a really simple ksp (only to test) in which I
> simply call another ksp.
> >> Here is the part of the code (and attached is the complete code):
> >> PetscErrorCode ierr;
> >> ierr = KSPSetFromOptions(ksp);CHKERRQ(ierr);
> >> KSPSetType(ksp,KSPFGMRES);
> >> ierr = KSPSolve(ksp, NULL, NULL); CHKERRQ(ierr);
> >> PetscFunctionReturn(0);
> >>
> >>
> >> For example with the ksp example 15:
> >> time $PETSC_DIR/arch-linux2-c-debug/bin/mpirun -np 3 ./ex15 -m 200
> -n 200 -ksp_type test -ksp_rtol 1e-10 -pc_type mg -ksp_monitor
> >> there is no problem, I obtain the good result.
> >>
> >> But with the ksp example 45 (based on dmda):
> >> time $PETSC_DIR/arch-linux2-c-debug/bin/mpirun -np 1 ./ex45
> -ksp_type test -ksp_rtol 1e-10 -pc_type mg -ksp_monitor
> >>
> >> It crashes with that (I didn't change anything in ex45)
> >> [0]PETSC ERROR: Object is in wrong state
> >> [0]PETSC ERROR: Vec is locked read only, argument # 1
> >> [0]PETSC ERROR: See http://www.mcs.anl.gov/petsc/documentation/faq.html
> for trouble shooting.
> >> [0]PETSC ERROR: Petsc Release Version 3.6.0, Jun, 09, 2015
> >> [0]PETSC ERROR: ./ex45 on a arch-linux2-c-debug named extinction by
> couturie Wed Jul 15 18:45:39 2015
> >> [0]PETSC ERROR: Configure options --download-openmpi --download-hypre
> --download-f2cblaslapack --with-fc=0
> >> [0]PETSC ERROR: #1 VecGetArray() line 1646 in
> /home/couturie/work/petsc-3.6.0/src/vec/vec/interface/rvector.c
> >> [0]PETSC ERROR: #2 VecGetArray3d() line 2411 in
> /home/couturie/work/petsc-3.6.0/src/vec/vec/interface/rvector.c
> >> [0]PETSC ERROR: #3 DMDAVecGetArray() line 77 in
> /home/couturie/work/petsc-3.6.0/src/dm/impls/da/dagetarray.c
> >> [0]PETSC ERROR: #4 ComputeRHS() line 84 in
> /home/couturie/work/petsc-3.6.0/src/ksp/ksp/examples/tutorials/ex45.c
> >> [0]PETSC ERROR: #5 KSPSetUp() line 277 in
> /home/couturie/work/petsc-3.6.0/src/ksp/ksp/interface/itfunc.c
> >> [0]PETSC ERROR: #6 KSPSolve() line 546 in
> /home/couturie/work/petsc-3.6.0/src/ksp/ksp/interface/itfunc.c
> >> [0]PETSC ERROR: #7 KSPSolve_TEST() line 43 in
> /home/couturie/work/petsc-3.6.0/src/ksp/ksp/impls/gmres/test.c
> >> [0]PETSC ERROR: #8 KSPSolve() line 604 in
> /home/couturie/work/petsc-3.6.0/src/ksp/ksp/interface/itfunc.c
> >> [0]PETSC ERROR: #9 main() line 51 in
> /home/couturie/work/petsc-3.6.0/src/ksp/ksp/examples/tutorials/ex45.c
> >>
> >>
> >> I don't understand since many examples in snes are good with my code
> and they use dmda.
> >> Do any of you have an idea of the problem?
> >>
> >> Of course, I have a bad solution that consists in creating a new ksp
> but I would like to do without, if possible.
> >>
> >> Thank you in advance,
> >> Regards,
> >>
> >> Raphaël Couturier
> >>
> >>
> >>
> >> --
> >> 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
> >
> > <tsirm.c>
>
>
--
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-dev/attachments/20150715/bffd5780/attachment.html>
More information about the petsc-dev
mailing list