[happy_clp]

hi,i have a problem in using petsc to compute 3D FEA.my linear systems always diverge when i use others solution methods except for lu.as you konw, there are many solution methods in petsc, such as bcgs,gmres,cg......,but they really don't adapt to my situation. And when using lu as pc or ksp, memory becomes the bottleneck. anybody can help me? thank you very much!

[mr_fluent]

make sure that what you are solving has diagonal dominance in matrix.
it seems that your matrix does not have this property, LU methods are somewhat less sensitive to this issue.

[happy_clp]

you mean to LU less sensitive to diagonal dominance in matrix? I cannot clearly understand your meaning. i cannot make sure my matrix is diagonal dominance, but i guess it is, beause FEM's matrix is general diagonal dominance.if diagonal dominance, whether or not the solution is easier? or reverse.

[mr_fluent]

Originally Posted by happy_clp

you mean to LU less sensitive to diagonal dominance in matrix? I cannot clearly understand your meaning. i cannot make sure my matrix is diagonal dominance, but i guess it is, beause FEM's matrix is general diagonal dominance.if diagonal dominance, whether or not the solution is easier? or reverse.

the reason i said about diagonal dominance is because you said nothing else is working. If matrix is diagonally dominant, even a normal method like gauss siedel would show some convergence. (though it would be slow for large cases).

[happy_clp]

thank you very much!
my work is about 3D dynamic FEA,and use PETSc to solve Ax=b.but now i find when i use direct method(namely LU),without sufficient memory.When I select iterative method, the result is always diverge,such as sor, cg, gmres etc..so i am confusing.

[mr_fluent]

is there any way for you to print the matrix and then check if ap > sum(ali).

[happy_clp]

Thank you!
I don't care if ap > sum(ali) and i also don't know what is ap or sum(ali)?But,i guess i can do it by matlab.But my current problem is how to solve the scalable sparse matrix?because direct method has not enough compute memory,iterative method can not converge(i tried many iterative methods of PETSc).
whether or not you have some experiences in PETSc?
Thank you again!

[mr_fluent]

Originally Posted by happy_clp

Thank you!
I don't care if ap > sum(ali) and i also don't know what is ap or sum(ali)?But,i guess i can do it by matlab.But my current problem is how to solve the scalable sparse matrix?because direct method has not enough compute memory,iterative method can not converge(i tried many iterative methods of PETSc).
whether or not you have some experiences in PETSc?
Thank you again!

if you have to solve it only once then use lsqr method, it would not diverge no matter what you through at it.

if you have written in little bit more detail about equations you are solving it might be possible to something.

[happy_clp]

Thanks a lot! I will try lsqr tomorrow! Before, I learned bicgstab or GMRES is a better iterative method for asymmetric matrix,but they were useless to my matrix. CG or PCG is very good for symmetric positive definite matrix, however, my matrix is symmetric and indefinite.Today, to my surprised, I found each diagonal entity of my matrix approximately equals to absolute value of sum of each row all entities. The determinant is infinite, and condition also is very great.Now, I think the best method is find for a preconditioner.Do you think? thank you very much!thank your help. you are welcome if you need more matrix information!

[happy_clp]

lsqr is very slower, iteration is up to 70000.

[mr_fluent]

this is why use it only if you wish to solve it once. it you have to solve it again and again.

[happy_clp]

can you speak chinese?

[mr_fluent]

Originally Posted by happy_clp

can you speak chinese?

i can only speak, english, hindi, portuguese, japanese.

i can not speak chinise, though after few years when i am proficient in japanese i have desire to learn it.

[jed]

Originally Posted by happy_clp

Thanks a lot! I will try lsqr tomorrow! Before, I learned bicgstab or GMRES is a better iterative method for asymmetric matrix,but they were useless to my matrix. CG or PCG is very good for symmetric positive definite matrix, however, my matrix is symmetric and indefinite.Today, to my surprised, I found each diagonal entity of my matrix approximately equals to absolute value of sum of each row all entities. The determinant is infinite, and condition also is very great.Now, I think the best method is find for a preconditioner.

It seems like you are confusing the preconditioner with the Krylov method. As the first step, you should find a preconditioner that works with GMRES and

-ksp_gmres_restart 300

(or some other large number, to eliminate restarts as a reason for not converging).

Since your matrix is indefinite, it's quite likely that no black-box preconditioner will work well. One method that is fairly reliable if your matrix isn't terribly conditioned is sparse approximate inverse. You can select this with

-pc_type hypre -pc_hypre_type parasails

if you built PETSc with Hypre or

-pc_type spai

if built with SPAI. Unfortunately, the inverse of anything interesting is dense, so sparse approximations are usually quite poor.

