Задачка для разминки мозга

← →
Rouse_ © (2014-12-25 21:12) [0]

Пришлось тут на днях решать для студента (первого курса) задачку.
Дословно выглядит вот так:

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 vfrcbvev 8 байт).

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

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

В итоге для первого приближения накидал такое решение:

// pk04_dop.cpp : Defines the entry point for the console application. // // вообще-то это олимпиадная задача, я хз зачем ее студенту первого курса дают, но раз задали - значит задали #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; }

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

← →
Kerk © (2014-12-25 21:26) [1]

> Вроде ничего сложного, однако-ж результат суммирования по
> модулю явно не уложится ни в один из поддерживаемых типов
> (int64 vfrcbvev 8 байт).

Не уверен.
Максимальное значение Int64 сильно больше, чем 10^12. При этом ни один из остатков не может превысить 10^2 / 2.

Эта задача скорее на знание чего-нибудь хитрого из теории чисел, имхо.

← →
Kerk © (2014-12-25 21:26) [2]

> не может превысить 10^12 / 2.

единицу потерял

← →
Rouse_ © (2014-12-25 21:27) [3]

Вот так будет на дельфи, если у кого затруднения:

program Project1; {$APPTYPE CONSOLE} {$R *.res} uses Windows, System.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); HalfResult := HalfResult or (Int64(SFlag) 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); HalfResult := HalfResult or (Byte(SFlag) 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.

← →
Rouse_ © (2014-12-25 21:28) [4]

> Не уверен.
> Максимальное значение Int64 сильно больше, чем 10^12. При
> этом ни один из остатков не может превысить 10^2 / 2.

Так этож входное значение, а нужно просчитать сумму всех чисел по модулю yfxbyfz jn Т до нуля. Она там за глаза выходит за Int64

← →
Inovet © (2014-12-25 21:35) [5]

А зачем это делать побитно? Пожно же сразу складывать младшие половины и потом старшие с учётом переноса. Перенос заполучить или из регистра состояния, или как-нибудь сравнением исходных чисел сделать.

← →
Kilkennycat © (2014-12-25 21:41) [6]

...N%1...
a % что означает?

← →
Inovet © (2014-12-25 21:42) [7]

Ну и это... в цикле гнать всю сумму как-то тупо. Здесь должно решаться математически, типа получается сумма целых от 1 до N, а это чё там такое N*(N-1)/2 что ли? Как-то просто слишком, а ну так длинное же и не влазит никуда.

← →
Inovet © (2014-12-25 21:44) [8]

> [6] Kilkennycat © (25.12.14 21:41)
> a % что означает?

Это деление по модулю

← →
Rouse_ © (2014-12-25 21:53) [9]

> a % что означает?

аналог MOD

← →
Rouse_ © (2014-12-25 21:55) [10]

> а это чё там такое N*(N-1)/2

Откуда ты это взял? Нет там такого. Там есть диапазон от 0 до 10 в 12-ой степени, в пределах которого мыжет быть N

← →
Kerk © (2014-12-25 22:13) [11]

> Так этож входное значение, а нужно просчитать сумму всех
> чисел по модулю yfxbyfz jn Т до нуля. Она там за глаза выходит
> за Int64

Максимальное значение суммы какое получается? Я в четверг вечером не готов посчитать :)

В целом, надо с другой стороны заходить. Олимпиады так в лоб не работают. Я запустил твою считалку полчаса назад, она до сих пор считает.

Например, можно разбить ряд на две части.

Одна часть будет
N%1 + N%2+ N%3 + ... + N%(N / 2)

Вторая часть
N%(N / 2 + 1) + N%(N / 2 + 1) + ... + N%N

Сумма второй части считается легко. Это банальная арифметическая прогрессия X*(X+1)/2, где X = ((N - 1) div 2).

То есть вторую часть можно вычислить вообще мгновенно. Над первой часть нужно дальше думать, искать способы. Возможно, тем же самым способом индуктивно можно пойти. Но я не готов сейчас глубже копать. Может быть, завтра.

← →
MBo © (2014-12-25 22:14) [12]

>Откуда ты это взял? Нет там такого
Похожее есть. Сумма второй половины слагаемых Sh находится нетрудно без сложения.
M = (N - 1) div 2
Sh = M * (M + 1) / 2

← →
Inovet © (2014-12-25 22:15) [13]

> [10] Rouse_ © (25.12.14 21:55)
> Нет там такого.

Ну значит я не понял вот это
N%1+ N%2+ N%3+...+ N%N

← →
Rouse_ © (2014-12-25 22:19) [14]

> Максимальное значение суммы какое получается? Я в четверг
> вечером не готов посчитать :)
>
> В целом, надо с другой стороны заходить. Олимпиады так в
> лоб не работают. Я запустил твою считалку полчаса назад,
> она до сих пор считает.

Ага, у меня тоже, даже оптимизированный вариант, но уже значение вышло за диапазон 8 байт :)

> MBo © (25.12.14 22:14) [12]
> >Откуда ты это взял? Нет там такого
> Похожее есть. Сумма второй половины слагаемых Sh находится
> нетрудно без сложения.
> M = (N - 1) div 2
> Sh = M * (M + 1) / 2

