UOJ Logo mcfxmcfx的博客

博客

快速乘法取模的一种奇怪实现

2021-08-20 17:02:51 By mcfxmcfx

常见的取模优化有 Montgomery multiplicationBarrett reduction 等。当模数固定时,编译器会帮我们进行 Barrett reduction 的优化,而当模数不固定时,这两种都需要手动实现。而手动实现的速度总会比编译器慢(比如模数啥的必须占一个寄存器位置)。

那么能不能在程序中内置编译器呢?即使能,也多半超过了码长限制。但是我们可以白嫖编译器的成果,让编译器对着某个模数优化,我们之后再改掉代码里的模数。

实现如下:

const int P=1000000007;

#ifdef _WIN32
#include<windows.h>
bool set_rwx(void*addr)
{
    char tmp[8];
    return VirtualProtect(addr, 0x1000, PAGE_EXECUTE_READWRITE, (PDWORD)tmp);
}
#else
#include<sys/mman.h>
bool set_rwx(void*addr)
{
    return !mprotect(addr, 0x1000, PROT_READ | PROT_WRITE | PROT_EXEC);
}
#endif

namespace mod_jit
{
    struct pointer{
        union{uint64_t*ul;uint32_t*ui;uint16_t*us;uint64_t v;};
        pointer(){}
        template<typename T>pointer(T x){ul=(uint64_t*)x;}
    };

    pointer range_l,range_r;

    bool add_segment(uint64_t addr,uint64_t elfsz,uint64_t memsz)
    {
        if(elfsz+0x10000<memsz)return 1; // skip global arrays
        uint64_t tmp=range_l.v+addr+elfsz;
        if(range_r.v<tmp)range_r.v=tmp;
        return 0;
    }

#ifdef _WIN32
    void get_range()
    {
        pointer t=get_range;
        t.v&=~0xfffull;
        while(*t.ui!=0x905a4d)t.v-=0x1000;
        range_l=t;
        t.v+=*(t.ui+15);
        int num=*(t.us+3);
        t.v+=0x18+*(t.us+10);
        for(;num--;t.v+=0x28){
            if(*(t.ui+9)>>25&1)continue;
            if(add_segment(*(t.ui+3),*(t.ui+4),*(t.ui+2)))break;
        }
    }
#else
    void get_range()
    {
        pointer t=get_range;
        t.v&=~0xfffull;
        while(*t.ui!=0x464c457f)t.v-=0x1000;
        int num=*(t.us+28);
        range_l=t;
        t.v+=*(t.ul+4);
        for(;num--;t.v+=0x38){
            if(add_segment(*(t.ul+2),*(t.ul+4),*(t.ul+5)))break;
        }
    }
#endif

    void get_change(int P,std::function<void(uint32_t)>add_int,std::function<void(uint64_t)>add_ll)
    {
        assert(P>(1<<29)&&P<(1<<30));
        auto tint=[&](uint32_t x){add_int(x);add_int(-x);};
        auto tll=[&](uint64_t x){add_ll(x);add_ll(-x);};
        for(int i=-3;i<=3;i++)tint(P+i);
        tint(((((1ull<<30)-P)<<32)+P-1)/P);
        for(int i=59;i<=61;i++)tint(((1ull<<i)+P-1)/P);
        for(int i=91;i<=93;i++)tll(((__uint128_t(1)<<i)+P-1)/P);
    }

    int vc,tc,pl;
    std::vector<std::pair<int,pointer>>p;

    void record_int(uint32_t x)
    {
        pointer a=range_l,b=range_r;
        b.v-=3;
        while(a.v!=b.v)
        {
            if(*a.ui==x)p.push_back(std::make_pair(vc,a));
            a.v++;
        }
        vc++;
    }

    void record_ll(uint64_t x)
    {
        pointer a=range_l,b=range_r;
        b.v-=7;
        while(a.v!=b.v)
        {
            if(*a.ul==x)p.push_back(std::make_pair(vc,a));
            a.v++;
        }
        vc++;
    }

    void change_int(uint32_t x)
    {
        for(;pl<p.size()&&p[pl].first==tc;pl++)
            *p[pl].second.ui=x;
        tc++;
    }

    void change_ll(uint64_t x)
    {
        for(;pl<p.size()&&p[pl].first==tc;pl++)
            *p[pl].second.ul=x;
        tc++;
    }

    uint64_t anti_optimize_zero()
    {
        uint64_t res=0;
        for(int i=0;i<0x100;i++){
            res=res+55555^114514;
            res^=res<<9;
        }
        return res^11608655508041216768ull;
    }

    struct init{init()
    {
        get_range();
        get_change(P+anti_optimize_zero(),record_int,record_ll);
        std::set<uint64_t>map_area;
        for(auto&a:p)map_area.insert(a.second.v&~0xfffull);
        for(uint64_t x:map_area)assert(set_rwx(*(void**)&x));
    }}init_;

    void change_mod(int P)
    {
        tc=pl=0;
        get_change(P,change_int,change_ll);
    }
}

