# INTCOMP: A Maple package for calculating the complexity of an integer
# (see http://mathworld.wolfram.com/IntegerComplexity.html for definitions)

with(numtheory):

small_divisors := proc(n)
  local AD, SD, d, ub;
  AD := divisors(n) minus {1};
  SD := {};
  ub := floor(sqrt(n));
  for d in AD do
    if d <= ub then
      SD := SD union {d};
    fi;
  od;
  return(SD);
end:

comp := proc(n)
  option remember;
  local SD, i, e;
  if n = 1 then
    return(1);
  fi;
  SD := [op(small_divisors(n))];
  return(min(
    seq(comp(SD[i]) + comp(n/SD[i]), i=1..nops(SD)),
#    comp(n-1) + 1
    seq(comp(e) + comp(n - e), e=1..floor(n/2))
  ));
end:

plot_seq := proc(s)
  local i, xunit, yunit, n, max_entry, height, width;
  height := 1.5;
  width := 5;
  n := nops(s);
  # We round max_entry up to the nearest multiple of 5.
  max_entry := max(op(s));
  # TODO: adjust the two below to allow for labels.
  xunit := evalf((width - 0.25) / (n+1));
  yunit := evalf((height - 0.25) / (max_entry + 1));

  printf("Copy and paste the following into the LaTeX file,\n");
  printf("together with \"\\usepackage{pst-plot}\"\n");
  printf("\n\n\n\n\n");
  printf("\\savedata{\\plotdata}[{");
  for i from 1 to n-1 do
    printf("{%a,%a},", i, s[i]);
  od;
  printf("{%a,%a}}]\n", n, s[n]);
  printf("\\psset{xunit=%gin, yunit=%gin}\n", xunit, yunit);
  printf("\\begin{pspicture}(0,0)(%a,%a)\n", width/xunit, height/yunit);
  printf("\\psaxes[dy=5,Dy=5,dx=100,Dx=100,showorigin=false](0,0)(%a,%a)\n", n-1, max_entry);
  printf("\\rput[c](0,%g){$c(m)$}\n", max_entry + 0.125/yunit);
  printf("\\rput[l](%g,0){$m$}\n", n + 0.125/xunit);
  printf("\\dataplot[plotstyle=dots,dotstyle=*,dotsize=3\\psxunit]{\\plotdata}\n");
  printf("\\end{pspicture}\n");
  printf("\n\n\n\n");
end: