@jfitzsimons 

Old post, but since I had the same question, here is a way, with Maple 2021.

First, define a few variables.

restart:
t:=Pi/17:
a1:=cos(3*t)+cos(5*t)+cos(7*t)+cos(11*t):
a2:=cos(t)+cos(9*t)+cos(13*t)+cos(15*t):
b1:=cos(3*t)+cos(5*t):
b2:=cos(7*t)+cos(11*t):
b3:=cos(t)+cos(13*t):
b4:=cos(9*t)+cos(15*t):
c1:=cos(t):
c2:=cos(13*t):

They are designed such that a1=b1+b2, a2=b3+b3 and b3=c1+c2. If we can compute the pairwise products as well as a1+a2, we are done, as we can solve easily the equation x+y=S, x*y=P: x and y are the roots of t^2-S*t+P.

First check the signs, they will be useful to identify x and y in the above trinomial roots.

andmap(is,[a1>0,a2<0,b1>0,b2<0,b3>0,b4<0,c1>0,c2<0]);

true

Therefore, the variables come in positive/negative pairs.

Now, we need a few values. Note that in this case it's easier for Maple to simplify the exponential form.

S:=e->simplify(convert(combine(e,trig),exp)):

S(a1+a2);
1/2

S(a1*a2);
-1

S(b1*b2);
-1/4

S(b3*b4);
-1/4

S(c1*c2+b1/2);
0

Now define a function to solve the equations x+y=S, x*y=P. You may have to check that for S>0, P<0, the negative solution is the second one.

T:=op@unapply([solve(x^2-s*x+p=0,x)],s,p);

We are almost done.

A1,A2:=T(1/2,-1):
B1,B2:=T(A1,-1/4):
B3,B4:=T(A2,-1/4):
C1,C2:=T(B3,-B1/2):

Numerical check, the values should be close to zero:

evalf[50]~([a1-A1,a2-A2,b1-B1,b2-B2,b3-B3,b4-B4,c1-C1,c2-C2]);

And finally the value of cos(pi/17):

simplify(C1);

1/16-1/16*17^(1/2)+1/16*(34-2*17^(1/2))^(1/2)
+1/16*((-2*17^(1/2)+2)*(34-2*17^(1/2))^(1/2)+12*17^(1/2)+16*(34+2*17^(1/2))^(1/2)+68)^(1/2)

optimize

You need a "predicate type", see satisfies:

TypeTools:-AddType('my_sin',specfunc(sin)^And(rational,satisfies(x->x>1)));
type(sin(x)^(1/2),'my_sin');
type(sin(x)^(3/2),'my_sin');

That does not really solve your problem, but:

assume(x::real):
is(ceil(x)=-floor(-x));
                              true

simplify(floor(-x)+ceil(x));
                      floor(-x) + ceil(x)

is(%=0);
                              true

If you are working with very large sets of combinations, it's usually a better idea to iterate through them one at a time, to keep memory usage low.

Maple can do that with combinat:-nextcomb, Iterator:-Combination or combstruct:-iterstructs.

You may also devise your own functions to achieve this. You can find algorithms for instance in Knuth's TAOCP volume 4A.

Here is an example with the three approaches:

with(combinat):
s:=firstcomb(5,2):
do
  print(s);
  s:=nextcomb(s,5)
until s=FAIL:

with(Iterator):
for a in Combination(5,2) do print(a) od:

with(combstruct):
it:=iterstructs('Combination'(5),size=2):
while not finished(it) do print(nextstruct(it)) od:

Try:

applyrule(exp(k::algebraic*a)=p^k,exp(-a));

First, you can do infolevel[IntegrationTools]:=3 to check what is printed by int. You will at least see that Maple needs the sign of n+1, which should not be a surprise.

To get an answer with or without assumptions:

int(x^n*exp(-x),x=0..infinity,method=MeijerG);

GAMMA(n+1)

And of course, int works with an assumption:

assume(n::nonnegint):
int(x^n*exp(-x),x=0..infinity);

factorial(n)

I see at least two problems:

* gamma is not assigned (you only had gamma = 0.1; not gamma := 0.1;). But you can't assign gamma as it's protected, unless you do first local gamma;.

* In the expression for Eq, you seem to be using square brackets as if they meant parentheses. They do not. To modify the order of evaluation, use parentheses only. Not square brackets, not curly brackets.

See How do I solve shooting method in MaplePrimes and Classroom Tips and Techniques: Numeric Solution of a Two-Point BVP in the Application Center.

Try this:

with(plots):
display(arrow([0,0],[4,4.1],
color=red,shape=double_arrow,width=0.03,border=false,head_width=0.15,head_length=0.2),
arrow([4,4.1],[1,7.9],difference=true,
color=blue,shape=double_arrow,width=0.03,border=false,head_width=0.15,head_length=0.2),
arrow([0,0],[1,7.9],
color=green,shape=double_arrow,width=0.03,border=false,head_width=0.15,head_length=0.2));

Note that the output does not match your definition:

a:=proc(n) option remember;
if n<2 then n else n-2+add(a(k),k=1..n-1) fi end proc:
['a(n)'$n=0..20];

 [0, 1, 1, 3, 7, 15, 31, 63, 127, 255, 511, 1023, 2047, 4095, 
   8191, 16383, 32767, 65535, 131071, 262143, 524287]

However, it's easy to get the same output. Note that n-2 is replaced by n-1 in the formula for a(n):

a:=proc(n) option remember;
if n<2 then n else n-1+add(a(k),k=1..n-1) fi end proc:
['a(n)'$n=0..20];

 [0, 1, 2, 5, 11, 23, 47, 95, 191, 383, 767, 1535, 3071, 6143, 
   12287, 24575, 49151, 98303, 196607, 393215, 786431]

However, even now, 1535 does not match your 1536.