Comments 56
Если обобщать до k пропавших, то может возникнуть проблема с большими порядками для вычислений. Есть альтернатива: ipsit.bu.edu/documents/ieee-it3-web.ps
Наверняка можно сочинить алгоритм для сведения этой системы к полиному одной переменной (фактически, ищем точку пересечения гиперповерхностей первого тире катого порядка, где к — количество изъятий), а затем найти для этого полинома корни методом Лобачевского.
Поэтому решение с битовыми индексами выглядит приятнее, а главное, легче обобщается.
Если подумать, то это вполне элегантное решение.
Интересно, не будет ли тогда проще внаглую решить систему из варианта решения со степенями обычным градиентным методом?
Относительно применимости градиентного метода к этой почти-Вандермондовой матрице — надо оценки выписывать и эксперименты проводить, иначе так среди допущений и останемся. Корень кстати у такой системы один (потому как она около-Вандермондова), значит проблем со сходимостью скорее всего не будет.
Меня очень смущает то, что система из задачи со степенями очень похожа на систему из задачи интерполяции (является ее минором по сути), как-то одна задача должна во вторую выливаться с одинаковым итогом — полиномом
По хорошему, её бы надо перевести в базисные симметрические полиномы. Если для этого есть явные формулы — хорошо. Иначе, при больших k на их вывод потребуется слишком много времени и памяти.
А если значения симметрических полиномов будут известны — останется щёлкнуть пальцами и найти корни одного многочлена k-й степени от одной переменной :) При хорошем раскладе, при вычислении коэффициентов не придётся даже выходить из кольца целых чисел.
Базисные симметрические полиномы задаются полиномиальными коэффициентами из 0 и 1, то есть нужно поискать возможность из наборов полиномиальных коэффициентов, соответствующих строкам системы «собрать» эти самые базисные. Интуиция мне подсказывает, что это можно сделать, на первых порах — пусть даже лобовым перебором, а потом поискать закономерности и т. п.
— берём случайное число a из Zp
— считаем многочлен Q=((x-a)(p-1)/2-1) mod R
— находим C=gcd(Q,R).
С вероятностью 1-21-deg( R ) многочлен C будет отличен от 1 и от R. Так что мы разложили R на два множителя: R=C*(R/C). Продолжаем процесс, пока не получим линейные множители.
class Test {
static uint P=4294967291;
static uint Add(uint a,uint b) {
return (uint)(((ulong)a+b)%P);
}
static uint Sub(uint a,uint b) {
return (uint)(((ulong)a+P-b)%P);
}
static uint Mul(uint a,uint b) {
return (uint)(((ulong)a*b)%P);
}
static uint Inv(uint a) {
uint b=P,c=1,d=0;
bool neg=false;
while(a!=1) {
uint x=b%a,y=b/a;
b=a; a=x;
x=d+c*y; d=c; c=x;
neg=!neg;
}
return neg ? Sub(0,c) : c;
}
static int K;
static uint[] PowSum; // [K+1]
static uint[] Polynom;
static void AddPowers(uint n,bool add) {
uint r=1;
for(int i=0;i<=K;i++) {
PowSum[i]=add ? Add(PowSum[i],r) : Sub(PowSum[i],r);
r=Mul(r,n);
}
}
static void NewtonConversion() {
Polynom=new uint[K+1];
Polynom[K]=1;
for(int s=1;s<=K;s++) {
uint res=0;
for(int a=1;a<=s;a++) {
if(a%2==0) res=Sub(res,Mul(Polynom[K-s+a],PowSum[a]));
else res=Add(res,Mul(Polynom[K-s+a],PowSum[a]));
}
Polynom[K-s]=Mul(res,Inv((uint)s));
}
}
static uint[] MulMod(uint[] pol1,uint[] pol2,uint[] mod) {
int D=pol1.Length;
uint[]res=new uint[D];
for(int i=D-1;;i--) {
uint a=pol2[i];
for(int j=0;j<D;j++) res[j]=Add(res[j],Mul(a,pol1[j]));
if(i==0) break;
a=res[D-1];
for(int j=D-1;j>0;j--) res[j]=Sub(res[j-1],Mul(mod[j],a));
res[0]=Sub(0,Mul(mod[0],a));
}
return res;
}
static void CalcPow(uint[] pol,uint pow,uint[] mod) {
uint maxpow=1u<<31;
while(maxpow>=pow) maxpow>>=1;
uint[] res=(uint[])pol.Clone();
while(maxpow>0){
res=MulMod(res,res,mod);
if((pow&maxpow)!=0) res=MulMod(res,pol,mod);
maxpow>>=1;
}
Array.Copy(res,pol,res.Length);
}
static int Deg(uint[] pol) {
for(int i=pol.Length-1;i>=0;i--) if(pol[i]!=0) return i;
return -1;
}
static uint[] CalcGCD(uint[] pol1,uint[] pol2) {
for(;;) {
int d1=Deg(pol1),d2=Deg(pol2);
if(d2<0) return pol1;
while(d1>=d2) {
uint a=Mul(pol1[d1],Inv(pol2[d2]));
if(Mul(a,pol2[d2])!=pol1[d1]) {
break;
}
pol1[d1]=0;
for(int i=0;i<d2;i++) pol1[i+d1-d2]=Sub(pol1[i+d1-d2],Mul(a,pol2[i]));
d1--;
}
uint[] p=pol1; pol1=pol2; pol2=p;
}
}
static void DivPol(int idx0,int pow,uint[] pol,int deg) {
uint b=Inv(pol[deg]);
for(int i=0;i<deg;i++) {
pol[i]=Mul(pol[i],b);
Polynom[idx0+pow-deg+i]=Sub(Polynom[idx0+pow-deg+i],pol[i]);
}
pol[deg]=1;
for(int j=pow-1;j>=deg;j--) {
b=Polynom[idx0+j];
int dr=idx0+j-deg;
for(int i=0;i<deg;i++) Polynom[dr+i]=Sub(Polynom[dr+i],Mul(b,pol[i]));
}
}
static void Factor(int idx0,int pow,uint shift,List<uint> ans) {
while(pow>1) {
uint[] pol=new uint[pow];
uint[] pol2=new uint[pow+1];
pol[1]=1; pol[0]=shift;
for(int i=0;i<pow;i++) pol2[i]=Polynom[i+idx0];
pol2[pow]=1;
CalcPow(pol,(P-1)/2,pol2);
if(Deg(pol)<0){ shift++; continue; }
pol=CalcGCD(pol,pol2);
int deg=Deg(pol);
if(deg==0) { shift++; continue; }
DivPol(idx0,pow,pol,deg);
for(int i=0;i<deg;i++) Polynom[i+idx0]=pol[i];
shift++;
Factor(idx0+deg,pow-deg,shift,ans);
pow=deg;
}
ans.Add(Polynom[idx0]);
}
static List<uint> Missing(IEnumerable<uint> seq,int k) {
K=k;
PowSum=new uint[K+1];
uint n=1;
foreach(uint a in seq) {
AddPowers(n++,true);
AddPowers(a,false);
}
for(int i=0;i<K;i++) AddPowers(n++,true);
NewtonConversion();
List<uint> res=new List<uint>();
Factor(0,K,0,res);
return res;
}
static void Main(string[] args) {
uint[] test3_1= { 1,5,7,2,6,8,9 };
List<uint> res=Missing(test3_1,3);
foreach(uint a in res) Console.Write("{0} ",a);
Console.WriteLine();
}
}
Со второй вершиной пути не очень понятно. Если стартовая вершина лежит на нескольких циклах, то вторую вершину можно определить по разному. Если она одна, сходу ничего в голову не лезет. Подумаю.
Пусть у нас есть алгоритм, который использует t бит памяти.
Мы строим граф на основе циклической перестановки (такой, что p(x) != x для всех x). Граф двудольный, состоит из долей A и B. В каждой доле n вершин. Сначала проводим ребра из вершины a_i в вершину b_{p_i} для каждого i. Скормим эти ребра алгоритму и посмотрим на его память.
Я утверждаю, что из памяти можно извлечь всю перестановку. Если это так, то для разных перестановок состояния памяти должны быть разные. Значит, должно быть, что 2^t >= (n — 1)! или t = Omega(n log n).
Почему можно извлечь перестановку? Зафиксируем некоторое значение j и проведем ребра из b_i в a_i для вcех i != j. Скормим их алгоритму. Мы получили Эйлеров граф, в котором вторая вершина будет как раз b_{p_j}. Алгоритм её должен выдать. Понятно, что такую процедуру можно запустить для всякого j и восстановить перестановку.
P.S. Примерно такой же конструкцией можно показать, что и вероятностные алгоритмы тоже не помогут.
Кстати, вот девятая задача действительно какая-то головоломная. Я знаю, как показывается нижняя оценка на подсчет числа инверсий (тоже Omega(n)). Но вот именно четность не получается пробить.
возьмём два момента: после анализа n/3 чисел и после анализа 2*n/3. Предположим, что в каждый из этих моментов у нас есть только k бит информации. И докажем, что мы проиграли :)
Вопрос — можно ли, сохраняя более одного бита на промежуточных этапах, тем не менее, уменьшить объём нужной информации?
Для n=3 существует конечный автомат с тремя состояниями (т.е. 1.6 бита), который может вычислить чётность всех 6 перестановок.
Ребра, которые сначала алгоритм обработает чтобы получить интересную память, будут иметь вид a_i -> b_{p_i}. Чтобы узнать, куда ползет элемент j, нужно добавить в граф ребра b_i -> a_i для всех i != j.
А с точками можно обойтись и без полной сортировки пар. Достаточно перебрать все точки, для каждой из них отсортировать направления на остальные точки, и потом проверить, есть ли одинаковые. Время то же — N^2*log(N), а памяти нужно гораздо меньше. Да и код получится несколько проще.
Но такой подход требует N+1 бита (можно обойтись N, элемент с индексом 0 нам не нужен). Не думаю, что это минимум. Я бы предположил, что хватит O(sqrt(N)), но это просто гадание ни на чём.
Задача 8 может решаться в один проход "подсчетом", отличающимся тем что считать нужно только однобитовый флаг присутствия (значит рабочая память будет M/8 байт) и тем, что при первой же необходимости прирастить счет значения которое уже равно 1 мы заканчиваем)
Если из условия "Задача 8. В последовательности записано M+1 целое неотрицательное число" предполагать что сама последовательность влезает в память компьютера. (Хотя это может быть не так)
(Пишу потому что) Хотя иначе мы приходим все же к вашей же сложности, но с более практическим пониманием идеи диапазонов. Мы выделяем столько счетчиков "парности" сколько можем по памяти, а после считаем присутствие для значений, N старших бит которого указывают соответствующий групповой счетчик. (Эта классика и вызывает у меня стойкую ассоциацию с сортировкой подсчетом)
Оценка по памяти: 1Гб сможет адресовать 8_388_608 диапазонов по 23 старшим битам. Последующими проходами мы просто игнорируем все не интересующее за пределами выбранного диапазона (идентифицированного по старшим битам) и точно так же индексируем значения внутри него по следующим битам, повторяя пока биты не закончатся. Т.о. например за три прохода даже с 512Мб памяти мы уже сможем решать это для 64битных значений M.
Счетчиками парности я их назвал наводя на идею о том, что если завести не однобитовые а больше флаги, то можно искать таким образом и наибольшие/наименьшие значения встречающиеся более/менее чем N раз. Например с 4Гб памяти за 3 прохода можно найти наименьшее/большее 64битное число встречающееся более 2 или 3 раз. (что на больших разреженных числах выглядит довольно полезным инструментом выявления каких-то корреляций. Практическое приложение - анализ коллизий хэш функций/ или случайных генераторов на домашнем компе)
log=async(...arg)=>new Promise(r=>{ //async log
console.log(...arg); setTimeout(r);
})
async function someCollisionOfUi32(iterable){
const map=new Uint32Array(65536); //256K
// поиск группы с кол-вом членов >65535 (минимум 1 повторился):
let i=0, collision= Math.floor(Math.random()*0x100000000);
console.log('predefined collision was:', collision);
for(x of iterable(collision, 123)){
if(!(++i&0x3FFFFFF))await log((i/0x100000000).toFixed(2));
if(++map[x>>>=16] > 65536){
map.fill(0);
// поиск числа из группы повторившегося более раза
for(y of iterable(collision, 123)){
if(!(++i&0x3FFFFFF))await log((i/0x100000000).toFixed(2));
if((y>>>16)==x&&map[y&65535]++)return y; //гарантированно случиться
}
break;
};
}
return -1; //дубли отсутствуют
}
function*TestSeq(mix, seed=Math.floor(Math.random()*0x100000000)){
let x=new Uint32Array(1);
x[0]=seed;
for(let n=0;n<0x100000000;n++){
x[0]=x[0]*69069 +1; // Ui32 modulo by overflow
yield x[0];
}
yield mix; yield mix;
return;
}
// console.log([...TestSeq(4)])
console.time("collise");
console.log('someCollisionOfUi32:',await someCollisionOfUi32(TestSeq));
console.timeEnd("collise"); //~219s
Попробовал решение 8й задачи на всех ui32 числах. С тестом как водиться провозился дольше чем с алгоритмом. В итоге:
- тест на JS и работает из консоли браузера
- лог асинхронный (что позволяет вкладке не виснуть намертво до конца выполнения) выводит процент готовности прохода в консоль (всего их два)
- потребление памяти 256Кб
- скорость обхода 2**k+2 чисел ~219s (3 с половиной минуты на ryzen 9)
- тестируется LCG PRNG последовательность всех ui32 чисел + в конец ее добавляется дважды случайно взятое число (оно и будет давать коллизию)
Дополнения:
- в принципе 2**32 чисел это всего то 2**29 = 536 870 912 байт индекса, что можно побитово замапить флагами четности и решить за один проход вдвое быстрее, но потратив 512Мбайт вместо 256Кбайт. Можно оценить разницу по памяти. Впрочем работая с ui64 надо именно так поступать, потому что битность кратно уменьшает кол-во проходов. Я потестил ui32 версию потому что хотелось строгий тест провести с наличием всех чисел. А для ui64 их в 2*32 раз больше... В принципе адаптировать из этого теперь элементарно. Но на полной последовательности просто забейте это дома считать...
- на самом деле для хранения переполнения 16битным счетчикам не хватало бы 1го значения. Пробовал версию которая бьет диапазоны по 65535 чисел - тогда памяти требуется даже 128K + 8б, но из-за делений индексов и взятия остатков по модулю не степени 2ки - где-то вдвое падает производительность. Можно отдельные флаги переполнения ввести для групп потратив +8К и по скорости будет почти также. Но будет более запутанный код.
- забавно сразу же прикладной фактор показал себя... пока тестил пришлось взять линейный конгруэнтный генератор, взял первый попавшийся из вики, и бах алг что то выбрасывает раньше, гляжу опять в вики в генераторе написано, что не все биты в нем используются (т.е. нет полной конгруэнтности, значит есть дубли) - и алг уже это показывал сразу. Причем даже с случайным сидом в основном отбрасывал на нескольких конкретных числах. Если бы нужно было проанализировать этот генератор, это была бы ценная информация для оценки куда копать
За один проход