第六章 數學問題 -------- 6.11【同餘方程組】POJ1006 生理週期

同餘方程組：

　　先來看一道題目：有物不知其數，三三數之剩二；五五數之剩三；七七數之剩二。問物幾何?  然後我們可以做如下變換，設x為所求的數。

　　　x%3=2              x ≡ a1(%m1)  ①
　　　x%5=3  ===>  x ≡ a2(%m2)  ②
　　　x%7=2              x ≡ a3(%m3)

　　根據前面兩式可以得到  

　　　　x = a1+m1y1     (1) 
　　　x = a2+m2y2

　　兩式相減得到  m1y1 - m2y2 = a2 - a1  這是一個線性不定方程,可解出y1  ---> linearEquation(m1,-m2,a2-a1)   帶回(1),得特解x0 = a1+m1*y1 --> 得到通解表示式 x =x0 + k*lcm(m1,m2) 得一個新方程 x = x0 (mod lcm(m1,m2))

　　程式碼：

/**
 * 
 * @param a 餘陣列成的陣列
 * @param m 模組成的陣列
 * @return
 * @throws Exception
 */
public static long linearEquationGroup(Long[] a, Long[] m) throws Exception {
  int len = a.length;
  if (len == 0 && a[0] == 0)
    return m[0];

  for (int i = 1; i < len; i++) {
    // 這裡往前看是兩個方程
    long a2_a1 = a[i] - a[i - 1];
    long d = linearEquation(m[i - 1], -m[i], a2_a1);
    // 現在的x是y1,用y1求得一個特解
    long x0 = a[i - 1] + m[i - 1] * x;
    long lcm = m[i - 1] * m[i] / d;
    a[i] = (x0 % lcm + lcm) % lcm;// x0變成正數
    m[i] = lcm;
  }
  // 合併完之後,只有一個方程 : x = a[len-1] (% m[len-1])
  //long d = linearEquation(1, m[len-1], a[len-1]);
  return a[len - 1] % m[len - 1];
}

 題目：POJ1006 生理週期

　　程式碼：

import java.util.ArrayList;
import java.util.List;
import java.util.Scanner;

public class POJ1006 {

    public static void main(String[] args) throws Exception {
        Scanner scanner = new Scanner(System.in);
        int t = 1;
        List<Long[]> aList = new ArrayList<Long[]>();
        List<Long> dList = new ArrayList<Long>();
        while(scanner.hasNext()){
            Long []a = {scanner.nextLong(),scanner.nextLong(),scanner.nextLong()};
            Long d = scanner.nextLong();
            if (a[0]==-1&&a[1]==-1&&a[2]==-1&&d==-1) {
                break;
            }else {
                aList.add(a);
                dList.add(d);
            }
        }
        for (int i = 0; i < aList.size(); i++) {
            Long[] a = aList.get(i);
            long d = dList.get(i);
            Long[] m = { (long) 23, (long) 28, (long) 33 };
            long res = MyGcd.linearEquationGroup(a, m);
            while (res <= d) {
                res += 21252;
            }
            System.out.println("Case " + (t++) + ": the next triple peak occurs in " + (res - d) + " days.");
        }
    }
    
    private static class MyGcd {
        static long x;
        static long y;

        /**
         * 
         * @param a 餘陣列成的陣列
         * @param m 模組成的陣列
         * @return
         * @throws Exception
         */
        public static long linearEquationGroup(Long[] a, Long[] m) throws Exception {
            int len = a.length;
            if (len == 0 && a[0] == 0)
                return m[0];

            for (int i = 1; i < len; i++) {
                // 這裡往前看是兩個方程
                long a2_a1 = a[i] - a[i - 1];
                long d = linearEquation(m[i - 1], -m[i], a2_a1);
                // 現在的x是y1,用y1求得一個特解
                long x0 = a[i - 1] + m[i - 1] * x;
                long lcm = m[i - 1] * m[i] / d;
                a[i] = (x0 % lcm + lcm) % lcm;// x0變成正數
                m[i] = lcm;
            }
            // 合併完之後,只有一個方程 : x = a[len-1] (% m[len-1])
            //long d = linearEquation(1, m[len-1], a[len-1]);
            return a[len - 1] % m[len - 1];
        }

        public static long inverseElement(long a, long mo) throws Exception {

            long d = linearEquation(a, mo, 1);
            x = (x % mo + mo) % mo;
            return d;
        }

        public static long linearEquation(long a, long b, long m) throws Exception {
            long d = ext_gcd(a, b);
            // m不是gcd(a,b)的倍數,這個方程無解
            if (m % d != 0)
                throw new Exception("無解");
            long n = m / d;// 約一下,考慮m是d的倍數
            x *= n;
            y *= n;
            return d;
        }

        public static long ext_gcd(long a, long b) {

            if (b == 0) {
                x = 1;
                y = 0;
                return a;
            }
            long res = ext_gcd(b, a % b);
            long x1 = x;// 備份x
            x = y;// 更新x
            y = x1 - a / b * y;// 更新y
            return res;
        }

    }

}

　　結果：