POJ 2115- C Looooops(扩展欧几里德)

C Looooops

Time Limit:1000MS     Memory Limit:65536KB     64bit IO Format:%I64d & %I64u

Submit  Status  Practice  POJ 2115

Appoint description:  System Crawler  (2015-05-26)

Description

A Compiler Mystery: We are given a C-language style for loop of type 

for (variable = A; variable != B; variable += C)

 statement;

I.e., a loop which starts by setting variable to value A and while variable is not equal to B, repeats statement followed by increasing the variable by C. We want to know how many times does the statement get executed for particular values of A, B and C, assuming that all arithmetics is calculated in a k-bit unsigned integer type (with values 0 <= x < 2  ^k) modulo 2  ^k. 

Input

The input consists of several instances. Each instance is described by a single line with four integers A, B, C, k separated by a single space. The integer k (1 <= k <= 32) is the number of bits of the control variable of the loop and A, B, C (0 <= A, B, C < 2  ^k) are the parameters of the loop. 

The input is finished by a line containing four zeros. 

Output

The output consists of several lines corresponding to the instances on the input. The i-th line contains either the number of executions of the statement in the i-th instance (a single integer number) or the word FOREVER if the loop does not terminate. 

Sample Input

3 3 2 16
3 7 2 16
7 3 2 16
3 4 2 16
0 0 0 0

Sample Output

0
2
32766
FOREVER

题意：对于这个循环for (variable = A; variable != B; variable += C)。给出A,B,C和k(k表示变量是在k位机下的无符号整数),判断循环次数,不能终止输出"FOREVER".

思路：我们可以根据所给的循环和容易的得出公式a+c*x=b(mod 2^k),然而可以转换成 c*x=(b-a)mod(2^k)。令A=c，B=b-a，n=2^k，所以原式变成Ax=B (mod n)。

设线性模方程的一个解为x0 
条件①：有d = gcd(a, n) 
条件②：有d = ax1 + ny, 由扩展欧几里得(Egcd)得到x1的值 
条件③：b % d == 0 (有解的条件) 
则x0 = x1*(b/d);

证明:

因为：容易求得d = gcd (a, n), 则有d = ax1 + ny①(扩展欧几里得定理) 
方程①2边同时模n得：d % n == ax1 % n② 
又因为：b % d == 0, 即b是d的倍数; 
所以(b/d)必为整数; 
所以由②得： b % n == d*(b/d) % n == ax1*(b/d) % n == ax % n 
所以很容易可以看出x = x1*(b/d)是方程的一个整数解，得证

需要注意的是: 把模线性方程求得的特解转化为正数之后,要模 b/gcd(a,b) ,而不是b

#include <stdio.h>
#include <math.h>
#include <string.h>
#include <stdlib.h>
#include <iostream>
#include <sstream>
#include <algorithm>
#include <set>
#include <queue>
#include <stack>
#include <map>
using namespace std;
typedef long long LL;
const int inf=0x3f3f3f3f;
const double eps=1e-10;
const double pi= acos(-1.0);
LL exgcd(LL a,LL b,LL &x,LL &y)
{
    if(b==0){
        x=1;
        y=0;
        return a;
    }
    LL r=exgcd(b,a%b,x,y);
    LL t=x;
    x=y;
    y=t-(a/b)*y;
    return r;
}

int main()
{
    LL A,B,C,k,x,y;
    while(~scanf("%lld %lld %lld %lld",&A,&B,&C,&k)){
        if(!A&&!B&&!C&&!k) break;
        LL a=C;
        LL b=(B-A);
        LL n=1ll<<k;
        LL d=exgcd(a,n,x,y);
        if(b%d)
           puts("FOREVER");
        else{
            x=x*(b/d)%n+n;
            printf("%lld\n",x%(n/d));
            //对于无数个解形成的一群余数：周期个数是d，周期长度是n/d，也就是最小正整数解在n/d里
        }
    }
    return 0;
}