## Wednesday, June 18, 2008

### Invert: An External Matrix Inversion program for GAMS

Matrix Inversion

It is possible to perform a matrix inversion using a standard LP solver. See http://www.amsterdamoptimization.com/pdf/lineq.pdf for a detailed description and examples of how this can be achieved. Sometimes it is quicker however just to call an external program to calculate the inverse of a square matrix. The program uses LAPACK's DGESV subroutine.

Attach:invert.exe

Usage

The command line parameters are explained when running invert without parameters:

``C:\>invert ---------------------------------------- INVERT: matrix inversion Erwin Kalvelagen, Amsterdam Optimization ---------------------------------------- Usage    > invert gdxin i a gdxout inva where    gdxin  : name of gdxfile with matrix    i      : name of set used in matrix    a      : name of 2 dimensional parameter inside gdxin    gdxout : name of gdxfile for results (inverse matrix)    inva   : name of 2 dimensional parameter inside gdxout   Calculates the inverse of a non-singular matrix a(i,j) where i and j are aliased sets. inva will contain the inverse matrix.    C:\>``

Example 1

The input looks like:

``\$ontext    Finds the inverse of a matrix through an external program    Erwin Kalvelagen, march 2005    Reference: model gauss.gms from the model library              http://www.gams.com/modlib/libhtml/gauss.htm \$offtext set i  /i1*i3 /;alias (i,j); table a(i,j) 'original matrix'      i1     i2     i3i1    1      2      3i2    1      3      4i3    1      4      3; parameter inva(i,j) 'inverse of a'; execute_unload 'a.gdx',i,a;execute '=invert.exe a.gdx i a b.gdx inva';execute_load 'b.gdx',inva; display a,inva;``

The program has an updated GDXIO interface that both works with old 22.5 GAMS and earlier and with new GAMS 22.6 and later. The GAMS people broke the GDX interface so these versions are really incompatible. However, some of my clients were using a combination of older and newer GAMS systems, so I had to spend lots of time to write a GDX interface that both work with all these GAMS systems. Here is a log from a run with 22.6:

``--- Job invert1.gms Start 02/22/08 13:29:46GAMS Rev 149  Copyright (C) 1987-2007 GAMS Development. All rights reservedLicensee: Erwin Kalvelagen                               G070509/0001CE-WIN          GAMS Development Corporation                               DC4572--- Starting compilation--- invert1.gms(28) 3 Mb--- Starting execution: elapsed 0:00:00.016--- invert1.gms(25) 4 Mb ---------------------------------------- INVERT: matrix inversion Erwin Kalvelagen, Amsterdam Optimization ---------------------------------------- DLL:GDX Library   Dec 24, 2007 WIN.FR.NA 22.6 239.000.000.vis P3PC GDXin:a.gdx Input set:i Input parameter:a LAPACK DGESV, n=           3 GDXout:b.gdx Output parameter:inva Done--- invert1.gms(26) 4 Mb--- GDXin=c:\tmp\New Folder\b.gdx--- invert1.gms(28) 4 Mb*** Status: Normal completion--- Job invert1.gms Stop 02/22/08 13:29:46 elapsed 0:00:00.063``

Example 2

``\$ontext    Finds the inverse of a matrix through an external program    Erwin Kalvelagen, march 2005    Reference: model gauss.gms from the model library              http://www.gams.com/modlib/libhtml/gauss.htm \$offtext set i  /i1*i5/;alias (i,j); parameter a(i,j) 'Hilbert matrix';a(i,j) = 1/(ord(i)+ord(j)-1); parameter inva(i,j) 'inverse of a'; execute_unload 'a.gdx',i,a;execute '=invert.exe a.gdx i a b.gdx inva';execute_load 'b.gdx',inva; display a,inva;``

The output with GAMS 22.5 is:

``--- Job invert2.gms Start 02/22/08 13:32:23GAMS Rev 148  Copyright (C) 1987-2007 GAMS Development. All rights reservedLicensee: Erwin Kalvelagen                               G070509/0001CE-WIN          GAMS Development Corporation                               DC4572--- Starting compilation--- invert2.gms(24) 3 Mb--- Starting execution--- invert2.gms(21) 4 Mb ---------------------------------------- INVERT: matrix inversion Erwin Kalvelagen, Amsterdam Optimization ---------------------------------------- DLL:_GAMS_GDX_237_2007-01-09 GDXin:a.gdx Input set:i Input parameter:a LAPACK DGESV, n=           5 GDXout:b.gdx Output parameter:inva Done--- invert2.gms(22) 4 Mb--- GDXin=c:\tmp\New Folder\b.gdx--- invert2.gms(24) 4 Mb*** Status: Normal completion--- Job invert2.gms Stop 02/22/08 13:32:23 elapsed 0:00:00.078``

Notes
• For Lapack info see: http://www.netlib.org/lapack/
• The program can efficiently invert matrices up to (1000 x 1000). For bigger matrices, the dense inversion routine becomes slower and memory requirements become large.
• Source: Main program Attach:invert.f90