Merhaba,
Bu mesajda farklı integral teknikleri hakkında biraz bilgi paylaşayım istiyorum. Basit veya karmaşık bir çok farklı integral alma yöntemi var. Bu yöntemlerin daha basit olanlarından bahsetmek istiyorum. Daha karmaşık olanlarından bahsetmeyeceğim çünkü henüz anlamış değilim. Bahsedeceğim yöntemlerin ayrıca adaptif (uyarlanabilir) versiyonlarının algoritmalarını da paylaşacağım.
İlk olarak Monte Carlo entegrasyon yöntemiyle başlayabiliriz. Monte Carlo yönteminin genel algoritması şöyle özetlenebilir:
from math import pi
from random import uniform
def monte_carlo_integration(f, a, b, n):
return (b - a) / n * sum(f(uniform(a, b)) for _ in range(n))
Bu yöntemle pi sayısını kestirmeye çalışalım. Pi sayısının hesaplanmasıyla alakalı şöyle bir başlık da açılmış bu arada:
a = 0
b = 1
n = pow(10, 6)
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * monte_carlo_integration(f, a, b, n))
print(pi)
Çıktı:
3.14273191750064
3.141592653589793
Yukarda oluşturduğumuz monte carlo yönteminde, integralin hassasiyeti, üreteceğimiz nokta sayısına bağlıdır. Monte Carlo’nun adaptif versiyonunda ise integralin hassasiyeti, üretilecek nokta sayısına değil, bölümlere ayrılmış integrallere bağlıdır. Yukardaki monte carlo yöntemini, şöyle adaptif hale getirebiliriz:
from random import uniform
def adaptive_monte_carlo_integration(f, a, b, t):
def recursive_monte_carlo(f, a, b, t):
c = (a + b) / 2
left = (c - a) * f(uniform(a, c))
right = (b - c) * f(uniform(c, b))
area = left + right
if abs(area - (b - a) * f(uniform(a, b))) <= t:
return area
else:
return recursive_monte_carlo(f, a, c, t) + recursive_monte_carlo(f, c, b, t)
return recursive_monte_carlo(f, a, b, t)
Adaptif yöntemlerde, n
parametresi yerine (0, 1)
arasında olan ve tolerance anlamına gelen t
parametresini kullanıyoruz. Adaptif Monte Carlo yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
t = 1e-6
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * adaptive_monte_carlo_integration(f, a, b, t))
print(pi)
Çıktı:
3.1415897817240817
3.141592653589793
Adaptasyonu iki alana göre yapmak zorunda değiliz elbette. Üç alana göre de yapabilirsiniz.
def adaptive_monte_carlo__integration_v2(f, a, b, t):
def recursive_monte_carlo(f, a, b, t):
c = (2 * a + b) / 3
d = (2 * b + a) / 3
left = (c - a) * f(uniform(a, c))
mid = (d - c) * f(uniform(c, d))
right = (b - d) * f(uniform(d, b))
area = left + mid + right
if abs(area - (b - a) * f(uniform(a, b))) <= t:
return area
else:
return recursive_monte_carlo(f, a, c, t) + recursive_monte_carlo(f, c, d, t) + recursive_monte_carlo(f, d, b, t)
return recursive_monte_carlo(f, a, b, t)
İkinci Adaptif Monte Carlo yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
t = 1e-6
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * adaptive_monte_carlo_integration_v2(f, a, b, t))
print(pi)
Çıktı:
3.141592560560807
3.141592653589793
Başka bir entegrasyon yöntemini ele alalım. Örneğin, 0. dereceden polinom (diktörtgen) kullanarak, Riemann entegrasyonu uygulayabiliriz. Metodun algoritması şu şekildedir:
def riemann_integration(f, a, b, n):
return (h := (b - a) / n) * sum(f(a + i * h) for i in range(n))
Riemann yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
n = pow(10, 6)
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * riemann_integration(f, a, b, n))
print(pi)
Çıktı:
3.141594652413976
3.141592653589793
Şimdi de, Riemann’ın yöntemini adaptif hale getirelim:
def adaptive_riemann_integration(f, a, b, t):
def recursive_riemann(f, a, b, fa, fb, t):
c = (a + b) / 2
fc = f(c)
left = (c - a) * fa
right = (b - c) * fb
area = left + right
if abs(area - (b - a) * (fa + fb)) <= t:
return area
else:
return recursive_riemann(f, a, c, fa, fc, t) + recursive_riemann(f, c, b, fc, fb, t)
return recursive_riemann(f, a, b, f(a), f(b), t)
Adaptif Riemann yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
t = 1e-6
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * adaptive_riemann_integration(f, a, b, t))
print(pi)
Çıktı:
3.141592480560611
3.141592653589793
Birinci dereceden polinom kullanarak da entegrasyon yapabiliriz. Birinci dereceden polinom kullanabileceğimiz yöntemlerden birisi Trapezoid’tir. Onun da algoritması şu şekilde:
def trapezoid_integration(f, a, b, n):
return 1 / 2 * (h := (b - a) / n) * sum(
sum(j * f(a + (i + k) * h) for k, j in enumerate([1, 1])) for i in range(n)
)
Trapezoid yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
n = pow(10, 6)
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * trapezoid_integration(f, a, b, n))
print(pi)
Çıktı:
3.141592652413872
3.141592653589793
Pi sayısına daha fazla yaklaşmışız. Şimdi de adaptif bir trapezoid yöntemi oluşturalım:
def adaptive_trapezoid_integration(f, a, b, t):
def recursive_trapezoid(f, a, b, fa, fb, t):
c = (a + b) / 2
fc = f(c)
left = 1/2 * (c - a) * (fa + fc)
right = 1/2 * (b - c) * (fc + fb)
area = left + right
if abs(area - 1/2 * (b - a) * (fa + fb)) < t:
return area
else:
return recursive_trapezoid(f, a, c, fa, fc, t) + recursive_trapezoid(f, c, b, fc, fb, t)
return recursive_trapezoid(f, a, b, f(a), f(b), t)
Adaptif Trapezoid integrasyon yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
t = 1e-20
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * adaptive_trapezoid_integration(f, a, b, t))
print(pi)
3.1415926535897705
3.141592653589793
Pi’ye daha da yaklaşmışız.
Polinom derecesinin 2 olduğu durumda Milne kuralını uygulayabiliriz. Onun algoritması ise şöyle:
def milnes_rule_integration(f, a, b, n):
return 4 / 3 * (h := (b - a) / n) * sum(
sum(j * f(a + (i + k / 4) * h) for k, j in enumerate([2, -1, 2])) for i in range(0, n, 4)
)
Milne kuralıyla pi sayısını kestirmeye çalışalım:
a = 0
b = 1
n = pow(10, 6)
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * milnes_rule_integration(f, a, b, n))
print(pi)
Çıktı:
3.141599646107457
3.141592653589793
Milne’nin yöntemini adaptif hale getirelim:
def adaptive_milnes_rule_integration(f, a, b, t):
def recursive_miles_rule(f, a, b, fa, fb, t):
c = (a + b) / 2
fc = f(c)
d = (a + c) / 2
e = (c + b) / 2
fd = f(d)
fe = f(e)
left = 4/3 * (c - a) / 4 * (2 * fa + -1 * fd + 2 * fc)
right = 4/3 * (b - c) / 4 * (2 * fc + -1 * fe + 2 * fb)
area = left + right
if abs(area - 4/3 * (b - a) / 4 * (2 * fa + -1 * fc + 2 * fb)) <= t:
return area
else:
return recursive_miles_rule(f, a, c, fa, fc, t) + recursive_miles_rule(f, c, b, fc, fb, t)
return recursive_miles_rule(f, a, b, f(a), f(b), t)
Adaptif Milne kuralı integrasyon yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
t = 1e-20
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * adaptive_milnes_rule_integration(f, a, b, t))
print(pi)
3.1415926535897736
3.141592653589793
Bu yöntemle de pi’ye baya yaklaşmışız.
İkinci dereceden polinom kullanabileceğimiz başka bir entegrasyon yöntemi ise Simpson 1/3 kuralı. Onun algoritması da şöyle:
def simpson_1_3_integration(f, a, b, n):
return 1 / 3 * (h := (b - a) / n) * sum(
sum(j * f(a + (i + k) * h) for k, j in enumerate([1, 4, 1])) for i in range(0, n, 2)
)
Simpson 1/3 yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
n = pow(10, 6)
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * simpson_1_3_integration(f, a, b, n))
print(pi)
Çıktı:
3.1415926531305414
3.141592653589793
Simpson 1/3 kuralını adaptif hale getirelim:
def adaptive_simpson_1_3_integration(f, a, b, t):
def recursive_simpson_1_3(f, a, b, fa, fb, t):
c = (a + b) / 2
fc = f(c)
d = (a + c) / 2
e = (c + b) / 2
fd = f(d)
fe = f(e)
left = 1/3 * (c - a) / 2 * (fa + 4 * fd + fc)
right = 1/3 * (b - c) / 2 * (fc + 4 * fe + fb)
area = left + right
if abs(area - 1/3 * (b - a) / 2 * (fa + 4 * fc + fb)) <= t:
return area
else:
return recursive_simpson_1_3(f, a, c, fa, fc, t) + recursive_simpson_1_3(f, c, b, fc, fb, t)
return recursive_simpson_1_3(f, a, b, f(a), f(b), t)
Adaptif Simpson 1/3 integrasyon yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
t = 1e-20
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * adaptive_simpson_1_3_integration(f, a, b, t))
print(pi)
Çıktı:
3.141592653589793
3.141592653589793
Polinomun derecesini 3’e çıkardığımızda, Simpson 3/8 kuralıyla entegrasyon yapabiliriz. Onun algoritması da şöyle:
def simpson_3_8_integration(f, a, b, n):
return 3 / 8 * (h := (b - a) / n) * sum(
sum(j * f(a + (i + k / 3) * h) for k, j in enumerate([1, 3, 3, 1])) for i in range(0, n, 3)
)
Simpson 3/8 yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
n = pow(10, 6)
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * simpson_3_8_integration(f, a, b, n))
print(pi)
Çıktı:
3.141596654558362
3.141592653589793
Simpson 3/8 yöntemini adaptif hale getirelim:
def adaptive_simpson_3_8_integration(f, a, b, t):
def recursive_simpson_3_8(f, a, b, fa, fb, t):
c = (a + b) / 2
d = (2 * a + c) / 3
e = (2 * c + a) / 3
o = (2 * c + b) / 3
m = (2 * b + c) / 3
fc = f(c)
fd = f(d)
fe = f(e)
fo = f(o)
fm = f(m)
left = 3/8 * (c - a) / 3 * (fa + 3 * fd + 3 * fe + fc)
right = 3/8 * (b - c) / 3 * (fc + 3 * fo + 3 * fm + fb)
area = left + right
if abs(area - 3/8 * (b - a) / 3 * (fa + 3 * fe + 3 * fo + fb)) <= t:
return area
else:
return recursive_simpson_3_8(f, a, c, fa, fc, t) + recursive_simpson_3_8(f, c, b, fc, fb, t)
return recursive_simpson_3_8(f, a, b, f(a), f(b), t)
Adaptif Simpson 3/8 integrasyon yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
t = 1e-20
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * adaptive_simpson_3_8_integration(f, a, b, n))
print(pi)
Çıktı:
3.141592653589793
3.141592653589793
Bu da, iyi bir sonuç verdi.
Polinomun derecesini 4’e çıkartalım. 4. dereceden polinom kullanabileceğimiz entegrasyon yöntemi Boole’nin yöntemidir ve algoritması aşağıdaki gibidir:
def booles_rule_integration(f, a, b, n):
return 1 / 90 * (h := (b - a) / n) * sum(
sum(j * f(a + (i + k / 4) * h) for k, j in enumerate([7, 32, 12, 32, 7])) for i in range(n)
)
a = 0
b = 1
n = pow(10, 6)
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * booles_rule_integration(f, a, b, n))
print(pi)
Çıktı:
3.141592653539437
3.141592653589793
Şimdi de, Boole kuralını adaptif hale getirelim:
def adaptive_booles_rule_integration(f, a, b, t):
def recursive_booles_rule(f, a, b, fa, fb, t):
c = (a + b) / 2
d = (3 * a + c) / 4
e = (a + c) / 2
g = (3 * c + a) / 4
h = (3 * c + b) / 4
m = (c + b) / 2
o = (3 * b + c) / 4
fc = f(c)
fd = f(d)
fe = f(e)
fg = f(g)
fh = f(h)
fm = f(m)
fo = f(o)
left = 1/90 * (c - a) * (7 * fa + 32 * fd + 12 * fe + 32 * fg + 7 * fc)
right = 1/90 * (b - c) * (7 * fc + 32 * fh + 12 * fm + 32 * fo + 7 * fb)
area = left + right
if abs(area - 1/90 * (b - a) * (7 * fa + 32 * fe + 12 * fc + 32 * fm + 7 * fb)) <= t:
return area
else:
return recursive_booles_rule(f, a, c, fa, fc, t) + recursive_booles_rule(f, c, b, fc, fb, t)
return recursive_booles_rule(f, a, b, f(a), f(b), t)
Adaptif Boole integrasyon yöntemiyle pi sayısını kestirmeye çalışalım:
a = 0
b = 1
t = 1e-20
f = lambda x: pow(1 - pow(x, 2), .5)
print(4 * adaptive_booles_rule_integration(f, a, b, t))
print(pi)
Çıktı:
3.141592653589793
3.141592653589793
Daha başka entegrasyon yöntemleri de var elbette. Başta da dediğim gibi, bu algoritmalar daha karmaşık yapıdalar. Alan içerisine yerleştirilecek şeklin kaçıncı dereceden bir polinom olacağına bağlı olarak yöntem icat etmek mümkün. 1. dereceden bir polinomun katsayıları [1, 1], iken 2. dereceden bir polinomun katsayıları [1, 4, 1], 3. dereceden bir polinomun katsayıları [1, 3, 3, 1], 4. dereceden bir polinomun katsayıları ise [7, 32, 12, 32, 7] olur. Bu katsayılara Newton-Cotes katsayıları da deniyor.
Şimdilik paylaşacaklarımın sonuna geldik. Faydası olur diye umuyorum. Tekrar görüşmek üzere. Sağlıcakla kalın.