On solving sets of linear equations in Fortran vs Julia

Specific Domains Numerics
, , ,
2

From what you describe it seems indeed that something wrong is happening on the FORTRAN side. However, to rule out that your making an error in the process of extracting data from FORTRAN to Julia I’d suggest to take a look at two options:

  • Compile your FORTRAN code into a shared object and directly call it from Julia with ccall to get the system matrices, the right-hand side and the solution from FORTRAN and then you can check bit by bit from Julia.
  • Perhaps a bit easier: Use an established file format. For MatrixMarket for example there is a FORTRAN implementation and a Julia one and it’s also pretty simple to work with. Maybe that helps you to get some insight?

For using MKL to solve your problem, take a look at MKL.jl … they do all the hard work of rerouting Julia’s default linear algebra operations to the MKL.

With respect to quantifying the quality of a solution, looking at the residual Ax - b as you describe is exactly what I would do.

2 Likes