program Wrg1Tab;                  {Produce a TABle of the number of
                                   representations of n as a sum of s
                                   k-th powers}
{$N+,E+}
uses crt, printer;
{$I GetInput.i }

const
MaxN = 99999999999.0;             {10^11 - 1.  If k ò 3 then we store
                                   n^(1/k) k-th powers in an array.  Each
                                   entry occupies 8 bytes.  Thus we need
                                   8n^(1/3) < 64K}
MaxK = 10;
MaxS = 75;
PrinterLimit = 299;               {This allows printing of r(n) for
                                   0 ó n ó 2999.}

var
k,                                {power}
s,                                {number of summands}
l                                 {length of string}
  : integer;
n0,                               {table begins at 10n0}
ll,                               {lower limit of table}
ul,                               {upper limit of table}
t,                                {temporary variable}
ts                                {table size selected for printing}
  : comp;
r                                 {number of representations}
  : array[0..2999] of comp;
pwr                               {k-th powers, computed in advance}
  : array[0..4642] of comp;       {Powers are pre-computed only for k ò 3.
                                   Since ul < 10^11, we have x[i] ó 4641
                                   for all i.}
OriginalMode                      {original video mode -- will be restored}
  : word;
Ch                                {character read from keyboard}
  : char;
St,                               {prompt}
Stk                               {k expressed as a string}
  : string[80];
PrinterOn,                        {Is the printer responding?}
FuncKey,                          {Was the key pressed a function key?}
InputOk                           {Is the input ok?}
  : Boolean;

procedure MakeTable(ll, ul : comp);

var
i, j,                             {indices}
MaxPwr,                           {largest number whose k-th power is needed}
i0,                               {marked index,
                                   used in lexicographic ordering}
mini0,                            {least value of i0 encountered so far}
xpos                              {number of positive coordinates}
  : longint;
n,                                {number represented}
m,                                {multiplicity of the representation}
y                                 {x[s]^k}
  : comp;
x                                 {variables---s of them}
  : array[0..MaxS] of longint;
t                                 {temporary variable}
  : extended;

begin
  MaxPwr := 1 + round(exp(ln(ul)/k));
  if k > 2
  then
    begin
      pwr[0] := 0;
      for i := 1 to MaxPwr do pwr[i] := exp(k*ln(i))
    end;
  for i := 0 to 2999 do r[i] := 0;
  for i := 0 to s do x[i] := 0;
  if ll > 0
  then
    begin
      t := exp(ln(ll)/k);
      x[s] := round(t);
      if pwr[x[s]] < ll
      then x[s] := x[s]+1
    end;
  i0 := s;
  xpos := 1;                      {for n = 0 this is incorrect;          }
  mini0 := s;                     {for such n a correction is made later }
  repeat
    n := 0;
    if k = 2                      {form the sum of s k-th powers}
    then for i := 1 to s do n := n + sqr(x[i])
    else for i := 1 to s do n := n + pwr[x[i]];
    if (n >= ll) and (n <= ul)
    then
      begin
        m := 1;                   {compute the multiplicity with which}
        j := 1;                   {this representation should be counted}
        for i := 1 to s do
          begin
            if (i > 1) and (x[i] = x[i-1])
            then j := j + 1
            else j := 1;
            m := m*i/j
          end;
        if k = 2                  {if k = 2 we allow negative coordinates}
        then m := m*exp(xpos*ln(2));
        i := round(n-ll);         {index in array}
        r[i] := r[i] + m;
        x[s] := x[s] + 1
      end
    else
      begin
        i0 := s-1;
        while (x[i0] = x[i0+1]) and (i0 >= 1) do i0 := i0 - 1;
                   {find i0 so that x[i0] < x[i0+1] = ... = x[s]}
        if i0 < mini0
        then
          begin
            mini0 := i0;
            xpos := s + 1 - i0
          end;
        if (i0 <= s - k) and (KeyPressed)    {allow termination}
        then
          begin
            TextMode(OriginalMode);
            Halt
          end;
        if i0 > 0
        then
          begin
            x[i0] := x[i0]+1;
            for i := i0 + 1 to s do x[i]  := x[i0]
          end;
        n := 0;
        if k = 2
        then
          begin
            for i := 1 to s - 1 do n := n + sqr(x[i]);
            y := sqr(x[s])
          end
        else
          begin
            for i := 1 to s - 1 do n := n + pwr[x[i]];
            y := pwr[x[s]]
          end;
        if n + y < ll
        then
          begin
            t := exp(ln(ll-n)/k);
            x[s] := round(t);
            if pwr[x[s]] < ll - n
            then x[s] := x[s]+1
          end
      end
  until i0 = 0;
  if ll = 0                       {fix error in counting multiplicity of  }
  then r[0] := 1                  {representation of 0                    }
end;                              {of procedure MakeTable}

procedure Display;

var
i,                                {index of rows in table}
j                                 {index of columns in table}
  : integer;