Another poor-man's solution is to use a Schwarz preconditioner with exact (LU) solves in each subdomain. The theory doesn't guarantee that this will converge, but in practice it often works, although usually with suboptimal convergence rates. Use

-pc_type asm -sub_pc_type lu

and either run in parallel or choose

-pc_asm_blocks

large enough that you can afford the memory for direct factorization. Note that sometimes third-party direct solvers will be better, you should probably at least try

-sub_pc_factor_mat_solver_package mumps

and

-sub_pc_factor_mat_solver_package umfpack

To do better, you will need to use problem-specific information. Such methods are based on block-LDU factorization and the crucial part is finding a preconditioner for the Schur complement which is necessarily a dense matrix (hence unformable). For incompressible flow, one method is obtained by approximately commuting operators in the Schur complement, producing a shifted advection-diffusion operator defined in the pressure space. The other good method for preconditioning the Schur complement in incompressible flow is known as the least-squares commutator (LSC). A special case of LSC which ignores the effect of the mass matrix is known as BFBt and approximates the inverse of the Schur complement using Moore-Penrose pseudoinverses. In PETSc-3, you should look at PCBFBT (it requires you to provide the necessary structure). More generally, you can use PCFIELDSPLIT, but you will have to handle the approximate inverse of the Schur complement yourself. For more difficult problems, I just write my own shell preconditioners.

A good general reference for indefinite problems is

@article{benzi2005nss,
  title={{Numerical solution of saddle point problems}},
  author={Benzi, M. and Golub, G.H. and Liesen, J.},
  journal={Acta Numerica},
  volume={14},
  pages={1--137},
  year={2005},
  publisher={Cambridge Univ Press}
}

For incompressible flow, this recent paper is quite good

@article{elman2008tcp,
  title={{A taxonomy and comparison of parallel block multi-level preconditioners for the incompressible Navier-Stokes equations}},
  author={Elman, H.C. and Howle, V.E. and Shadid, J. and Shuttleworth, R. and Tuminaro, R.},
  journal={Journal of Computational Physics},
  volume={227},
  number={1},
  pages={1790--1808},
  year={2008},
  publisher={Academic Press}
}

[jed]

Also, the PETSc mailing list http://lists.mcs.anl.gov/pipermail/petsc-users/ is a very good place to ask questions like this. Provide as much detail as possible and exact error messages.

[happy_clp]

Thank you very much!
GMRES and other iterative methods all cann't obtain convergence to my problem, unless select lu as a pc or only use lu as a ksp,however, memory becomes new bottleneck again. PETSc team recommends MUMPS,but I find MUMPS cannot be compatible with fortran 77. My program need compute eigenvalues, so I utilize arpack package and it is only compatible with fortran 77.
I also want to use hypre,but I don't know how to use? Whether or not after I install hypre, I can utilize hypre by -pc_type hypre -pc_hypre_type parasails at once or I need modify my program yet? I am short of relevant
informations and experiences.Can you help me? My program is about vibration analysis of anisotropic material.I compute eigenvalue by arpack package.But in arpack, I need utilize own program to compute Ax=b,so I select PETSc package to do it.After several attempts, I only find, if convergence, ksp or pc must include lu.

[jed]

Originally Posted by happy_clp

Thank you very much!
GMRES and other iterative methods all cann't obtain convergence to my problem, unless select lu as a pc or only use lu as a ksp,however, memory becomes new bottleneck again. PETSc team recommends MUMPS,but I find MUMPS cannot be compatible with fortran 77. My program need compute eigenvalues, so I utilize arpack package and it is only compatible with fortran 77.
I also want to use hypre,but I don't know how to use? Whether or not after I install hypre, I can utilize hypre by -pc_type hypre -pc_hypre_type parasails at once or I need modify my program yet? I am short of relevant
informations and experiences.Can you help me? My program is about vibration analysis of anisotropic material.I compute eigenvalue by arpack package.But in arpack, I need utilize own program to compute Ax=b,so I select PETSc package to do it.After several attempts, I only find, if convergence, ksp or pc must include lu.

To repeat, your lack of convergence is not a failing of GMRES or any other Krylov iteration, it's use of an inappropriate preconditioner. There are no good black-box preconditioners for indefinite problems, sparse approximate inverse is generally terrible, it's just often better than nothing. To do better, you really need to use problem-specific information. Read the paper by Benzi and others for some general guidelines on this.

After building PETSc with Hypre, you can control it on the command line; no code modification is necessary. Run with -help to see what options are available with the current solver combination (piping through grep is helpful).

A more flexible and scalable eigensolver is SLEPc. It uses PETSc data structures and brings the same flexibility to eigenproblems. It provides a uniform interface to several third-party solvers including ARPACK.