# PETSc input format for linear solvers

I’m going through some considerable effort to translate one of my codes from MATLAB. It’s a type of finite element code and I haven’t implemented the solver yet but comparing CPU times for simply assembling the global system between MATLAB and C, the latter is more than two orders of magnitude faster...

Anyways, I have no idea how to use PETSc but will attempt to learn it in the next few weeks. Just reading over some code and documentation, it seems like it might be pretty hard to translate my code for PETSc.

I noticed that there’s an option to read a matrix from a “file” (see here). I’m just wondering if I go through with translating my code as I’m doing now, can I just save matrix/rhs to a .dat file or something then read from there? I don’t have all that much experience with C so I don’t know if this would be difficult to do or not, especially for larger systems like I’ll have.

In my code, I’ll have to solve a linear system at every Newton iteration (it’s a nonlinear system), and combine this with the fact that the problem is time-dependent. The only time in my entire code that I’ll need PETSc is for the linear solve (admittedly, it’s a lot of linear solves...).

If I’m going to have to learn the syntax for PETSc now, I’d just like to know so I’m not frustrated later.

On another note, I’m really surprised how much faster the code is in C so far. Maybe my MATLAB coding isn’t that efficient, but I’ve been keeping the for loops at a minimum where possible.

I would recommend looking at a time-dependent example because PETSc can provide a lot more diagnostics if you formulate at that level. For example, you could use a Rosenbrock method, an additive Runge-Kutta IMEX method, or others. Some involve Newton iteration, but that is not required and is not always the best approach. To use those methods, you can start by just providing your residual function, $f(t,u)$, as in $\dot u = f(t,u)$, or alternatively, $g$ in $g(t,u,\dot u) = 0$. PETSc can automatically compute a Jacobian via coloring, and if you later decide to assemble the Jacobian for efficiency, PETSc can check whether what you assembled was correct and show you the errors. If you roll your own time integration and your own Newton, these algorithmic variants and debugging features will not be available unless you write them yourself.
Your proposal to use the file system will be a huge bottleneck even in serial. It complicates the workflow greatly and will limit the methods that you can use. But if you just want to test a matrix once, you can use PetscBinaryWrite() from MATLAB and load the system and solve using src/ksp/ksp/examples/tutorials/ex10.c.
• My understanding from your other question is that you have a sparse matrix. Storing that in a dense format turns an $O(n)$ algorithm into $O(n^2)$ storage, which is crippling if you want to solve larger problems. For testing, I recommend using PetscBinaryWrite() to write the sparse MATLAB matrix, plus a vector, to the file. You can use -ksp_view_solution binary:solution.petsc with the example (ex10.c), and read that with PetscBinaryRead() from MATLAB. – Jed Brown Nov 20 '13 at 16:49
• ah, I see. I’m looking through some of the examples now. If I’m set on writing the assembly and such on my own for now, it seems that I could go through my code and translate the types for compatibility? i.e. int to PetscInt, double to PetscScalar, etc.? And same for implicit functions like malloc. I guess there’s nothing to it but to start coding a very small example. Also, I apologize for all these questions that must seem very silly; for someone just learning this gist of PETSc, the documentation is pretty dense. – Justin Dong Nov 20 '13 at 18:25