[Trilinos-Users] don't understand AztecOO with multiple processors
Michael A Heroux
maherou at sandia.gov
Sun Dec 16 23:52:32 MST 2007
Davood,
I don¹t think you want to specify the column map. Let it be done for you by
the matrix class functions. Simply remove that argument from the
Epetra_CrsGraph:
Epetra_CrsGraph Graph(Copy, RowMap, NumNzPtr, true) ;
Also, I would not set the true flag on the final argument until your code is
debug first.
Mike
On 12/16/07 10:53 PM, "Davood Ansari" <david.ansari at gmail.com> wrote:
> Hi All
>
> I have trouble using the AztecOO solver when I use more than one process.
> I think the problem originates form the fact that I don't exactly understand
> the parallel.
> Here is the details:
>
> I have managed to construct my FE matrix and to impose the bc's and loads.
> I have a serial code to which I compare the resulting matrix problem and the
> right hand side vector.
>
> I use Epetra's simple linearly divided maps to construct a RowMap that divided
> my rows
> (or simply DOF's) by the number of processes. At the same time I chose the
> column map to
> be a map that covers all DOFs for each process. Thus for each process is
> associated with
> a row block of the matrix.
> In this way the range map of the matrix is also a map that covers all DOFs for
> each process.
>
> Epetra_Map::Epetra_Map RowMap(newNoOfDOFs, 0, Comm) ;
> Epetra_Map::Epetra_Map ColMap(newNoOfDOFs, newNoOfDOFs, 0, Comm) ;
>
> int* NumNzPtr = (int*)(NumNz) ;
> Epetra_CrsGraph Graph(Copy, RowMap, ColMap, NumNzPtr, true) ;
> build_crs_graph (Graph, NumNz, rank, RowMap, ColMap) ;
>
> Epetra_FECrsMatrix A(Copy, Graph, true) ;
> Epetra_FEVector RHS(RowMap, true) ;
>
>
> For a small example (23 DOFs only) , I did check the matrix and the RHS and
> verified that they are accurate.
>
> Now I am using the AztecOO to solve my sample problem.
> Here I have to define a LHS (or the solutions) vector.
> So I did the following: (A is the matrix in A . LHS = RHS)
>
> Epetra_FEVector LHS(A.DomainMap(), false) ; // this shall contain the
> solution
> Epetra_LinearProblem Problem(&A, &LHS, &RHS) ;
>
>
> AztecOO Solver(Problem) ;
>
> Solver.SetAztecOption( AZ_solver, AZ_gmres) ;
> Solver.SetAztecOption( AZ_precond, AZ_none) ;
> Solver.Iterate(50,1E-9) ;
>
> The code works only when I have a single process. With more than one process,
> lots of
> weired things happen. For example this is the result I get when I use :
>
> cout << LHS ;
>
>
> MyPID GID Value
> 0 0 0
> 0 1 0
> 0 2 0.411514
> 0 3 0
> 0 4 0.104817
> 0 5 0.835807
> 0 6 0
> 0 7 0.332329
> 0 8 0.987688
> 0 9 0
> 0 10 0.500903
> 0 11 0.649448
> 0 12 0
> 0 13 0
> 0 14 0
> 0 15 0
> 0 16 0
> 0 17 0
> 0 18 0
> 0 19 0
> 0 20 0
> 0 21 0
> 0 22 0
> Epetra::MultiVector
> 1 0 0.997976
> 1 1 -0.239904
> 1 2 0
> 1 3 0
> 1 4 0.905788
> 1 5 -0.196422
> 1 6 0.482265
> 1 7 - 0.636227
> 1 8 -0.588409
> 1 9 0.707857
> 1 10 -0.282271
> 1 11 0
> 1 12 0
> 1 13 0
> 1 14 0
> 1 15 0
> 1 16 0
> 1 17 0
> 1 18 0
> 1 19 0
> 1 20 0
> 1 21 0
> 1 22 0
>
>
> While the correct solution is:
>
> 0.000000
> 0.000000
> 0.411514
> 0.000000
> 0.104817
> 0.835807
> 0.000000
> 0.332329
> 0.987688
> 0.000000
> 0.500903
> 0.649448
> 0.996917
> 0.772527
> 0.000000
> 0.000000
> 0.904827
> 0.433849
> 0.481754
> 0.916029
> -0.587785
> 0.707107
> 0.464568
>
>
> The solution seems to have been build as two separate vectors both of the
> global dof's size !!!!
> I guess some thing is terribly wrong with my methodology.
> Kindly comment.
>
> Davood
>
>
>
> _______________________________________________
> Trilinos-Users mailing list
> Trilinos-Users at software.sandia.gov
> http://software.sandia.gov/mailman/listinfo/trilinos-users
-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://software.sandia.gov/mailman/private/trilinos-users/attachments/20071217/aa8676dc/attachment.html
More information about the Trilinos-Users
mailing list