In my work, it’s more convenient to export matrix and RHS to PETSc to solve Ax=b. Because I want to develop a special efficient solving method for the problem based on PETSc.
I want to export a PETSc VECMPI type from FreeFem++, to use “VecLoad” to load a vector(eg. RHS) in the following work.
I try to export RHS according to the example, but the type is matseqdense or matmpidense.
real[int] bb;
real[int, int] rhsView(bb.n, 1);
rhsView(:,0) = bb;
ObjectView(rhsView, format = "binary", name = "./data/bb.dat");
You should be able to “develop a special efficient solving method for the problem based on PETSc” directly in FreeFEM, what is missing?
What you are doing is correct. In your code, instead of using the loaded MatSeqDense or MatMPIDense, just use MatDenseGetColumnVec — PETSc v3.17.1-371-gabf86ae5d7 documentation to get the corresponding VecSeq or VecMPI.
A few weeks ago I had a similar problem as @DuanXY and solved it as you suggested using MatDenseGetColumnVec.
In my case it was about testing the feasibility to solve a certain system with AmgX, which is unfortunately not included in PETSc, but there is a wrapper including a solveFromFiles example.
So maybe it is worth it to consider saving vectors in the standard binary vector storage format, because at least two people have tripped over this?
I forgot that you could also use a dummy system, e.g., with the identity matrix, call KSPSolve(), and use -ksp_view_rhs binary. This will get the job done as well.
The tricky part about saving a vector is that users have to be careful whether they save it in FreeFEM or PETSc numbering, but that could be doable (feel free to make a PR over at GitHub).
Awesome. I was unaware of that merge request. Thanks for the hints.
I am not sure I can do this myself, as I am a novice when it comes to PETSc and the freefem plugin
I guess that I would have to add Global.Add("ObjectView", "(", new PETSc::view< KN<PetscScalar>, 0 >); Global.Add("MatLoad", "(", new PETSc::view< KN<PetscScalar>, 1 >); and a third implementation of view_dispatched?
Doesn’t the tricky part you mentioned also apply to the existing implementation? If I save a real[int] rhs using a real[int, int] rhsView(rhs.n, 1), I also have to pay attention to the numbering and use ChangeNumbering if necessary, or am I missing something?
Doesn’t the tricky part you mentioned also apply to the existing implementation?
That is correct, but most PETSc routine that uses a PetscScalar[int, int] as input expects this input to be in PETSc numbering. On the contrary, if you use ^-1 * in or KSPSolve(, in, ), then in must either be in FreeFEM or PETSc numbering.
