Python’da Sembolik Hesaplama İle Grafik Çizimi

Şu sayfada adına Glasser fonksiyonu denilen ve $$G(x):=\int_0^x\sin(t\sin t)\;dt$$ olarak tanımlanan bir fonksiyon gördüm, ilgili sayfada bu fonksiyonun bazı ilginç özellikleri anlatılıyor ve Mathematica yazılımı ile çizilmiş grafikleri yer alıyor. Sayfada yer alan bazı geometrik özelliklerini daha yakından gözlemlemek için grafiğini ben de Python ile çizmeye çalıştım, bu yazımda bunu nasıl yaptığımı açıklayacağım.

Bildiğiniz gibi Python’da bir fonksiyonun grafiğini çizmek çoğu dilde olduğu gibi son derece basit bir fikre dayanır, istediğiniz aralıkta çok sayıda nokta belirlersiniz ve bu noktaların her birinde $(x, f(x))$ çiftlerini hesaplarsınız (bu işleme vektörizasyon denir). Sonra da düzlemde bu koordinatları düz doğru parçalarıyla birleştirirsiniz, eğer yeterince çok sayıda nokta alırsanız köşeler gözle görülmez ve düzgün bir grafiğe sahip olursunuz. Bu işlemler çok sık yapıldığı için NumPy ve Matplotlib gibi kütüphanelerde bunlar için hazır fonksiyonlar vardır, bu konuları şu sayfadaki ders notlarımda detaylıca açıklamıştım (ders notlarında Python 2 için anlatmıştım ama versiyon 3 için de ilgili notları kullanabilirsiniz.)

Fakat Glasser fonksiyonu için durum biraz daha karmaşık, vektörizasyon için her bir $x$ noktasında $G(x)$ değerini hesaplamanız gerekiyor. Bu integral değerlerini hesaplama işlemi elemanter fonksiyonlar cinsinden yapılamadığı için sembolik hesaplama araçlarından yardım almak zorundayız.

Benim aklıma ilk gelen SymPy oldu, aşağıdaki kodları çalıştırmadan önce yazının devamını okumanızı tavsiye ederim.

from sympy import symbols, sin, Integral
from sympy.plotting import plot
x, t = symbols('x t')

G = Integral(sin(t*sin(t)),(t,0,x))
plot(G, (x,0,10), show=True)

Yukarıdaki kodlar ile $[0,10]$ aralığında grafik çizdiriliyor ancak hesaplama oldukça yavaş oluyor, aralığı büyüttüğünüz zaman verimsizliği çok daha belirgin olarak gözlemleyebilirsiniz. Bu kodların çıktısı aşağıdaki grafik oldu.

Daha sonra alternatif bir yaklaşım olarak SciPy kütüphanesini kullanmayı düşündüm, bu kütüphanede SymPy’de olduğu gibi genel amaçlı integrasyon için kullanılabilen integrate.quad fonksiyonu var.

import numpy as np
from scipy import integrate
import matplotlib.pyplot as plt

def f(x):
    return np.sin(x*np.sin(x))

X = np.arange(0,10,0.001)

def F(x):
    return integrate.quad(f,0,x)

plt.plot(X, [F(x)[0] for x in X])

Bu programın da grafik üretmesi hatırı sayılır zaman aldı ama SymPy’ye göre daha hızlı olduğunu söyleyebilirim, aşağıdaki çıktıyı aldım ki SymPy ile aynı sonucu vermesi doğal.

Sonra aklıma birden nümerik integrasyon geldi, SciPy içinde bulunan integrate.cumtrapz fonksiyonunu kullanarak ilgili integrali standart bir hata payı ile yamuklar kuralını kullanarak yaklaşık olarak hesaplattım. Bu sayede vektörizasyon aşaması dramatik olarak hızlandı ve geniş aralıklarda bile çok hızlı sonuç alabildim.

import numpy as np
from scipy import integrate
import matplotlib.pyplot as plt

x = np.linspace(0, 20, num=1000)
y = np.sin(x*np.sin(x))
y_int = integrate.cumtrapz(y, x, initial=0)
plt.plot(x, y_int)
plt.show()