[petsc-users] Inverting a Large Matrix with Fortran Code

Jed Brown jed at jedbrown.org
Wed Aug 20 16:24:15 CDT 2014


Koki Imada <koki.imada at york.ac.uk> writes:

> Hello,
>
> I have a problem of the kind:
>
> G(i) = (invM.P).G(i+1) +invM.(d -R.b),
>
> where M, P and R are square tri-diagonal matrices, invM is the inverse of M, G(i) 
> and G(i+1) are the linearised solution vectors (given G(i+1), G(i) is determined 
> using above), and "d" and "b" are vectors.
>
> The problem is, M at grid point (i) depends on M(i-1), so I #cannot# reduce the 
> problem to:
>
> G(i) = invM(i).H(i+1), where H(i+1)=P.G + (d-R.b),
>
> because I first need to calculate alpha(i)=invM(i).P(i) for all i, before 
> calculating G(i) from G(i+1).

I can't tell exactly what you're getting at here, but inv(M)*P is a
dense matrix, so creating it explicitly would turn an O(n) algorithm
into at least O(n^2).

> So I would like to directly invert M, using MatLUFactor and MatMatSolve. The problem 
> here is, my code is written in Fortran 90, and there is a note on MatLUFactor manual 
> page saying:
> "fortran interface is not autogenerated as the f90 interface defintion cannot be 
> generated correctly [due to MatFactorInfo] "
>
> Does this mean the MatLUFactor routine cannot be used in Fortran? Or is there a 
> workaround?

It's called a "Developer Note" for a reason.  The Fortran interface
exists, but it is custom instead auto-generated.  That does not concern
users.
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 818 bytes
Desc: not available
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140820/0759ac96/attachment.pgp>


More information about the petsc-users mailing list