Sàng nguyên tố và một số ứng dụng

Sàng nguyên tố Eratosthenes

Hôm nay chúng ta sẽ tìm hiểu sàng nguyên tố Sieve of Eratosthenes

Giới thiệu thuật toán

Độ phức tạp thuật toán

Độ phức tạp của thuật toán này là O(nloglogn)

Tư tưởng

Đánh dấu tất cả các số là bội của số nguyên tố p từ p^2 ->n

Demo code (python)

def sieve(n):
    ''' sàng nguyên tố [0,n] '''
    
    danh_dau=[True]*(n+1) # giả sự lúc đầu đều có thể là snt
    
    can_n=int(n**0.5)+1 # = floor(sqrt(n))+1
    
    for i in range(2,can_n+1): # i= 2->can_n
        if danh_dau[i]: # i là số nguyên tố
            
            for j in range(i*i,n+1,i): # j=i*i, i*i+i , ...,n
                danh_dau[j]=False ## j khong là số nguyên tố
    
    
    primes=[]
    for i in range(2,n+1): #i= 2->n
        if danh_dau[i]: primes.append(i) #liệt kê lại số nguyên tố vào mảng mới
    return primes

print sieve(100)
"""
#output:
[2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97]
"""

Ứng dụng

Dưới đây là ứng dụng của sàng số nguyên tố, trong bài toán phân tích ra thừa số nguyên tố.

Code dưới đây thay mảng đánh dấu không phải là [True|False] mà là mảng đánh dấu số nguyên tố trước đó

def gen_sieve_table(n):
    
    ''' sàng nguyên tố [0,n] '''
    
    danh_dau=[0]*(n+1) # giả sự lúc đầu đều có thể là snt
    
    can_n=int(n**0.5)+1 # = ceil(sqrt(n))+1
    
    for i in range(2,can_n+1): # i= 2->can_n
        if danh_dau[i]==0: # i là số nguyên tố -> giá trị =0 [không có ước nguyên tố nhỏ hơn #1]
            
            for j in range(i*i,n+1,i): # j=i*i, i*i+i , ...,n
                danh_dau[j]=i ## j khong là số nguyên tố
    
    
    return danh_dau
    

def factor(n,sieve_table):
    if sieve_table[n]==0: ## là số nguyên tố trả lại ước là chính nó
        return [n] 
    else:
        d=sieve_table[n]  ## chứa 1 ước nguyên tố nhỏ nhất là d
        return [d] + factor(n//d,sieve_table)

sieve_table=gen_sieve_table(100000)

print factor(12345,sieve_table)
output:
[5, 3, 823]

Demo code (C programming language)

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#define maxn 1000000

int a[maxn+1];
int primes[maxn],primes_len;
void sieve(){
    int i,j;
    memset(a,sizeof(int)*(maxn+1),0);
    for(i=2;i*i<=maxn;++i){
        if(a[i]) continue; //i la hop so
        for(j=i*i; j<=maxn;j+=i){ // cac so la boi cua i tu i*i ->n
            a[j]=i;  // cai nay luu lai so nguyen to nho nhat ma j chia het cho i
        }
    }
    /// liet ke lai so nguyen to
    
    primes[0]=2;
    primes_len=1;
    for(i=3;i<=maxn;i+=2){
        if(a[i]) continue;
        primes[primes_len]=i;
        primes_len++;
    }
    
}

int ft[50],fn;

void recrusive_factor(int n){
    if(!a[n]){
        ft[fn]=n;
        fn++;
    }else{
        recrusive_factor(n/a[n]);
        ft[fn]=a[n];
        ++fn;
        
    }
}

void factor(int n){
    fn=0;
    recrusive_factor(n);
}

int main() {
    
    sieve();
    
    printf("so luong so nguyen to <= (%d) = %d\n",maxn,primes_len);
    int n=12345;
    printf("phan tich thua so nguyen to %d=",n);
    factor(n);
    int i;
    
    for(i=0;i<fn-1;++i){
        printf("%d*",ft[i]);
    }
    printf("%d\n",ft[fn-1]);
    return 0;
}
output:
so luong so nguyen to <= (1000000) = 78498
phan tich thua so nguyen to 12345=823*3*5
9 Likes

Bài này em viêt công phu vậy anh nghĩ giờ chưa có ai sửa được wiki đâu. Mà độ phức tạp ở bài này trông lạ nhỉ. Gio giải thích lại xem, lâu rồi anh không nghiên cứu thật toán nên nghĩ chậm quá.

Cách chứng minh độ phức tạp thì em cũng không nắm rõ lắm, trên wiki họ chứng minh rồi. Em chỉ nói cách cài đặt bởi rất đơn giản ~ 10 dòng code, mà lại nhiều ứng dụng.

1 Like

Anh thắc mắc cái này. Trông cái này khó hiểu quớ. Ví dụ như nlogn thì đơn giản. Nhưng nloglogn thì tính sao nhỉ, lặp 2 tầng?


Thông tin ngoài lề tí, mọi người đừng reply liên quan đến nội dung ở dưới, tránh loãng topic hay của Gió. Cảm ơn.

Anh phải trang bị lại kiến thức thuật toán mới được. Anh đi làm về embedded linux, các công việc chủ yếu là liên quan đến kiến trúc của Linux, filesystem, device driver, memory management, … Kiến thức càng sâu thì càng có lợi. Nhưng hiện giờ anh chỉ làm ở tầng cơ bản, tức là áp dụng lại những gì các cao thủ trên thế giới họ làm trước rồi. Anh chỉ đọc hiểu và làm theo.

Như trong kernel có nhiều thuật toán sắp xếp, tìm kiếm để giải quyết việc tìm kiếm và gọi đúng cái hàm phục vụ cho đúng loại Hardware. Tuy nhiên các thuật toán này đã được viết sẵn hết rồi nên anh cũng không quan tâm mấy, tập trung hẳn vào việc áp dụng nó. Trong tương lai, để phát triển nghề nghiệp thì biết sâu về giải thuật rất có lợi. Có điều “tương lai” đó hơi xa =))

