[petsc-dev] problem when creating a new ksp that calls another one
Raphaël Couturier
raphael.couturier at univ-fcomte.fr
Wed Jul 15 14:50:02 CDT 2015
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
> <mailto: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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20150715/8160339e/attachment.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: tsirm.c
Type: text/x-csrc
Size: 7232 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-dev/attachments/20150715/8160339e/attachment.bin>
More information about the petsc-dev
mailing list