gistfile1.txt

namespace dspsv_func{

const int PRINT_DEBUG_STUFF = 0;

function dspsv_exec(var A,var b,int uplo)
{

	var pla = loadlib("linalg_group");
        var lapack = loadlib('liblapack.so');
        if (lapack == null || !lapack)
                die("Cannot find liblapack. Search for the correct paths");
        
	int m,n,nrhs,ldb,info;
	m=A.rows;
	n=A.cols;
	
	nrhs=b.cols;
	int ipiv[n];
	int j;
	for(j=0;j<n;++j)
		ipiv[j]=0;
	ldb=max(1,n);
say(b);
	var AP;
	AP=create_ap(A,uplo);
say(AP);

	var dspsv = dlfunc(lapack, "dspsv_", "vpppppppp");
        if(dspsv == null || !dspsv)
                die("Not DSPSV");
	
	dspsv(uplo,n,nrhs,AP,ipiv,b,ldb,info);
	
	if(info==0)
		return AP;
}

function create_ap(var A,var uplo)
{
	int m,n;
	m=A.rows;
        n=A.cols;
	var pla = loadlib("linalg_group");
	int ap_size;
        ap_size=((n*(n+1))/2);
	var ap=new 'NumMatrix2D';
	ap.resize(1,ap_size);
	int i,j;
	if(uplo == 85)
	{
		for(i=1;i<=m;++i)
			for(j=1;j<=n;++j)
				ap[(i + (j-1)*j/2)-1] = A[i-1,j-1];
	}
	else if(uplo == 76)
	{
		for(i=1;i<=m;++i)
			for(j=1;j<=n;++j)
				ap[(i + (j-1)*((2*n)-j)/2)-1] = A[i-1,j-1];
	}
	else
	{
		die("invalid input");
	}
	return ap;
}

function max(var a,var b)
{
        return a>b?a:b;
}

function debug(var matrix, string msg, var args [slurpy])
{
    if (PRINT_DEBUG_STUFF) {
        say(sprintf(msg, args));
        say(matrix);
    }
}

}



function main[main](var args)
{
	var pla = loadlib("linalg_group");
	var a = new 'NumMatrix2D';
    var b=new 'NumMatrix2D';
    var ap=new 'NumMatrix2D';
	int uplo;
	uplo=55;
	a.initialize_from_args(3, 3,
                                1.0, 0.0, 0.0,
                                0.0, 1.0, 0.0,
                                0.0, 0.0, 1.0);
	b.initialize_from_args(3, 1,
				1.0,
				2.0,
				3.0);
        using dspsv_func.dspsv_exec;
        ap=dspsv_exec(a,b,uplo);

	say("successful");
	say("AP=");
	say(ap);
	say("B=");
	say(b);
}