Ну что-то примерно похожее: 12499999999999999999999999,875

← →
Kerk © (2014-12-25 22:21) [15]

> X*(X+1)/2, где X = ((N - 1) div 2).

И даже при N = 10^12 сумма второй части получается намного меньше Int64 :)

← →
Inovet © (2014-12-25 22:29) [16]

> [14] Rouse_ © (25.12.14 22:19)
> Ну что-то примерно похожее: 12499999999999999999999999,875

Это для какого N?
калькулятор Windows для 10^12
500 000 000 000 500 000 000 000

← →
Inovet © (2014-12-25 22:32) [17]

124 999 999 999 750 000 000 000

← →
Kerk © (2014-12-25 22:33) [18]

> намного меньше Int64 :)

Неправильно посчитал, действительно может вылезти.

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

← →
jack128 © (2014-12-25 23:11) [19]

> Ну в любом случае, просто суммировать - не выход.

почему не выход?

← →
jack128 © (2014-12-25 23:18) [20]

> Насчет первой части тоже нужно что-то математическое придумать.

Насколько я понимаю, начиная с I := N/3+1 до N/2 - 1 остаток будет прогрессией от 1 до N/3 - 1

← →
Kerk © (2014-12-25 23:18) [21]

Долго потому что. Должна быть возможность решить за разумное время. Это же олимпиадная задачка. Такие задачи изначально не рассчитаны на то, чтоб их по часу считать.

← →
Rouse_ © (2014-12-25 23:30) [22]

Оптимизированный вариант на максимальном значении работает около 6 секунд, не думаю что это сильно долго :) финальное число пока не называю, бо интрига :)

← →
Rouse_ © (2014-12-25 23:31) [23]

Кстати жек, тыж алго решения знаешь, проверь у себя, у тебя такое же время выходит?

← →
Sha © (2014-12-26 10:24) [24]

> MBo © (25.12.14 22:14) [12]
> Сумма второй половины слагаемых Sh находится нетрудно без сложения.

а потом и второй половины от первой половины )

← →
MBo © (2014-12-26 10:28) [25]

Сделал решение за асимптотическое время Sqrt(N) (т.е. порядка нескольких миллионов операций для 10^12, считается мгновенно).
А вот арифметику для 128 бит, вероятно, реализовал с ошибками, отлаживать пока некогда.
Для N = 1000 получается 176919, что совпадает со значением, вычисленным в рамках простой арифметики.
для N = 1000000000000 получается результат со старшим int64-словом $257A=9594 ~ 176*10^21. Младший и старший байты $E0 и $25

← →
Rouse_ © (2014-12-26 12:10) [26]

А как делал?
Через рассчет нижнего и верхнего лимитов арифметических прогрессий?

← →
MBo © (2014-12-26 12:28) [27]

Через сумму квадратов с чередующимися знаками.
Поправка по моим числам (упустил в обном месте возможность переполнения, не уверен, что ещё чего-то не осталось):
старшее слово $2578=9592 ~ 176*10^21. Младший и старший байты $1E и $25

← →
MBo © (2014-12-26 13:14) [28]

Еще одна итерация ;)
старшее слово $2598=9624 ~ 177*10^21. Младший и старший байты $1E и $25

Во как - я по другому ускорял... как доделаешь - покажи, можно для начала в приват rouse79@yandex.ru

> Вроде ничего сложного, однако-ж результат суммирования по
> модулю явно не уложится ни в один из поддерживаемых типов
> (int64 vfrcbvev 8 байт).

Составитель читал книгу Вирта "Описание языка паскаль" и сподобился заменить формулу расчёта факториала своей :-)

> Плохиш © (26.12.14 14:19) [30]
> Составитель читал книгу Вирта "Описание языка паскаль" и
> сподобился заменить формулу расчёта факториала своей :-)

Развивай мысль :)

Угу вечером гляну.
Мой вариант данное число вычисляет за
Z := 10000000000000 div 900000;
используя две вершины каждого из пиков.
т.е. цикл крутится порядка полутора секунд, но еще не до конца допилен (это третий вариант решения).
Твой будет четвертым получается :)

Глянул - суматор для 128 бит у тебя верный, один в один можно сказать. А вот остальной подход пока никто из нашего отдела не понял с квадратами, правда у нас корпоратив и сложно сходу в код вникнуть... :)

- Я потерял нить рассуждения?
Собеседник отвечает, глядя в сторону остальных участников,
- По-моему, он потерял не только нить.

> Inovet © (26.12.14 19:35) [35]

Я отвечал MBo :)
А так-то вообще мой подход (уход от цикла через рассчет арифметической прогрессии, порционно по 900 тысяч элеменов) с посказками доступен тут: http://alexander-bagel.blogspot.ru/2014/12/blog-post.html
Почему именно 900 тысяч? Ну округлил. А так база была - максимальная сумма деления по модулю от 10 в 12-ой которая при суммировании влезет в 8 байт.

