четверг, 25 декабря 2014 г.

Предновогодняя задачка для разминки мозга

Пришлось тут на днях решать для студента (первого курса) задачку.
Дословно выглядит вот так:
2014 ACM-ICPC China Hunan Invitational Programming Contest
There is a simple problem. Given a number N.
You are going to calculate N%1+ N%2+ N%3+...+ N%N.
Input: The length N(1<=N<=10^12).
Output: Answer.
Sample input 5
Sample output 4
Вроде ничего сложного, однако-ж результат суммирования по модулю явно не уложится ни в один из поддерживаемых типов (int64 максимум 8 байт).

Сразу уточню: задачка олимпиадная, причем свежего 14-го года :)


Конечно странно что такие вещи задают первокурснику, но...

В итоге для обкатки накидал такое решение на двух полусумматорах - банальная битовая арифметика:

#include "stdafx.h"

// хранилище для оооочень большого числа
unsigned char UInt128[16];

void printUInt128(){
 printf("result = 0x");
 for (int i = 15; i >= 0; i--)
  printf("%hhX", UInt128[i]);
 printf("\n");
}

// полный битовый сумматор реализующий таблицу истинности: http://ivatv.narod.ru/zifrovaja_texnika/1_04.htm
bool fullAdder(bool a, bool b, bool p0, bool *p){
 // первый полусумматор
 bool bRes = false;
 *p = false;
 if ((a && !b) || (b && !a))
  bRes = true;
 if (a && b)
  *p = true;
 if (!p0)
  return bRes;
 // второй полусумматор
 *p = true;
 bRes = !bRes;
 if (!a && !b)
  *p = false;
 return bRes;
}

// сумматор на битовых операциях
void add(long long x){
 bool pFlag = false;
 unsigned long long halfResult = 0; 
 unsigned long long bigValue;
 int i;
 bool aBit, bBit, sFlag;

 // получаем указатель на массив, содержащий большое число
 unsigned long long *p;
 p = (unsigned long long*)&UInt128[0];

 // берем младшие 8 байт
 bigValue = (unsigned long long)(*p);
 for (i = 0; i < 64; i++){
  // и побитно, посредством полного сумматора складываем два числа
  aBit = ((unsigned long long)1 << i & x) > 0;
  bBit = ((unsigned long long)1 << i & bigValue) > 0;
  sFlag = fullAdder(aBit, bBit, pFlag, &pFlag);
  if (sFlag)
   halfResult |= (unsigned long long)1 << i;
 };
 // результат помещаем обратно
 *p = halfResult;

 // если нет переноса бит от предыдущей операции, то можно выходить
 if (!pFlag)
  return;

 halfResult = 0;

 // сдвигаемся на 8 байт
 p++;

 // берем старшие 8 байт
 bigValue = (unsigned long long)(*p); 
 for (i = 0; i < 64; i++){
  // увеличиваем значение опираясь на бит переноса
  bBit = ((unsigned long long)1 << i & bigValue) > 0;
  sFlag = fullAdder(false, bBit, pFlag, &pFlag);
  if (sFlag)
   halfResult |= (unsigned long long)1 << i;
 };
 // результат помещаем обратно
 *p = halfResult;
}


int _tmain(int argc, _TCHAR* argv[])
{
 unsigned long long a, i;

 printf("enter value: ");
 scanf("%lld", &a);

 printf("calculating...\n");
 i = a;
 while(i > 0){
  add(a % i);
  i--;
 };
 printUInt128();

 getchar();
 getchar();

 return 0;
}

Либо, чтобы было проще понять, вариант на дельфи, но (увы) без коментариев:

program Project1;
 
{$APPTYPE CONSOLE}
 
{$R *.res}
 
uses
  Windows,
  SysUtils,
  Math;
 
type
  Int128 = array [0..15] of Byte;
 
var
  VeriBigInteger: Int128;
 
