Gauss:=proc(n) local w,xi,i,j,firsteq,eqns,sols,m; global res;
  m:=trunc(n/2); w:=array(0..m); xi:=array(1..m);
  firsteq:=2*sum(w[i],i=1..m);
  if irem(n,2)=1 then firsteq:=firsteq+w[0]  fi;
  eqns:={2=firsteq};
  for j from 2 by 2 to 2*(n-1) do
    eqns:=eqns union{int(x^j,x=-1..1)=2*sum(w[i]*xi[i]^j,i=1..m)};
  od;
  if irem(n,2)=1 then sols:={w[0]}
    else sols:={}
  fi;
  for j from 1 to m do
    sols:=sols union {w[j],xi[j]};
  od;
  res:=solve(eqns,sols);
  evalf(res);
end: