[POI 21 Stage 3] Solar Panels

作法：

題目簡單來說就是對於一個矩形的區域，問這個區域裡 x 坐標根 y 坐標的 gcd 的最大值是多少。考慮一個樸素的算法：枚舉一個數 d 當 gcd ，看看這個區域裡面是否有一個點滿足他的 x 坐標和 y 坐標都是 d 的倍數，接下來要改進這個算法。

假設矩形的左下角叫 ( x0 , y0 )，右上角叫 ( x1 , y1 ) ，那麼當左下角和 d 確定之後，最有可能落在詢問矩形內的兩坐標都是 d 的倍數的數也就確定了，他的 x 坐標就是 >= x0 的第一個 d 的倍數， y 坐標同理。以下記 L1 = x1 - x0 ， L2 = y1 - y0 。因為>= x0 的第一個 d 的倍數和 x 的差會等於 d - ( ( x0 - 1 ) % d ) - 1，所以這個數 <= x1 的充要條件就是 d - ( ( x0 - 1 ) % d ) - 1 <= L1　，或是　( x0 - 1 ) % d >= d - L1 - 1 。y 坐標同理有 ( y0 - 1 ) % d >= d - L2 - 1，所以 d 是一個可行的 gcd 若且唯若上面兩條式子成立。

再來要用到一個神奇的取餘數的性質，就是當 a / b = q 且 a % b = r 時，如果 b 增加 1 並且保持他們的商 q 不變的話，那麼餘數就會減少 q 。因為 a = b * q +r = ( b + 1 ) * q + ( r - q )。為了使用這個性質，我們還需要知道 b 最大可以到多大( 還能夠保持和 a 相除的商是 q  )。設 b 的最大值為 x 好了，因為 a / x = q，所以必須要有 q * x <= a < ( q +1 ) * x ，可以推出 x = a / q 。所以 b 可以到達的最大數就是 a / q = a / ( a / b ) 。這東西也在 POI 14 Stage 1 Queries ( HOJ 22 駭客 ) 出現過，可以參考看看這裡。

這樣一來，記 d2 = min( ( x0-1 ) / (( x0-1 )/d) , ( y0-1 ) / ((y0-1)/d) )，就會得到在 d ~ d2 之間上面那兩條式子的左邊的變動會是線性的。具體來說，設 x0-1 = d * q1 + r1 ，那麼當 d 增加為 d +x 的時候，新的商還是一樣是 q1 ，但新的餘數就會變成 r1 - x * q1 ，而我們希望的是 ( x0 - 1 ) % (d+x) >= (d+x) - L1 - 1 ，其中 d +x 越大越好，現在已經知道左式等於 r1 - x * q1 了，所以移項一下得到 x 必須滿足 ( q1 +1 ) x <= r1 + L1 +1 - d ，同理設 y0-1 = d * q2 + r2 ，那麼 x 也必須滿足 ( q2 +1 ) x <= r2 + L2 +1 - d ，再加上先天的條件 x <= d2 - d ，就可以得到在 d ~ d2 之間最佳的答案是多少了，用他去更新答案即可。

另外在實作上要注意的小細節是，在計算 d2 的時候可能有一邊會變成除以 0 ，這時只要特判讓它變成一個很大的數就好了，反正等一下取 min 就會把它弄掉了 ( 不可能兩數都變成除以0 )。

code :

#include<bits/stdc++.h>
#define INF 1100000000
using namespace std;
 
int divi(int x,int y)
{
    return y==0 ? INF : x/y ;
}
 
main()
{
    int T ; scanf("%d",&T) ;
    while(T--)
    {
        int x0,L1,y0,L2 ;
        scanf("%d%d%d%d",&x0,&L1,&y0,&L2) ;
        L1-=x0 ; L2-=y0 ;
        int ans=1 ;
        for(int d=1;d<=x0 || d<=y0;)
        {
            int d1=min( divi(x0-1,(x0-1)/d) , divi(y0-1,(y0-1)/d) ) ;
            int q1=(x0-1)/d , r1=(x0-1)%d ;
            int q2=(y0-1)/d , r2=(y0-1)%d ;
            if(r1+L1+1>=d && r2+L2+1>=d)
            {
                int x=d1-d ;
                x=min(x, (r1+L1+1-d)/(q1+1) ) ;
                x=min(x, (r2+L2+1-d)/(q2+1) ) ;
                ans=max(ans,d+x) ;
            }
            d=d1+1 ;
        }
        printf("%d\n",ans) ;
    }
}