procedure PrintInt128;
var
  I: Integer;
begin
  Write('0x');
  for I := 15 downto 0 do
    Write(IntToHex(VeriBigInteger[I], 2));
  Writeln;
end;
 
function FullAdder(A, B, P0: Boolean; var P: Boolean): Boolean;
begin
  Result := False;
  P := False;
  if A and not B then
    Result := True;
  if B and not A then
    Result := True;
  if A and B then
    P := True;
  if not P0 then Exit;
  P := True;
  Result := not Result;
  if not A and not B then
    P := False;
end;
 
procedure Add(Value: Int64);
var
 I: Integer;
 ABit, BBit, SFlag, PFlag: Boolean;
 HalfResult, BigValue: Int64;
begin
 PFlag := False;
 HalfResult := 0;
 BigValue := PInt64(@VeriBigInteger[0])^;
 for I := 0 to 63 do
 begin
   ABit := (Int64(1) shl I and Value) > 0;
   BBit := (Int64(1) shl I and BigValue) > 0;
   SFlag := FullAdder(ABit, BBit, PFlag, PFlag);
   if SFlag then
      HalfResult := HalfResult or (Int64(1) shl I);
 end;
 PInt64(@VeriBigInteger[0])^ := HalfResult;
 
 HalfResult := 0;
 BigValue := PInt64(@VeriBigInteger[8])^;
 for I := 0 to 63 do
 begin
   BBit := (1 shl I and BigValue) > 0;
   SFlag := FullAdder(False, BBit, PFlag, PFlag);
   if SFlag then
      HalfResult := HalfResult or (Int64(1) shl I);
end;
 
 PInt64(@VeriBigInteger[8])^ := HalfResult;
end;
 
function TstNative(X: int64): int64;
var
  I: Int64;
begin
  I := X;
  Result := 0;
  while I > 0 do
  begin
    Result := Result + (X mod I);
    Dec(I);
  end;
end;
 
procedure TstOld(X: int64);
var
  I: Int64;
begin
  ZeroMemory(@VeriBigInteger[0], SizeOf(Int128));
  I := X;
  while I > 0 do
  begin
    Add(X mod I);
    Dec(I);
  end;
  PrintInt128;
end;
 
var
  N: int64;
begin
  try
  repeat
    Readln(N);
    Writeln('0x', IntToHex(TstNative(N), 32));
    TstOld(N);
  until N = 0;
  except
    on E: Exception do
      Writeln(E.ClassName, ': ', E.Message);
  end;
  readln;
end.

Сразу скажу - ЭТО работает, причем работает правильно (можно даже использовать для проверки), однако есть более грамотный и более быстрый вариант решения в 17 строчек кода, вместо процедур Add и FullAdder (на дельфи или на си - кому как удобнее) :)

Задача - попробовать найти данный вариант решения или предложить еще более оптимальный вариант :)

С наступающим :)

UPD: дам подсказку, график прохода по модулю выглядит вот так:


Определение верхнего лимита прогрессии:

X := Limit div (Z + 1) - (Z - (Limit mod (Z + 1)));
Где:
X - предел положительной вершины
Limit - значение числа N из оригинальной задачи
Z - текущее прирастание прогрессии (т.е. при Z = 1, прогрессия 0, 1, 2... при Z = 5, прогрессия 0, 5, 10... )

т.е. грубо идем от 0 до 1000 - пики выглядят:

499, 332, 247, 196, 165 и т.д....

