Post date: Jul 5, 2011 11:48:01 AM
Tuyến tính (linear): Tiếng Hán-Việt, "tuyến" có nghĩa là "thẳng". Một phương trình tuyến tính (hay còn gọi là phương trình bậc một hay phương trình bậc nhất) là một phương trình đại số có dạng:
f(x) = ax + b
Đồng dư (congruential): "Cho số nguyên dương n, hai số nguyên a,b được gọi là đồng dư theo mô-đun n nếu chúng cho cùng số dư khi chia cho n (hay là a-b chia hết cho n)." Khi nhắc đến đồng dư, đồng nghĩa với nhắc đến phép mod
100 mod 3 = 1
4 mod 3 = 1
Ta nói, 100 và 4 là đồng dư theo module 3.
Hàm đồng dư tuyến tính (Linear congruential equation): có dạng f(x) = (ax + b) mod m
f(x) = (2x + 1) mod 10
Sinh số ngẫu nhiên theo phương pháp đồng dư tuyến tính được D.Lehner đưa ra vào 1951, với công thức
a[0] = seed
a[i] = (a[i-1]*b + 1) mod m;
Trong đó, seed, b, m là các hằng số. Thông thường khi dùng các thư viện random trong java, hay C, ta truyền vào tham số seed.
Nhược điểm:
Giả sử cần tính a*b không tràn:
a = a1 << k + a0 ( ví dụ 28 = (1<<4) + 12)
b = b1 << k + b0
a * b = (a1 * b1) << (2 * k) + (a1 * b0 + a0 * b1) << k + a0 * b0
Giả sử a là một số sign 32bit, để loại bỏ thành phần dấu, ta ngầm định hiểu 1 số unsigned 30bits. chọn k = 15. Thành phần (a1*b1)<<(2*k) gây tràn số -> bỏ qua
a*b = (a1*b0 + a0*b1)<<15 + a0*b0
Để tránh trường hợp (a1*b0 + a0*b1)<<15 gây tràn số, ta dùng thủ thuật sau
a * b = ((a1 * b0 + a0 * b1) & 0xEFFF)) << 15 + a0*b0
Với
a1 = (a>>15)
a0 = a & 0xEFFF
b1 = b>>15
b0 = b&0xEFFF
Gọi hàm mult dùng cho việc nhân 2 số không tràn, được viết như sau:
function mult(a, b)
{
a1 = a >> 15;
a0 = a& 0xEFFF b1 = b >> 15;
b0 = b & 0xEFFF return (((a1 * b0 + b0 * b1) & 0xEFFF) << 15) + a0 * b0;
}
Lúc này, a[i] được viết lại thành
a[i] = (mult(a[i-1],b) + 1) mod m
Để tính random trong 1 khoảng 0-> (r-1)
a[i] = (mult(a[i-1], b) + 1) mod m
randomint = (a*r)/m
Sinh số ngẫu nhiên theo phương pháp đồng dư cộng (Additive congruential)
Phương pháp này phát triển từ phương pháp đồng dư tuyến tính dựa trên ý tưởng "thanh ghi dịch chuyển hồi tiếp tuyến tính - linear-feedback shift register". Với đồng dư cộng, sẽ có n số random được tính trước theo phương pháp đồng dư tuyến tính, các số tiếp theo được tổ hợp (cộng) từ các số đã tính.
a[0] = 0;
a[i] = a[i] = (a[i-1]*b + 1) mod m; với i < NUM
a[j] = (a[j - b] + a[j - c]) mod m; với i >= NUM
class CongruentialRandom {
public final static int m_iMod_Shift = 30;
public final static long m_iMod_Shift_Mask = 0x3FFFFFFF L;
public final static int m_iHaftMod_Shift = 15;
public final static long m_iHaftMod_Shift_Mask = 0x7FFF L;
public CongruentialRandom() {}
public long Mult(long a, long b) {
long a1 = (long) a >> m_iHaftMod_Shift;
long a0 = a & m_iHaftMod_Shift_Mask;
long b1 = (long) b >> m_iHaftMod_Shift;
long b0 = (long) b & m_iHaftMod_Shift_Mask;
long re = ((((a1 * b0 + a0 * b1) & m_iHaftMod_Shift_Mask) << m_iHaftMod_Shift) + a0 * b0) & m_iMod_Shift_Mask;
return re;
}
public int NextInt(int max) {
return 0;
}
}
class LinearCongruentialRandom extends CongruentialRandom {
public final static long k_iLinear = ((long) 1 << 28) + 21;
long m_iVal;
public LinearCongruentialRandom(long seed) {
super();
m_iVal = seed;
}
public int NextInt(int max) {
m_iVal = (Mult(m_iVal, k_iLinear) + 1) & m_iMod_Shift_Mask;
long re = (m_iVal * max) >> m_iMod_Shift;
return (int) re;
}
}
class AdditiveCongruentialRandom extends CongruentialRandom {
public final static long k_iLinear = 31;
public final static int k_iSizeOfCache = 55;
public final static int k_iMaxOfCache = k_iSizeOfCache - 1;
public final static int k_iMidOfCache = (k_iSizeOfCache >> 1) + 1;
long m_iVal[] = new long[k_iSizeOfCache];
int m_iGenCounter = 0;
public AdditiveCongruentialRandom(long seed) {
super();
m_iVal[0] = seed;
for (int j = 1; j < k_iSizeOfCache; j++) {
m_iVal[j] = (Mult(k_iLinear, m_iVal[j - 1]) + 1) & m_iMod_Shift_Mask;
}
m_iGenCounter = k_iSizeOfCache - 1;
}
public int NextInt(int max) {
m_iGenCounter = (m_iGenCounter + 1) % k_iSizeOfCache;
m_iVal[m_iGenCounter] = (m_iVal[(m_iGenCounter + k_iMidOfCache) % k_iSizeOfCache] + m_iVal[(m_iGenCounter + k_iMaxOfCache) % k_iSizeOfCache]) & m_iMod_Shift_Mask;
long re = (m_iVal[m_iGenCounter] * max) >> m_iMod_Shift;
return (int) re;
}
}