real vector range1(real scalar n, real scalar i) {
	if (i < 1 | i > n) {
		return(1::n)
	} else if (i == 1) {
		return(2::n)
	} else if (i == n) {
		return(1::n-1)
	} else {
		return(1::i-1\i+1::n)
	}
}

real matrix submat(real matrix a, real scalar i, real scalar j) {
	return(a[range1(rows(a), i), range1(cols(a), j)])
}

real scalar sumrec(real matrix a, real scalar x) {
	real scalar n, s, p
	n = rows(a)
	if (n==1) return(a[1,1])
	s = 0
	p = 1
	for (i=1; i<=n; i++) {
		s = s+p*a[i,1]*sumrec(submat(a, i, 1), x)
		p = p*x
	}
	return(s)
}