17 комментариев:

  1. Мы кембриджев не кончали, в верхней математике не сильно шарим - мы по простому

    type
    Int128 = packed record
    case integer of
    0: (i128: array [0..15] of byte);
    1: (s0,s1: int64);
    end;
    var
    R: int128;
    i,N: int64;
    Over: boolean;
    begin
    N:=Round(power(10,12));
    R.s1:=0;
    R.s0:=0;
    i:=1;
    Over:=false;
    while i <= N do begin
    R.s0:=R.s0+(N mod i);
    Inc(i);
    if (R.s0 < 0) and (not Over) then Over:=true;
    else if (R.s0 >= 0) and Over then begin
    Over:=false;
    Inc(R.s1);
    end;
    end;
    end;

    Веселого нового года :)

    ОтветитьУдалить
    Ответы
    1. А я все думаю - зачем люди строки нумеруют?
      Вот и мне подарок к новому году , теперь знаю - чтобы глупая машина структуру кода не съедала.

      1 type
      2 Int128 = packed record
      3 case integer of
      4 0: (i128: array [0..15] of byte);
      5 1: (s0,s1: int64);
      6 end;
      7
      8 var
      9 R: int128;
      10 i,N: int64;
      11 Over: boolean;
      12
      13 begin
      14 i:=$7fffffffffffffff;
      15 N:=Round(power(10,12));
      16 R.s1:=0;
      17 R.s0:=0;
      18 i:=1;
      19 Over:=false;
      20 while i <= N do begin
      21 R.s0:=R.s0+N mod i;
      22 Inc(i);
      23 if (R.s0 < 0) and (not Over) then begin
      24 Over:=true;
      25 end
      26 else if (R.s0 >= 0) and Over then begin
      27 Over:=false;
      28 Inc(R.s1);
      29 end;
      30 end;
      31 end;

      Удалить
  2. Ответы
    1. В комментариях подсветка кода и разметка не поддерживается :)
      По самому коду, в принципе да. как-то так и надо, только ждать долго, около трех дней, поэтому нужно что-то более эффективное :)

      Удалить
  3. Модуль и сложение никуда не денешь.
    Более эффективное - только на ассемблере

    ОтветитьУдалить
    Ответы
    1. Конечно денешь, в этом то собственно и есть суть задачи. Более того скажу, что самый быстрый из известных мне вариантов решения решает данную задачу для N = 10^12 за 73 миллисекунды :)
      1 января решение будет опубликовано.

      Удалить
    2. суть олимпиадных задач и состоит в знаниях и применениях алгоритмов для быстрого решения, а не решения в лоб. Так что естественно решаться задача будет не по данной в условии формуле прямым счетом.

      Удалить
  4. Так где же решение?

    ОтветитьУдалить
  5. You are going to calculate N%1+ N%2+ N%3+...+ N%N - а это разве не арифметическая прогрессия?

    ОтветитьУдалить
    Ответы
    1. Не совсем, грубо говоря это набор арифметических прогрессий.

      Удалить
  6. Скажите, где можно найти вторую часть ответа на задачку номер 2?

    ОтветитьУдалить
    Ответы
    1. А я ее не публиковал, но там в кратце трюк заключается в восстановлении pwndOrg и снятия флага DCX_INUSE у структуры tagDCE при вызове GetDCEx с флагом DCX_INTERSECT

      Удалить
    2. А нужно сделать, чтобы было всем понятно, а то не совсем понятно. :)

      Удалить
    3. Статью писать лень, вот тут я основные позиции вкратце накидал, сишный код взят из исходников w2k.
      http://rouse.drkb.ru/tmp/canvas.zip

      Удалить
    4. Думаю не только мне интересно.
      Да и как-то странно выглядит, что нет 2-й части на ответ. Не хорошо это.
      А заметочку можно за неделю осилить для блога, чтобы он жил.

      Посмотрел код, честно говоря не понял откуда в такой простой задачке появились структуры и т.д. :)

      Удалить
    5. Структуры - это скажем там внутренности GDI. А статьи пока что не будет, дело в том что когда я готовил вторую часть ответа и перепроверял себя наткнулся на один не понятный мне самому момент, без изучения которого ответ был бы не полным. А на его изучение требуется время, которого и так мало, поэтому решил все это отложить в сторонку. Как созрею - так сразу опубликую все целиком :)

      Удалить