用法:粘上这坨,然后用这个固定模数 P 写代码。要改模数的时候调用 mod_jit::change_mod

速度参考:https://loj.ac/s/1233497 (原始提交:https://loj.ac/s/1233405 https://loj.ac/s/1170461

这坨代码的主要流程是,先找到程序基地址,然后读取每个段,得到整个程序的地址区间,接下来在程序中找到编译器可能优化出来的常数,把对应内存设为可修改的,每次 change_mod 时就修改他们。

然而这份代码也只是个概念验证,想粘个板子直接用,目前还不太现实。

有一个主要问题,这个找常数的过程找到的不一定真的是常数,他有可能是指令的一部分恰好和常数相同(虽然概率非常小)。为了解决这个问题,就需要把程序按指令划分开,虽然能做到(比如参考这个),但是成本会高得多。

另一个问题是,编译器优化出的常数,可能性非常多,很难枚举完。在 get_change 函数里面可以看到许多带 P 的表达式,他们勉强能覆盖大多数情况,但是显然不能覆盖所有情况。另外这里还要求 229<P<230,这是因为编译器根据模数不同,生成的移位指令也不同,而想要写代码找出这些移位指令是比较困难的。(另一种解决方法是,对于每个 2k2k+1 之间的 P 都生成一份代码,这倒也许还不错)

挑战图同构(大雾)

2021-08-02 18:10:28 By mcfxmcfx

前几天在 EI 群看到 EtaoinWu 在 UOJ 搞了一道编码题 #679. 无标号图编码,感觉很有意思,于是就来做了做。

根据 OEIS,不同的图个数约为 2n(n1)/2n!,而这个结果也非常的精确。写一下代码可以发现,在 n=70 时,这个估计和实际值的差距已经不超过实际值的 1016

n=100 时,2n(n1)/2n!1.17624425。根据信息熵,能编码的 01 串长度最多只能是 4425,而此时分数为 170。

引入

首先考虑把点定一个顺序。一种简单的方法是,把 0K1 连成独立集,然后 i+K 往前面 j 连边当且仅当 gi 的第 j 位是 1,其中 gi[0,2K)。 令 fi=j[gj=i],我们需要保证 fi{0,1},并且改变 0K1 的顺序会使得 f 不同(这样才能还原出 0K1 的顺序)。假设 f 是在程序中预先写死的。 解码时可以先找到独立集,然后根据后面点的连边情况还原整个序列。假设 K=7,7 个点形成独立集的概率非常低,而 f 也满足的概率就更低了。

这个做法一共使用了 762+793=672bits 来区分顺序,剩下还有 93922=4278bits 可以存储实际序列。

算法一

7 实际上也太大了。7 元独立集本身就占用了许多边,而很多数只需要 6bits 确定,不需要 7bits。我们可以让 K=5,然后令 fi=j[gj=i],并且移除 fi{0,1} 的限制。可以让 f 中存在一个 fi=1,一个 fi=2,并且所有 fi8fi=1 的点标号可以唯一确定。设他为 x,我们可以用是否向 x 连边区分出 fi=2 的点。接下来所有点都可以用和这三个点的连边情况来区分了。胡乱实现一个可以得到 93 分

算法二

这种排序方法仍然浪费了大量 bit。他的本质是,用 n 个 bit,可以把大小为 n 的集合分为两个集合。当 n 较大时,这种方法的浪费不大(如 2log10+20log202.5bits),但是 n 小的时候会造成非常大的浪费(如 n=2 需要 2bits,然而实际信息量只有 1bit)。

假设我们按照前面的方法确定了几个点的标号(把这几个点叫做 B,最前面 K 个叫做 A),接下来问题的就是,如何更快地把剩下的集合排序。对于 fi 对应的集合,可以强制要求里面每个点到 B 的连边情况都不同,那么它花掉了 |B||fi|bits,但是可以贡献 log(2|B|fi)bits。经过计算可以发现,当 B 较大时,净损失几乎就是 logfi!,也就是几乎没有额外开销。实现上,可以在标程里抄一个 BinomTable

在此基础上,还有一些优化。把每个集合贡献的 log(2|B|fi)bits 乘起来,作为一个高精度整数考虑,可以多赚回来一些 bit。把 K 换成 4,可以节省一些 bit。

这些优化全部应用上可以得到 138 分

算法三

刚才的优化始终没有绕开必须有 fi=1fi=2 作为 B 的种子这个条件。为了绕开,我们可以要求 B 的内部是一个不对称图(asymmetric graph),这样就可以直接在内部确定编号。

假设 |B|=88 个点的不对称图有 3696 个。这样净开销就只有 872log(36968!)=0.85bits,非常节省。应用这个可以得到 143 分

算法四

前面的过程一直假设 f 是固定的,而这也损失了大量的信息。如果把 f 变成动态的,f 本身也可以编码大量信息。我当时的实现可以得到 158 分。(但是这代码有个 bug,一直到 170 分的提交才发现并修复,所以应该还可以再得更高一点的分)

算法五