begin
  TextColor(15);
  TextBackground(0);
  ClrScr;
  if k = 2
  then
    begin
      Write('                     NUMBER OF REPRESENTATIONS AS A SUM OF ');
      WriteLn(s:1, ' SQUARES')
    end;
  if k = 3
  then
    begin
      Write('                       NUMBER OF REPRESENTATIONS AS A SUM OF ');
      WriteLn(s:1, ' CUBES')
    end;
  if k > 3
  then
    begin
      Write('                     NUMBER OF REPRESENTATIONS AS A SUM OF ');
      WriteLn(s:1, ' ',k:1, '-TH POWERS')
    end;
  WriteLn;
  Write('       n     10n   10n+1  10n+2  10n+3  10n+4  10n+5  10n+6  ');
  Write('10n+7  10n+8  10n+9');
  for i := 0 to 19 do
    begin
      Write((n0+i):10:0);
      TextColor(14);
      TextBackground(1);
      for j := 0 to 9 do
        begin
          t := r[round(10*(n0+i)+j-ll)];
          if t > 999999
          then
            begin
              str(t:7:0, St);
              l := length(St);
              St := ' ..' +  copy(St, l - 3, 4);
              Write(St)
            end
          else Write(r[round(10*(n0+i)+j-ll)]:7:0)
        end;
      TextColor(15);
      TextBackground(0)
    end;
  WriteLn;
  GoToXY(1, 25);
  TextBackground(7);
  ClrEoL;
  if n0 > 0
  then TextColor(4)
  else TextColor(15);
  Write('     PgUp');
  if n0 < 9999999980.0
  then TextColor(4)
  else TextColor(15);
  Write('        PgDn');
  TextColor(4);
  Write('        s        k        N');
  if PrinterOn
  then
    begin
      TextColor(4);
      Write('        P');
      TextColor(0);
      Write('rint')
    end
  else
    begin
      TextColor(15);
      Write('        Print')
    end;
  TextColor(4);
  Write('        Esc');
  GoToXY(1, 25);                  {hide the cursor}
  TextColor(7);
  Write(' ');
  GoToXY(1, 25);
  TextColor(0)
end;                              {of procedure Display}

procedure PrinTab;                {print the table}

var
i, j                              {indices}
  : integer;

begin
  WriteLn(lst);
  if k = 2
  then
    begin
      Write(lst, '                        WARING''S PROBLEM');
      WriteLn(lst, ' FOR SUMS OF SQUARES')
    end;
  if k = 3
  then
    begin
      Write(lst, '                            WARING''S PROBLEM');
      WriteLn(lst, ' FOR SUMS OF CUBES')
    end;
  if k > 3
  then
    begin
      Write(lst, '                        WARING''S PROBLEM FOR SUMS OF ');
      WriteLn(lst, k:1, '-TH POWERS')
    end;
  WriteLn(lst);
  Write(lst, '     k     10k   10k+1  10k+2  10k+3  10k+4  10k+5  10k+6  ');
  WriteLn(lst, '10k+7  10k+8  10k+9');
  for i := 0 to round(ts) do
    begin
      Write(lst, i:6);
      for j := 0 to 9 do Write(lst, r[10*i+j]:7:0);
      WriteLn(lst)
    end;
  WriteLn(lst)
end;                              {of procedure PrinTab}