Функции 128-битной арифметики, которые мне были нужны (не обрабатываются ситуации, которые исключены в данной задаче):
TUInt128 = record class operator Add(a, b: TUInt128): TUInt128; class operator Subtract(a, b: TUInt128): TUInt128; class operator Implicit(Value: UInt64): TUInt128; class function Sqr128(Value: TUInt128): TUInt128; static; case Integer of 0: (Lo, Hi: UInt64); 1: (LoLo, LoHi, HiLo, HiHi: Cardinal); 2: (Bytes: array [0 .. 15] of Byte); end; { TUInt128 } class operator TUInt128.Add(a, b: TUInt128): TUInt128; begin if a.Lo > UInt64($FFFFFFFFFFFFFFFF) - b.Lo then Inc(a.Hi); Result.Lo := a.Lo + b.Lo; Result.Hi := a.Hi + b.Hi; end; class operator TUInt128.Implicit(Value: UInt64): TUInt128; begin Result.Hi := 0; Result.Lo := Value; end; class function TUInt128.Sqr128(Value: TUInt128): TUInt128; var ad: TUInt128; begin Result.Hi := UInt64(Value.LoHi) * Value.LoHi; Result.Lo := UInt64(Value.LoLo) * Value.LoLo; ad.Lo := UInt64(Value.LoHi) * Value.LoLo; ad.Hi := ad.LoHi; ad.Lo := ad.Lo shl 32; Result := Result + ad; Result := Result + ad; end; class operator TUInt128.Subtract(a, b: TUInt128): TUInt128; begin if a.Lo < b.Lo then Dec(a.Hi); Result.Lo := a.Lo - b.Lo; Result.Hi := a.Hi - b.Hi; end;

Вот начало ряда сумм S(N) и модули, из которых эти суммы складываются.
N S mods 1 0 0 2 0 0 0 3 1 0 1 0 4 1 0 0 1 0 5 4 0 1 2 1 0 6 3 0 0 0 2 1 0 7 8 0 1 1 3 2 1 0 8 8 0 0 2 0 3 2 1 0 9 12 0 1 0 1 4 3 2 1 0 10 13 0 0 1 2 0 4 3 2 1 0 11 22 0 1 2 3 1 5 4 3 2 1 0 12 17 0 0 0 0 2 0 5 4 3 2 1 0 13 28 0 1 1 1 3 1 6 5 4 3 2 1 0
Видно, что правая часть представляет собой арифм. последовательность, сумму которой я уже приводил, а вот левая похитрее. Компактного математического выражения не нашел, но эмпирическое исследование привело к следующему:
главный член суммы резко растет на каждом втором (нечетном) числе, рост квадратичный, и это квадрат выражения (N-1) div 2
Q1 = Sqr((N-1) div 2)
Теперь, если вывести разницу Q1 - S(N), то можно увидеть, что она чаще всего скачет на каждом третьем числе, а сдвиг от начала будет уже 3
Q2 = -Sqr((N-3) div 3)
Далее аналогично
Q3 = Sqr((N-6) div 4)
Q4 = -Sqr((N-10) div 5)
и т.д. Числа сдвига - "треугольные", вида k(k-1)/2, где k - знаменатель.

Итого с использованием длинной арифметики получается функция, которая отличается от такой же для короткой арифметики только типом результата и заменой Sqr на TUint128.Sqr128 (если Multiply реализовать, то в обоих случаях можно одинаково X*X)

function SumOfMods128(N: Int64): TUInt128; var Minus, Divider: UInt64; Even: Boolean; begin Result := 0; Minus := 1; Divider := 2; Even := True; while Minus < N do begin if Even then Result := Result + TUint128.Sqr128((N - Minus) div Divider) else Result := Result - TUint128.Sqr128((N - Minus) div Divider); Minus := Minus + Divider; Divider := Divider + 1; Even := not Even; end; end;

> Rouse_ © (26.12.14 14:27) [31]
>
>
> > Плохиш © (26.12.14 14:19) [30]
> > Составитель читал книгу Вирта "Описание языка паскаль"
> и
> > сподобился заменить формулу расчёта факториала своей :
> -)
>
> Развивай мысль :)

Там был приведён код расчёта факториала для любых чисел с выводом результата :-)
Результат хранился в массиве строк или символов, точно не помню. Почти 20 лет прошло.

Если без больших чисел, то как-то так:

function LoSum(m, n: int64): int64; begin; Result:=0; while m>0 do begin Result:=Result+n mod m; m:=m-1; end; end; function HiSum(m, n: int64): int64; var i: integer; k: int64; begin; Result:=(n-m) * (n-m+1) shr 1; i:=1; while true do begin i:=i+1; k:=n div i; if k<m then exit; k:=(k+m) * (k-m+1) shr 1; Result:=Result-k; end; end; procedure TForm1.Button1Click(Sender: TObject); var m, n, sum: int64; begin; n:=1000; //sum:=LoSum(n,n); m:=trunc(sqrt(n+0.0)); sum:=LoSum(m,n)+HiSum(m+1,n); Memo1.Lines.Add(Format("%d %d",[n, sum])); end;

Задачка для разминки мозга Найти похожие ветки