upto(n, {q = 3}) = {vl = vector(n, i, List()); v = [1, 1, 1]; qterms = q; 
	while(v!=0, 
		s = sum(i = 1, 3, v[i]^2); 
		if(s <= n, 
			if(!!vl[s], 
				listput(vl[s], v); 	
				if(#vl[s] > qterms,
					vl[s] = 0; 
				)
			); 
			v = nxtVl(v); 
		,
			v = nxtVs(v); 
		)
	); 
	printStringforms(vl);
	res = select(x -> #x == qterms, vl, 1); 
	res
}


nxtVs(v) = {
forstep(i = #v - 1, 1, -1, 
	if(v[i] + 2 <= v[#v], 
		v[i]++; for(j = i + 1, #v, v[j] = v[i]); return(v)
	) 
); return([0]); 
}

nxtVl(v) = {
	v[#v]++; v
}

printStringforms(vl) = { 
	for(i = 1, #vl, 
		if(#vl[i] == qterms, 
			print(stringform(i, vl[i]))
		)
	)
}

stringform(n, l) = {
	res = concat(["", n]); 
	for(i = 1, #l,
		res = concat([res, " = ", l[i][1], "^2"]); 
		for(j = 2, qterms, 
			res = concat([res, " + ", l[i][j], "^2"])
		)
	); res 
	
}