3 Likes

Cho hỏi, phần sàng nguyên tố: sao nói set cho nó là ngto, giờ lại hợp số?

memset(a,sizeof(int)*(maxn+1),0);
    for(i=2;i*i<=maxn;++i){
        if(a[i]) continue; //i la hop so
        for(j=i*i; j<=maxn;j+=i){ // cac so la boi cua i tu i*i ->n
            a[j]=i;  // cai nay luu lai so nguyen to nho nhat ma j chia het cho i
        }
    }

Lúc đầu a[i]=0 cho tất cả các số, là giả sử nó là snt. về sau phần j chạy để loại bỏ dần. Vì a[i]!=0 tức là có 1 ước ngto 1<a[i]<i nên không xét

vậy code đó chổ nào là bỏ dần thế?

Như vậy phải ko?

  const long long maxn=100+5;
    bool isPrime[maxn];
    void Prime(){
    	for(long long i=2;i<=maxn;i++){
    		isPrime[i]=true;
    	}
    	for(long long i=2;i<=maxn;i++){
    		if(isPrime[i]) 
    			for(long long j=i;j*i<=maxn;j++){
    				isPrime[j*i]=false;
    			}
    	}
    	
    }

i chỉ cần chạy tới căn maxn, tức là i*i<maxn.
Chú ý dk < maxn vì mảng có maxn phần tử

Chào a, e là thành viên mới. Em mới nghiên cứu python còn vài điều chưa hiểu a có thể giúp dùm e nha.danh_dau=[True]*(n+1) em chưa hiểu gì mấy? là list hả anh.

Thiếu segmented sieve :slight_smile: khắc phục được nhược điểm sử dụng quá nhiều mem bằng cách chỉ sử dụng một vùng nhớ nhỏ hơn kích thước của mem (OS) page và tốc độ cũng cao hơn.

Đầu tiên ta có tổng 1 + 1/2 + 1/3 + 1/5 + 1/7 + … + 1/lastprime(n) = 1 + O(loglog(n)) https://arxiv.org/abs/math/0504289 . Với phiên bản đầy đủ dễ thấy tổng số lần đánh dấu là n(1 + 1/2 + 1/3 + 1/5 + …), hay đpt là n*O(loglogn) với tất cả các thao tác đều là O(1).


Phần trừ ra có sai số rất lớn nên rất vô nghĩa.

======================================

Với segmented sieve thì các tham số của bài là d và n. Để sàng mỗi đoạn [m-d+1…m] ta dùng ngay các kết quả từ trước, nhỏ hơn sqrt(m+d-1). Mỗi đoạn mất O(d*loglog(m)). Cộng tất cả các đoạn thì đpt là O(d*loglog(d*(2d)*(3d)*...)) = O(d*loglog(d^(n/d) * (n/d)!) = O(d*log(n/d*log(d) + n/d*log(n/d)) = O(n*loglogn)

3 Likes
def SieveOfEratosthenes(n):
    prime=list(range(2,n+1))
    i=2
    while i*i<=n:
        for j in range(i*i,n+1,i):
            if j in prime: 
                prime.remove(j)
        i+=1
    return prime

Code này không hay chính vì thao tác remove là O(n).

3 Likes
83% thành viên diễn đàn không hỏi bài tập, còn bạn thì sao?