I think that in order to achieve complete clarity in finding all real roots of the equation  a=1-(1+x)*exp(-x)  for any real  a  is useful to make the simplest study of the function defined of the right-hand side of the equation:

f:=x->1-(1+x)*exp(-x):

g:=diff(f(x), x);

map(solve, [g=0, g<0, g>0]);

minimize(f(x), x=-infinity..0, location),  maximize(f(x), x=-infinity..0, location);

minimize(f(x), x=0..infinity, location),  maximize(f(x), x=0..infinity, location);

plot([f(x), 1], x=-1.5..5, thickness=[2,1], linestyle=[1,2], scaling=constrained);

We can see (and the plot confirms) that the function  f  is strictly decreasing in the domain  RealRange(-infinity,0) and takes all the values from  RealRange(0,infinity)  and  f  is strictly increasing in the domain  RealRange(0,infinity) and takes all the values from  RealRange(0, Open(1)) . 

On the basis of this information written procedure  SS , which finds all the real roots of the equation with a given accuracy (the numbers of digits):

restart;

SS:=proc(a::realcons, n::posint:=10)

local eq, b;

Digits:=n;

b:=evalf(a);

eq:=a=1-(1+x)*exp(-x);

if b<0 then return `No real solutions` elif

b=0. then return x=0 elif

b>0 and b<1 then return x[1]=fsolve(eq, x=-infinity..0), x[2]=fsolve(eq, x=0..infinity) else return x=fsolve(eq, x=-infinity..0) fi;

end proc:

Examples of use:

SS(-2);  SS(0);  SS(1/Pi);  SS(1);  SS(3, 20);

This system can be solved instantly (and we obtain all solutions) by  solve  command:

restart;

eq[1]:= 0.223569*c_1+2.35589*c_2*c_1^2+0.002356*c_1*c_2^2:

eq[2]:= 1.277899*c_1*c_3-2.350023*c_2*c_3^2+7.5856*c_3*c_2^2:

eq[3]:= 3.225989*c_1^2-2.35589*c_3*c_1^2-7.28356*c_3*c_2^3:

solve({eq[1], eq[2], eq[3]}, {c_1, c_2, c_3});

Your example can be easily generalized to any natural number, if you use the binary number system. In it  any number  S such that  1<= S <2^n  can be recorded as the sequence of no more than  n  of  "0"  and  "1" .

The first procedure for a given number  S  returns the list of the number of dollars in each bag that your condition is true. The second procedure returns the corresponding partition for any number.

MinBags:=proc(S)

[seq(2^n, n=0..ilog[2](S))];

end proc:

Parts:=proc(s)

local L;

L:=convert(s, base, 2);

convert(s,symbol)=convert([seq(`if`(L[i]<>0,convert(L[i]*2^(i-1),symbol),NULL), i=1..nops(L))], `+`);

end proc:

Examples.

Your example:

MinBags(15);

for n from 1 to 15 do

Parts(n);

od;

Another example:

MinBags(137);

for n from 100 to 110 do

Parts(n);

od;

You can easily solve the problem by converting the output into the list, and reversing it:

A:=<1,1,1,1; 1,2,3,4; 4,3,2,1>;

LinearAlgebra[NullSpace](A);

ListTools[Reverse](convert(%, list));

L:={1,2,3,4,5,6,7,8}:

P:=combinat[choose](L, 3):

k:=0:

for p in P do

a:=p[1]; b:=p[2]; c:=p[3];

if a+b>c and a+c>b and b+c>a then k:=k+1 fi;

od:

k;

                        22

Addition:  faster code for large sets.

Example for  L:={$ 1..500}:

restart;

ts:=time():

k:=0:

for a from 1 to 498 do

for b from a+1 to 499 do

for c from b+1 to 500 do

if a+b>c and a+c>b and b+c>a then k:=k+1 fi;

od: od: od:

k;

time()-ts;

                                 10323125

                                    23.875

The first code fails in this example.

solves the problem for  specific parameters values.

Example:

Optimization[QPSolve](2*x+5*y+3*x^2+3*x*y+2*y^2-x*z+4*y*z+2*z^2, {x+y+z=1, 2*x-3*y+z=5});

         [-4.18080357142857, [x = 0.410714285714286, y = -0.897321428571428, z = 1.48660714285714]]

ben maths  You wrote about the mission to give a better code. Here is the shorter and faster code: 

Simp38:=proc(f,x0,xn,n)

local  h:=(xn-x0)/n,  x:=i->x0+i*h;

3*h/8*(f(x0)+f(xn)+3*add(f(x(i)),i=1..n-1,3)+3*add(f(x(i)),i=2..n-1,3)+2*add(f(x(i)),i=3..n-1,3));

end proc:

Test:

ts:=time():

simp38(x->x^5-5.15*x^2+8.55*x-4.05045,1,5,3000);

int(x^5-5.15*x^2+8.55*x-4.05045,x=1..5);  

time()-ts;

ts:=time():

Simp38(x->x^5-5.15*x^2+8.55*x-4.05045,1,5,3000);

int(x^5-5.15*x^2+8.55*x-4.05045,x=1..5);

time()-ts; 

plots[inequal](not (x>4 and y>4), x=-10..10, y=-10..10, color=yellow);

In older versions of Maple (for example in Maple 12)  ListTools:-FindMaximalElement  command does not work. Workaround is   ListTools:-Search  or   ListTools:-SearchAll  commands:

L := [1,2,3,7,6,5,4]:

m:=max(L);

pos:=ListTools:-Search(m,L);

                        m := 7

                       pos := 4

limit(n^(1/n), n=infinity);

limit(exp(n)/n^4, n=infinity);

                               1

                           infinity

For visual purposes only.

Example:

restart;

applyrule(-cos(x::anything)=cos(x+pi), sin(t)-cos(t));

subs(pi=Pi, %);

%;

For a correct perception  add the option  view=[-5..0, 0..1]  into Carl's code after  axes=boxed :

For  x>0  this equation is equivalent to the equation  add(ln(x+j), j=0..2015)=0 , so we have:

Digits:=20:

fsolve(add(ln(x+j), j=0..2015), x=0..infinity);

              0.86678004277349545751*10^(-5785)

It is obvious that this is the only positive root.

restart;

sys:=diff(fi1(t),t,t)=(m0-m)/5, diff(fi2(t),t,t)=m/50:

m0:=200:  m:=(fi1(t)-fi2(t))*32.2+(diff(fi1(t),t)-diff(fi2(t),t))*10:

dsol2:=dsolve({sys,fi1(0)=0,fi2(0)=0,D(fi1)(0)=0,D(fi2)(0)=0}):

fi1:=unapply(rhs(dsol2[1]), t): fi2:=unapply(rhs(dsol2[2]), t):

t1:=fsolve(D(fi1)(t)=10);  t2:=fsolve(D(fi2)(t)=10);

plot([piecewise(t>=0 and t<=t1, D(fi1)(t), 10), piecewise(t>=0 and t<=t2

, D(fi2)(t), 10)], t=0..10, color=[red,green], thickness=2); 

These are the very simple problems. See help on Maple commands:

1.  dsolve

2.  plots[fieldplot]   and   plots[fieldplot3d]

3.  plot3d  - 3rd variant in Calling Sequence