在算法四的基础上,经过调参可以发现,如果放宽了 f 的范围,那么一个图可能会有多种解码方式。假设在这些方式中选择 hash 值最小的。

经过多次调参,发现可以让算法四将 100 个点的图和长为 4427 的 01 序列建立联系,并且约有 30% 的序列编码再解码可以得到自身(实际上上界是 1.176/4=0.294,也就是说几乎每个图都可以编码了)。

于是我们可以考虑这个新问题:有一个 F:[0,24427){0,1},约 30% 的 F(x) 是 1,并且是随机分布的,需要找一个 G:[0,24425)[0,24427),以及 G1,使得 F(G(x))=1

假设我们建了一棵 4425bits 的 trie 树,每个叶节点有 04y 满足 F(y) 为 1,有一个需要求 G 的值 x。每个节点内部先尽量把 xy 匹配,然后匹配不完的再往父亲丢。

这个 trie 树当然不用实际建出,可以从叶子开始往上一层一层搜索,直到找到匹配位置。

显然最后所有的 x 都能找到匹配,只是时间复杂度可能会爆,感性理解一下可以发现大概有 exp(i) 的概率在 exp(i) 的时间内跑出来。

最后我的实现在题目要求的时限内错误率略高于 2%,多次提交可以得到理论上界 170 分

(这份代码实现上有一些小区别,但是大致思想如此)

碎碎念

看到这题时,我首先想到了之前看到的一道题,大概是把一个长为 100 的串编码到 100 个点的图中,但是他还限制了边数也不超过 100。那题可以用一条链来做,这给了我把点编号的启发。

这题到达了理论上界并不意味着图同构被解决了。假设有两个图 A,B,其中 A 是随机图,B 要么和 A 同构,要么是另一个随机图。那么将度数序列排序后比较可以有几乎为 1 的正确率。这题其实就类似这种情况。

这题的对偶版本——将图编码成 01 串,看起来也已经解决了(并且几乎达到了下界)。但是有个大问题:对于 01 串编码的情况,即使串本身不是随机的,也可以随机异或一些东西,把他变随机;但是给你一个精心构造的图,并不能把他变得随机。不过,这其实也并不是问题,因为要是 corner case 都能做,那就直接可以判定图同构了。

部分交互题的通用 hack 方法

2020-12-05 08:57:59 By mcfxmcfx

背景

对于一部分交互性题目,他们的考点是调用某个函数的次数。对于某些题目,如果询问次数足够大,那么暴力算法(或者比正解简单的多的做法)也可以通过。

假设交互器不是 adaptive 的(即数据是预先确定的),那么我们可以考虑调用询问函数,然后再把询问次数改回去。这其实很容易做到,比如在程序中内嵌一个 qemu 之类的模拟器,把询问函数模拟执行一遍,并记录下他改动过的内存,事后再改回去。

实现

qemu 体量太大,比较难修改,并且修改完可能也压不进代码长度限制。所以需要找一个比较简单的 x86_64 模拟器。

我经过一番搜索,找到了 https://github.com/blbogdan/x64vcpu 这个项目,他的体量就比较小。去除 FPU 和 SSE 等一些用处不大的东西(这里其实埋了个坑,之后会讲到),再去掉 asm(为了满足 uoj 的要求),得到了 https://github.com/mcfx0/x64vcpu 这样一个东西。

接下来就可以去 hack 交互题了,比如 https://uoj.ac/submission/441760

不过对于其他一些题目,比如 https://uoj.ac/problem/165 ,他的交互次数非常多,暴力本身的运行时间就已经接近时限,这种情况下如果还是每次询问都模拟,显然会 TLE。所以可以在最开始模拟一次询问,获取他修改的内存位置,最后再把这些内存改回去,于是询问次数也被改了回去(比如 https://uoj.ac/submission/441773 )。

至于前面说到的坑,其实是 libc 的 memset 用到了 SSE,这导致询问函数中用到 memset,这份代码就失效了。(但是只要把 SSE 的实现加回来,就没有这个问题了)

可能的防范方法

  • 让交互器变得 adaptive(比如 https://uoj.ac/problem/486 )。不过这并不意味着完全防住了本方法。比如也许可以在最开始进行若干次随机询问,使得数据固定,之后再使用本方法。

  • 算询问次数使用一些累加之外的手段,比如第 i 次询问把 a[i] 设为 1。这样做可以防住只在最后把内存改回去(因为中间询问改了哪是完全不知道的),但是不能防范每次都暴力模拟。

  • 在 grader 中使用各种奇怪的指令集。治标不治本。

  • 使用 IO 交互。

  • 让代码段不可读取。其实我不太清楚这样会不会影响程序运行,有机会的话可以试试。不过这样搞的话,评测机那边可能也需要相应做出一些调整。

  • 不在网络赛出这样的交互题。线下赛应该不可能有人背的住这样的板子。

12.07 更新:修了个 bug,改了改 api,可以见 github,也可以见 https://uoj.ac/submission/441827

mcfxmcfx Avatar