begin                             {main body}
  OriginalMode := LastMode;
  TextBackground(1);
  TextColor(14);
  ClrScr;
  GoToXY(6, 2);
  Write('                      WARING''S PROBLEM');
  GoToXY(6,4);
  WriteLn('Will construct a table of the number r(n) of representations');
  GoToXY(6, 5);
  WriteLn('of n as a sum of s k-th powers.  For k > 2 we take r(n) to be');
  GoToXY(6, 6);
  WriteLn('the number of ordered s-tuples (x_1, ..., x_s) of non-negative');
  GoToXY(6, 7);
  WriteLn('integers such that');
  GoToXY(6, 10);
  Write('                  x_1');
  GoToXY(WhereX, 9);
  Write('k');
  GoToXY(WhereX, 10);
  Write(' + ùùù + x_s');
  GoToXY(WhereX, 9);
  Write('k');
  GoToXY(WhereX, 10);
  WriteLn(' = n.');
  WriteLn;
  GoToXY(6, 12);
  WriteLn('In the case k = 2 we allow the x_i to be arbirary integers,');
  GoToXY(6, 13);
  WriteLn('possibly negative (because with this convention the r(n) have');
  GoToXY(6, 14);
  WriteLn('additional properties that otherwise would be missed.)');
  WriteLn;
  WriteLn('     WARNING: The time required to calculate these tables is');
  WriteLn('     roughly proportional to x^((s-1)/k), which is quite large');
  WriteLn('     if s is larger than k.');
  WriteLn;
  str(MaxS:1, St);
  St := '                               (1 ó s ó ' + St + ')';
  s := round(GetInput(6, WhereY, 'Choose the number of summands, s = ',
             St, 1, MaxS));
  WriteLn;
  str(MaxK, St);
  St := '                     (2 ó k ó ' + St + ')';
  k := round(GetInput(6, WhereY, 'Select the exponent, k = ',
                      St, 2, MaxK));
  WriteLn;
  Write('          Computing valuesÄÄÄpress any key to terminate . . . ');
  ll := 0;
  ul := 199;
  MakeTable(ll, ul);
  Ch := 'F';
  FuncKey := False;
  n0 := 0;
  Repeat
    {$I-}
    Write(lst,#0);
    if IOResult = 0
    then PrinterOn := True
    else PrinterOn := False;
    {$I+}
    if (not FuncKey) and (UpCase(Ch) = 'S')
    then
      begin
        str(MaxS:1, St);
        St := '                     (1 ó s ó ' + St + ')';
        s := round(GetInput(1, 24, 'Enter number of summands, s = ',
                            St, 1, MaxS));
        n0 := 0;
        ll := 0;
        ul := 199;
        ClrEoL;
        Write('          Computing valuesÄÄÄpress');
        Write(' any key to terminate . . . ');
        MakeTable(ll, ul)
      end;
    if (not FuncKey) and (UpCase(Ch) = 'K')
    then
      begin
        str(Maxk:1, St);
        St := '          (2 ó k ó ' + St + ')';
        k := round(GetInput(1, 24, 'Enter exponent k = ', St, 2, Maxk));
        n0 := 0;
        ll := 0;
        ul := 199;
        ClrEoL;
        Write('          Computing valuesÄÄÄpress');
        Write(' any key to terminate . . . ');
        MakeTable(ll, ul)
      end;
    if (not FuncKey) and (UpCase(Ch) = 'N')
    then
      begin
        n0 := GetInput(1, 24, 'Start table at 10N, where  N = ',
        '                      (0 ó N ó 9999999980)', 0, 9999999980.0);
        if (10*n0 < ll) or (10*n0 + 199 > ul)
        then                      {need to compute new table of values}
          begin
            if n0 <= 186
            then
              begin
                ll := 5*n0;
                ul := 20*n0+200
              end
            else
              begin
                ll := 10*n0 - 1500;
                ul := ll + 2999
              end;
            if ul > 99999999999.0
            then
              begin
                ll := 99999997000.0;
                ul := ll + 2999
              end;
            ClrEoL;
            Write('          Computing valuesÄÄÄpress');
            Write(' any key to terminate . . . ');
            MakeTable(ll, ul)
          end
      end;
    if (not FuncKey) and (UpCase(Ch) = 'P')
    then
      begin
        ll := 0;
        ts := GetInput(1, 24, 'Will print r(0) ... r(10N+9) where N = ',
              '                             (19 ó N ó 299)', 19, 299);
        ul := 9 + 10*ts;
        GoToXY(1, 24);
        ClrEoL;
        GoToXY(1,25);
        ClrEoL;
        Write('          Computing valuesÄÄÄpress');
        Write(' any key to terminate . . . ');
        MakeTable(ll, ul);
        GoToXY(1, 25);
        ClrEoL;
        Write('         Printing table of values for 0 ó n ó ', ul:1:0);
        Write('. . . ');
        PrinTab
      end;
    if (FuncKey) and (Ch = #73)     {PgUp}
    then
      begin
        n0 := n0 - 20;
        if n0 < 0 then n0 := 0;
        if 10*n0 < ll
        then                      {need to compute new table of values}
          begin
            ll := 5*n0;
            ul := 20*n0+200;
            if ul - ll > 2999
            then
              begin
                ll := 10*n0 - 1500;
                ul := ll + 2999
              end;
            if ul > 99999999999.0
            then
              begin
                ll := 99999997000.0;
                ul := ll + 2999
              end;
            ClrEoL;
            Write('          Computing valuesÄÄÄpress');
            Write(' any key to terminate . . . ');
            MakeTable(ll, ul)
          end
      end;
    if (FuncKey) and (Ch = #81)     {PgDn}
    then
      begin
        n0 := n0 + 20;
        if n0 > 9999999980.0 then n0 := 9999999980.0;
        if 10*n0 + 199 > ul
        then                      {need to compute new table of values}
          begin
            ll := 5*n0;
            ul := 20*n0+200;
            if ul - ll > 2999
            then
              begin
                ll := 10*n0 - 1500;
                ul := ll + 2999
              end;
            if ul > 99999999999.0
            then
              begin
                ll := 99999997000.0;
                ul := ll + 2999
              end;
            ClrEoL;
            Write('          Computing valuesÄÄÄpress');
            Write(' any key to terminate . . . ');
            MakeTable(ll, ul)
          end
      end;
    Display;
    repeat
      Ch := ReadKey;
      if Ch <> #0
      then FuncKey := False
      else
        begin
          FuncKey := True;
          Ch := ReadKey
        end;
      InputOk := False;
      if FuncKey and (Ch = #73) and (n0 > 0)
      then InputOk := True;
      if FuncKey and (Ch = #81) and (n0 < 9999999980.0)
      then InputOk := True;
      if (not FuncKey) and (Ch IN [#27, 'k', 'K', 'n', 'N', 's', 'S'])
      then InputOk := true;
      if PrinterOn and (not FuncKey) and (UpCase(Ch) = 'P')
      then InputOk := True
    until InputOk
  until (not FuncKey) and (Ch = #27);       {Esc was pressed}
  TextMode(OriginalMode)
end.