### script to invert a matrix (stored in GDX) in R
if (! require(gdxrrw)) stop ("gdxrrw package is not available")
## check for and use the optional GAMS sysdir arg
ca <- commandArgs(trailingOnly=TRUE)
if (length(ca) > 0) {
print (paste('GAMS sysdir passed to invert1.r:',ca[1]))
igdx(ca[1])
}
if (0 == igdx(silent=TRUE)) stop ("the gdx shared library has not been loaded")
fnIn <- 'a.gdx'
fnOut <- 'c.gdx'
# we want to read the matrix as square and dense, essentially as A(i,i),
# so we filter both dimensions. This way, we get exactly the structure we want
# regardless of the universe in the GDX or zero rows/cols in the matrix
i <- rgdx(fnIn,list(name='i'))
uu <- list(i$uels[[1]],i$uels[[1]])
a <- rgdx(fnIn,list(name='a',form='full',uels=uu))
if (a$dim != 2) stop ('matrix read in from GDX should have dimension 2')
# inverse matrix is very similar to A: same size, etc
ainv <- a
ainv$val <- solve(a$val)
ainv$ts <- 'inverse of a computed in R via solve(a)'
ainv$name <- 'inva'
# introduce some noise to exercise the check for a bad inverse
# ainv$val[1,1] <- ainv$val[1,1] * (1 + 1e-8)
wgdx(fnOut, ainv)
print (paste('R script all done: inverse of A written to',fnOut))