Dear community,
I am trying to solve linear elasticity equation in strong form using meshless methods. More precisely i use the DCPSE interpolating operators for the approximation of the shape functions and the construction of the square matrices. The method works fine for 2D problems and currently I am trying to proceed in 3D cases.
My limitation is that so far I use direct solvers (UMFPack) and in the 3D problem I will have memory overflow problems so I will need to shift to iterative solvers. My problem is that even for the 2D case I haven't find any iterative solver from the many I tried to reach to a correct solution. I guess it is due to the preconditioning that must be applied to the matrix.
I would appreciate any insight and experience in the solution of this type of matrices. For help I am attaching the sparsity pattern of my matrix for the 2D scenario.
Best regards,
Konstantinos.
There are some things you can try without resorting to iterative methods. You might start with reducing the bandwidth using a reverse cuthill mckee (RCM) technique. If that does not reduce the size of your problem enough to run in memory, you can try adaptive cross approximation (ACA). However, it may take you some time to learn RCM and ACA to the point that you can apply them.
Best of luck.
Dear Gregory,
very interesting your answer indeed. Thank you and i will look into it.
Best regards,
Konstantinos
Konstantinos,
I have a few more thoughts for you regarding ACA and RCM.
The ACA method can determine an effective pre-conditioner for a variety of iterative solvers. I believe it is easier to use ACA in this way than it is to couple it with a hierarchical direct solver.
In using an iterative (in core) method to solve the equation A*X=B for X, the ACA method determines a hierarchical matrix Ah which approximates A. Here, A and Ah are n x n matrices while B and X are n x m matrices.
An effective pre-conditioner M can be made up of the "inadmissible" entries of Ah within a user prescribed band. The entries of M can be stored in memory using the LAPACK banded storage scheme. The matrix M can be LU factored (in place) using a function provided in the Math Kernal Library (MKL). Various iterative methods (e.g. bicg, bicgstab, cgs, gmres, lsqr, qmr) each require the solution to M*x = b and possibly M'*x=b. These may also determined through use of an MKL function. See the ?gbtrf and ?gbtrs sections in the MKL manual (mklman.pdf).
To access the MKL library under Windows, you will need to write a mex interface (mex is a Matlab command). This will create a DLL that can be called just like a Matlab ".m" file, only it runs much faster. Similarly for a Linux OS.
As you can see, implementing all of this will take you some time, not the least of which is understanding how to create the hierarchical matrix Ah and identifying its inadmissible entries.
On the other hand, the RCM algorithm has an easier learning curve, is easier to implement, and does not require a hierarchical direct solver, (just a regular direct solver). RCM should reduce memory usage to some extent. However, the only way I know of to determine whether it is enough for your problem to run in memory, is to try it.
You might also look at the conference Paper: Fast Algebraic Methods in Computational Electromagnetics, which you can download from my ResearchGate page. It will give you some background information on ACA and hierarchical methods along with several references.
Regards, Greg.
Dear Greg,
I apologize for not answering later. I am grateful for the extended explication you provided here and for your remarks. I am looking into both alternatives currently. Your paper is indeed very helpfull. Thank you. I will inform you for any progress on my problem and which method will be more appropriate.
Best regards,
Konstantinos.
- Badges
- Science topic