[Trilinos-Users] Problem with AztecOO (with AZ_reuse in nested
loops)
mlx82 at fastwebnet.it
mlx82 at fastwebnet.it
Fri Feb 29 09:03:51 MST 2008
Good evening,
I'm using Trilinos 8.0.5 (without MPI) to solve a Stokes problem. I'm
solving the system with the Schur complement and two nested loops.
As it is a "difficult" problem I need also some preconditioners.
In particoular I'm using two different ILUT preconditioners one for
the
inner loop and the other for the external one.
Everything is already implemented with AztecOO package (and it is
really
farster than IFPACK). I want to compute the preconditioner
for the inner loop only once, and reuse it for every iteration of
external loop.
To be as clear as possible I'm introducing here some notation; suppose
that the linear system is:
K * x = b
where:
K = [A B
C 0]
and:
x = [ V P]^T;
b = [bV bP]^T;
The Schur problem is:
(-C * A^-1 * B) P = bP - C * A^-1 * bV
The resolution of the system is made by two steps:
1) Compute the RHS for the external loop: this is done solving a first
time a linear problem with A matrix and it's ILUT preconditioner (that
is computed now for the first time).
2) Compute P with the two nested loops: this time I want to reuse the
solver (and the preconditioner) used in step 1 for the inner loop. So
I
don't want to recompute the preconditioner.
Unfortunately after the first step I can't reuse the preconditioner of
matrix A in the inner loop. This is the error message that appears:
Error: Did not find previous factorization (requested
by setting options[AZ_pre_calc] to AZ_reuse).
To find this factorization, the following
parameters must match the previous factorization:
1) Total number of unknowns.
2) options[AZ_overlap]
3) options[AZ_scaling]
4) options[AZ_precond]
5) options[AZ_reorder]
6) options[AZ_type_overlap]
7) options[AZ_subdomain_solve]
8) options[AZ_graph_fill]
9) params[AZ_ilut_fill]
10) params[AZ_drop]
11) data_org[AZ_name]
XXXP4302 '*! 1.0277XXX -7926 208
( -7926, 208 32352016) ==> P4302 '*! 1.0277
( -914901, 269872 20123904) ==> vblock in gmres0
( -914901, 1632 19968448) ==> ptrs in gmres0
( -914901, 90760 15778480) ==> general in gmres0
( 2076, 80 32351840) ==> T332 '*! 1.0277
( 2076, 208 32351536) ==> P332 '*! 1.0277
( 2074, 34416 32097392) ==> x_reord P4302 '*! 1.0277
( 2074, 701688 14846224) ==> bindx P4302 '*! 1.0277
( 2074, 1403376 39390720) ==> val P4302 '*! 1.0277
( 2074, 17212 20028912) ==> iu P4302 '*! 1.0277
( 2074, 17212 20011600) ==> ordering P4302 '*! 1.0277
( 2074, 17212 19994288) ==> inv_ordering P4302 '*! 1.0277
( 2074, 224 32349920) ==> A_over P4302 '*! 1.0277
( 2074, 80 32349104) ==> T4302 '*! 1.0277
( 2074, 208 32348800) ==> P4302 '*! 1.0277
I have spent three days tryng to solve this problem without success.
The VERY strange thing is that if I try to solve again a problem with
matrix A before the beginning of step two (and external loop) the code
works, and reuse correctly the preconditioner. Also, if I re-compute
again
the preconditioner(after step one) during the first iteration of
external loop, and than reuse it for all the next iterations (inside
the external loop) the code works again. Unfortunately, with this
latter way I've not only to recompute the same preconditioner two
times, but also I've got a memory leakage! This is certanly introduced
by the
options keep_info = 1 in the internal loop. Valgrind gives me this
output:
==2038== searching for pointers to 9 not-freed blocks.
==2038== checked 3,851,536 bytes.
==2038==
==2038==
==2038== 2,192,376 bytes in 9 blocks are still reachable in loss
record 1 of 1
==2038== at 0x4C21D06: malloc (in /usr/lib64/valgrind/amd64-
linux/vgpreload_memcheck.so)
==2038== by 0x503E9A: AZ_allocate (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x505C20: AZ_manage_memory (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x510626: AZ_mk_context (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x51092B: AZ_oldsolve_setup (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x5117F8: AZ_oldsolve (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x51257A: AZ_iterate (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x4FA178: AztecOO::Iterate(int, double) (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x454FF7: shurMatrix::Apply(Epetra_MultiVector const&,
Epetra_MultiVector&) const (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x4F5FF0: Epetra_Aztec_operatorvec (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x503360: AZ_compute_residual (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038== by 0x5083EF: AZ_scale_true_residual (in
/home/mlx/Steam3D/Steam3D_MLX/bin/Steam3D)
==2038==
==2038== LEAK SUMMARY:
==2038== definitely lost: 0 bytes in 0 blocks.
==2038== possibly lost: 0 bytes in 0 blocks.
==2038== still reachable: 2,192,376 bytes in 9 blocks.
==2038== suppressed: 0 bytes in 0 blocks.
To better understood this leakage I've tryed to solve the same problem
recomputing always the preconditioner and setting AZ_keep_info = 0,
AZ_pre_calc = "calc". And the leakage disapear.
Below I've reported a short description of the code that I'm using
now.
At the beginning I'm setting some parameters for the loops:
//Setting parameters for internal loop
List_Solver_int = new Teuchos::ParameterList("Solver: internal loop");
List_Solver_int->set("solver", Sol_int_method);
List_Solver_int->set("conv", "r0");
List_Solver_int->set("tol", Sol_int_tol);
List_Solver_int->set("max_iter", Sol_int_ITmax);
List_Solver_int->set("kspace", Sol_int_Kril);
List_Solver_int->set("output", "last");
List_Solver_int->set("precond", "dom_decomp");
List_Solver_int->set("reorder", 1);
List_Solver_int->set("subdomain_solve", Prec_type);
List_Solver_int->set("overlap", 0);
List_Solver_int->set("type_overlap", "standard");
List_Solver_int->set("graph_fill", Prec_lof);
List_Solver_int->set("ilut_fill", Prec_lof);
List_Solver_int->set("drop", Prec_dt);
List_Solver_int->set("pre_calc", "calc");
List_Solver_int->set("keep_info",1);
//Setting parameters for external loop
List_Solver_ext = new Teuchos::ParameterList("Solver: external loop");
List_Solver_ext->set("solver", Sol_ext_method);
List_Solver_ext->set("conv", "rhs");
List_Solver_ext->set("tol", Sol_ext_tol);
List_Solver_ext->set("max_iter", Sol_ext_ITmax);
List_Solver_ext->set("kspace", Sol_ext_Kril);
List_Solver_ext->set("output", "none");
List_Solver_ext->set("precond", "dom_decomp");
List_Solver_ext->set("reorder", 1);
List_Solver_ext->set("subdomain_solve", Prec_type);
List_Solver_ext->set("overlap", 0);
List_Solver_ext->set("type_overlap", "standard");
List_Solver_ext->set("graph_fill", Prec_lof);
List_Solver_ext->set("ilut_fill", Prec_lof);
List_Solver_ext->set("drop", Prec_dt);
List_Solver_ext->set("pre_calc", "calc");
List_Solver_ext->set("keep_info", 0);
Than I have a C++ class "invA" to solve the general problem: A * LHS =
RHS. The constructor of this class is:
solver.SetUserMatrix(A, flag);
solver.SetLHS(X);
solver.SetRHS(RHS);
solver.SetParameters(List_Solver_int);
where A, X, RHS are addresses, and solver is a private member (AztecOO
solver).
Therefore defining once the class "invA" I create only one solver for
this kind of problem. The main idea for this class is to change RHS
values outside the class without changing it's address.
To perform the first step I call a member function of "invA" that
does:
solver.Iterate(maxIT,tol);
this compute automatically the preconditioner.
Than I call a member function of "invA" that uses:
solver.SetAztecOption(AZ_pre_calc, AZ_reuse);
so, for the future call of solver.Iterate(maxIT,tol), it should reuse
the preconditioner!
To perform the second step and compute the pressure P, I define an
epetra operator "schurMatrix" in which is present the function
member "apply" that contains the operation that computes (-C * A^-1 *
B) * P,
and pointers to "invA" class and it's member.
I start the external loop (Aztec00 Sys_solver, with List_Solver_ext
options and appropiate LHS and RHS), after setting Sys_solver.
SetPrecMatrix(Mass) where Mass is a Mass matrix for compute external
ILUT preconditioner. Inside apply there is a call to solver.Iterate
(maxIT,tol)
that is inside "invA". And here the program ends with that
error about AZ_reuse.
Thanks for any help.
Cristiano Malossi
More information about the Trilinos-Users
mailing list