Wednesday, October 8, 2008

การเขียนโปรแกรมหาจำนวนเฉพาะโดยใช้ Python

ไหนๆก็ไปโพสท์ไว้ใน blognone แล้ว ก็เลยเอามาแปะไว้นี่ด้วยเลยเผื่อเป็นประโยชน์

ย้ายไปเขียนบน LongSpine.com เนื่องจากโค้ดเละ

วิธีการหาจำนวนเฉพาะ ..แบบเขียนง่ายๆอ่านง่ายๆ


ไอเดียแรกคือการหาจำนวนเฉพาะใหม่จากจำนวนเฉพาะเก่าที่มีอยู่แล้ว โดยมีหลักคือจำนวนเฉพาะนั้นคือจำนวนที่ไม่มีจำนวนเฉพาะอื่นๆหารมันลงตัว

prime = [2]
for i in range(3, 100000):
flag = True # สมมุติว่าจำนวน i เป็นจำนวนเฉพาะ
for j in prime: # จำนวน j ใดๆที่เป็นจำนวนเฉพาะต้องหารจำนวนนี้ไม่ลงตัว
if (i % j == 0): # แต่ถ้าดันหารได้ลงตัว (มีเศษเป็น 0)
flag = False # เราก็พบว่าจำนวนนั้นไม่ใช่จำนวนเฉพาะ
break # ออกจากลูปทันที
if flag == True: # ถ้าพิสูจน์แล้วว่าเป็นจำนวนเฉพาะจริงๆ
prime.append(i) # ให้ใส่ค่านั้นเข้าไปในตัวแปร prime

ไอเดียที่สองคือลดจำนวนตัวที่ต้องวนหารลง เนื่องจากจำนวนเต็มใดๆ j สูงสุดที่จะหาจำนวนเต็มใดๆ i ได้นั้น จะต้องมีค่าไม่มากไปกว่ารากที่สองของ i

prime = [2]
for i in range(3, 100000):
flag = True # สมมุติว่าจำนวน i เป็นจำนวนเฉพาะ
for j in prime: # จำนวน j ใดๆที่เป็นจำนวนเฉพาะต้องหารจำนวนนี้ไม่ลงตัว
if (j ** 2 > i): # ถ้า j มีค่ามากกว่าตัวหารสุงสุดที่เป็นไปได้ (มีค่ามากกว่ารากที่สองของ i)
break # ให้ออกจากลูป
if (i % j == 0): # แต่ถ้าดันหารได้ลงตัว (มีเศษเป็น 0)
flag = False # เราก็พบว่าจำนวนนั้นไม่ใช่จำนวนเฉพาะ
break # ออกจากลูปทันที
if flag == True: # ถ้าพิสูจน์แล้วว่าเป็นจำนวนเฉพาะจริงๆ
prime.append(i) # ให้ใส่ค่านั้นเข้าไปในตัวแปร prime


ไอเดียที่สาม เรายังทำให้เร็วได้มากไปกว่านี้อีก โดยใช้วิธีซีฟของเอราทอสเทนีส ซึ่งใช้หน่วยความจำมาช่วยในการคำนวน (อ่านรายละเอียดได้ตามลิงค์) วิธีนี้ข้อเสียคือเปลืองหน่วยความจำมาก

import numpy

n = 100000 # หาจำนวนเฉพาะที่มีค่าน้อยกว่า 100000

primeArr = numpy.zeros(n) # สร้าง array ขนาด n
primeArr[0] = 1 # 0 ไม่ใช่จำนวนเฉพาะ
primeArr[1] = 1 # 1 ก็ไม่ใช่

for i in xrange(2, n):
if primeArr[i] == 0: # ถ้าจำนวนนั้นเป็นจำนวนเฉพาะ
for j in xrange(i * 2, n, i): # จำนวนอื่นๆที่จำนวนนั้นหารลงตัว
primeArr[j] = 1 # ..จะไม่ใช่จำนวนเฉพา

prime = []

for i in xrange(2, n):
if primeArr[i] == 0: # เก็บจำนวนเฉพาะไว่ในตัวแปร prime
prime.append(i)
ทดลองเปรียบเทียบความเร็วของทั้งสามตัวอย่าง

poomk@gemini:/tmp$ time ./prime1.py >0

real 0m19.622s
user 0m19.581s
sys 0m0.020s

poomk@gemini:/tmp$ time ./prime2.py >0

real 0m0.675s
user 0m0.676s
sys 0m0.000s

poomk@gemini:/tmp$ time ./prime3.py >0

real 0m0.535s
user 0m0.504s
จะเห็นว่าเร็วขึ้นตามลำดับ

ไอเดียที่สี่ เนื่องจากตัวอย่างที่สามนั้นไม่เหมาะกับการหาจำนวนเฉพาะที่มีขนาดมากๆ (มากๆๆๆๆ...ๆ) ทำให้เราไม่สามารถใช้วิธีของซีฟได้เพราะเรามีหน่วยความจำไม่เพียงพอ วิธีหนึ่งที่(ผมเคย)ใช้ก็คือการหาจำนวนเฉพาะตั้งต้นด้วยซีฟ(ตัวอย่างที่สาม) จากนั้นค่อยน้ำตัวเฉพาะตั้งต้นเหล่านั้นไปหาตัวเฉพาะอื่นๆต่อตามตัวอย่างที่สอง

import numpy

n = 1000000 # หาจำนวนเฉพาะที่มีค่าน้อยกว่า 1000000 (หนึ่งล้าน)

primeArr = numpy.zeros(n)
primeArr[0] = 1
primeArr[1] = 1

for i in xrange(2, n):
if primeArr[i] == 0:
for j in xrange(i * 2, n, i):
primeArr[j] = 1

prime = []

for i in xrange(2, n):
if primeArr[i] == 0:
prime.append(i)

for i in range(n, 2000000): # หาต่อจนถึงตัวที่ 2000000 (สองล้าน)
flag = True
for j in prime:
if (j ** 2 > i):
break
if (i % j == 0):
flag = False
break
if flag == True:
prime.append(i)
เมื่อจับเวลาก็จะเห็นความแตกต่าง
poomk@gemini:/tmp$ time ./prime2.py >0

real 0m27.190s
user 0m27.090s
sys 0m0.096s

poomk@gemini:/tmp$ time ./prime4.py >0

real 0m20.422s
user 0m20.369s
sys 0m0.040s

อนึ่งโปรแกรมตัวอย่างนี้เขียนขึ้นเพื่อความอ่านง่าย (และมักง่าย) มีจุดประสงค์เพื่อบอกแนวโน้มของความเร็วที่เพิ่มขึ้นของแต่ละไอเดีย ซึ่งหากเขียนไปใช้งานจริงๆอาจจะมีลูกเล่นที่ทำให้ทำเวลาได้ดีขึ้นกว่านี้อีก

No comments: