30/09/2018, 16:27

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

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

ĐPT 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

###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 phân tích ra thừa số ngyên tố:

####nếu thay mảng dá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]

C code

#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
",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
",ft[fn-1]);
    return 0;
}
output:
so luong so nguyen to <= (1000000) = 78498
phan tich thua so nguyen to 12345=823*3*5
Nguyễn Minh Dũng viết 18:30 ngày 30/09/2018

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á.

Gió viết 18:39 ngày 30/09/2018

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.

Nguyễn Minh Dũng viết 18:43 ngày 30/09/2018

O(nloglogn)

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 =))

Đạt Đỗ viết 18:31 ngày 30/09/2018

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
        }
    }
Gió viết 18:37 ngày 30/09/2018

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

Đạt Đỗ viết 18:31 ngày 30/09/2018

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

Đạt Đỗ viết 18:38 ngày 30/09/2018

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;
    			}
    	}
    	
    }
Gió viết 18:30 ngày 30/09/2018

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ử

Đoàn Vũ Luân viết 18:38 ngày 30/09/2018

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.

rogp10 viết 18:36 ngày 30/09/2018

Thiếu segmented sieve 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 chấp nhận, không chứng minh tổng 1 + 1/2 + 1/3 + 1/5 + 1/7 + … + 1/lastprime(n) = 1 + O(loglog(n)) (định lí Mertens). 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).

math.stackexchange.com
Eric Naslund

How does $ \sum_{p<x} p^{-s} $ grow asymptotically for $ \text{Re}(s) < 1 $?

number-theory, asymptotics, analytic-number-theory
answered by Eric Naslund on 04:36PM - 04 Jul 11

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(dloglog(m)). Cộng tất cả các đoạn thì đpt là O(d*loglog(d*(2d)*(3d)*...)) = O(dloglog(d^(n/d) * (n/d)!) = O(dlog(n/dlog(d) + n/dlog(n/d)) = O(nloglogn).

Bài liên quan
0