Ödev: 2

FİZ220 - Bilgisayar Programlama II | 10/05/2020

Matrisler, Operatörler

Son gönderim tarihi: 20 Mayıs Perşembe, 23:59
Gönderim şekli: FIZ220_Odev_02Grup#.ipynb isimli jupyter ipynb formatında dosyayı ödev sayfasından göndermek suretiyle
Gönderecek kişi: Grup temsilcisi

Dr. Emre S. Taşcı, emre.tasci@hacettepe.edu.tr
Fizik Mühendisliği Bölümü
Hacettepe Üniversitesi

Matrisler, operatörler

1. Soru: Dönüş operatörü

1. x-ekseninde, x=1, y=0 noktasında başlıyoruz. (1000x2)'lik bir noktalar matrisi tanımlıyoruz.
2. $\theta \in [0,360)$ aralığında rasgele bir $\theta$ değeri seçip ($\theta$'nın tam sayı olması gerekmiyor), noktamızı o kadar döndürüyoruz.
3. Geldiği noktanın x ve y değerlerini noktalar matrisinde kaydediyoruz: matrisin her bir satırı, o elde gelinen koordinatı gösteriyor.
4. 2-3 adımlarını 1000 kere tekrarlıyoruz (noktanın bir önceki elde geldiği noktadan devam ederek)
5. noktalar matrisinin x ve y değerlerini tuttuğumuz 0 ve 1 sütunlarının ortalamalarını ekrana yazdırıyoruz
6. 1-5 adımlarını 30 kere tekrarlıyoruz.

Çözüm

import numpy as np

#6
for j in range(30):
    
    #1
    P = np.array([1,0])
    noktalar = np.empty((1000,2),float)
    
    #4
    for i in range(1000):
        #2
        theta = np.random.rand()*(2*np.pi)
        R = np.array([[np.cos(theta),-np.sin(theta)],[np.sin(theta),np.cos(theta)]])
        Q = np.dot(R,P)
        
        #3
        noktalar[i,:] = Q
        
        # Donmeye su anda bulundugumuz noktadan
        # devam etmek icin P noktasinin konumunu
        # Q'nunki olarak guncelliyoruz
        P = Q
        
    #5
    print (np.mean(noktalar,axis=0))

[-0.01214762 -0.03941625]
[-0.01455028  0.0126116 ]
[ 0.05055427 -0.02501256]
[0.04292054 0.03274785]
[-0.00951909 -0.00797527]
[-0.02275414  0.04071829]
[ 0.05948083 -0.03073184]
[-0.00592383  0.03216359]
[0.04318133 0.00341541]
[-0.01709916  0.02287009]
[ 0.01488003 -0.02979313]
[ 0.01337292 -0.02128784]
[ 0.02727008 -0.03598392]
[0.02444532 0.00101701]
[-0.00587281  0.04122044]
[-0.0194059   0.00888909]
[-0.00286285 -0.01168498]
[-0.02369914  0.01144   ]
[-0.03302397 -0.00444223]
[-0.02185576 -0.02884876]
[-0.01302587  0.05396949]
[-0.00421027 -0.00735289]
[-0.04255289 -0.00144056]
[-0.00931617  0.0060035 ]
[-0.01382919  0.0035068 ]
[0.03184999 0.0200361 ]
[-0.01408975 -0.05669687]
[ 0.03798385 -0.02955132]
[-0.00975376  0.00644513]
[-0.01711323 -0.0203906 ]

Bu alıştırmadaki amaç, düzgün dağılımın gayet düzgün çalıştığını göstermekti. Döndürmeler sonucunda ulaşılan yeni noktalar koordinat sistemine o kadar güzel yayılmış ki, ortalamaları hemen hemen (0,0) çıkıyor. Her seferinde aynı (1,0) noktasını rastgele bir açıyla çeviriyor olsaydık belki bu normal gelecekti (dönüş açılarının $[0,360^o)$ aralığında olmasından ötürü) ama burada her seferinde kaldığımız yerden yeni bir dönüş yapıyoruz, ona rağmen, gayet güzel, sağlıklı, düzgün bir dağılıma sahip oluyoruz. Hazır ortalamalardan bahsediyorken, bir de "Merkezi Limit Teoremi" (Central Limit Theorem) adı verilen, "dağılımların dağılımları" ile ilgilenen müthiş bir teorem var ki, ona da FİZ371 - Bilimsel ve Teknik Hesaplamalar dersinde ucundan değiniyoruz. 8)

2. Soru: 65.419 (veya -114.581) özdeğerse, 145.544 ne? Üveydeğer mi? 8~(

$$A=\begin{pmatrix}1&2\\3&4\end{pmatrix}$$

matrisini ele alalım: bu matrisin iki tane özdeğeri, onlara karşılık gelen de iki tane özvektörü var. Derste öğrendiğimiz eig() fonksiyonu ile bunları bulalım:

import numpy as np
A = np.matrix([[1,2],[3,4]])
[ozdegerler, ozvektorler] = np.linalg.eig(A)
print(ozdegerler,"\n",ozvektorler)

[-0.37228132  5.37228132] 
 [[-0.82456484 -0.41597356]
 [ 0.56576746 -0.90937671]]

yani:

$$\lambda_1 = -0.372\leftrightarrow u_1 = \begin{bmatrix}-0.825\\0.566\end{bmatrix}$$$$\lambda_2 = 5.372\leftrightarrow u_2 = \begin{bmatrix}-0.416\\-0.909\end{bmatrix}$$

1. Her bir özvektörün x-ekseni ile yaptığı açıyı hesaplayın. (145.544 ve -114.581 derece bulacaksınız - buraya kadar bir sıkıntı yok 8)

Çözüm

import numpy as np

A = np.matrix([[1,2],[3,4]])
[ozdegerler, ozvektorler] = np.linalg.eig(A)

#0
lambda_1 = ozvektorler[:,0]
lambda_2 = ozvektorler[:,1]

aci_1 = np.rad2deg(np.arctan2(lambda_1[1],lambda_1[0]))[0,0]
aci_2 = np.rad2deg(np.arctan2(lambda_2[1],lambda_2[0]))[0,0]

print("Açı 1: {:8.3f} derece".format(aci_1))
print("Açı 2: {:8.3f} derece".format(aci_2))

Açı 1:  145.544 derece
Açı 2: -114.581 derece

Özvektörleri yine derste öğrendiğimiz "vura vura buldurma" metoduyla bulmaya çalışalım:

1. Bileşenleri ${x,y}\in [-1,1]$ aralığında olan iki boyutlu rasgele bir k vektörü oluşturun (bileşenler tamsayı olmak zorunda değil)
2. Bu vektörün yatayla yaptığı açıyı hesaplayın
3. Vektörün üzerine A matrisini etki ettirin, ortaya çıkan vektörün yatayla yaptığı açıyı hesaplayın -- bu açıyla, bir önceki hesapladığınız açı arasındaki fark 0.0001 dereceden küçük olana kadar bu adımı tekrarlayın (kaç adımda bu küçük farka ulaştığınızı ekrana yazdırın). Yani k vektörüne matrisimizi etki ettirdik, diyelim ki k' vektörünü elde ettik, aralarındaki açı 0.0001 dereceden büyükse, matrisi k' vektörüne vurdurup, bu sefer k'' vektörünü bulduk. Şimdi açı kontrolünü k' ile k'' arasında yapıyoruz: yani hep şimdi hesapladığımız açı ile bir önceki hesapladığımız açı arasında oluyor bu kontrol.
4. En son elde ettiğiniz vektörün boyunu 1 olacak şekilde ayarlayın, yönü aynı kalsın (anahtar terim: "normalize etmek" 8) Açıyı ve vektörü yazdırın.
5. 1-4 arasını 20 kere tekrarlayın. Bulduğunuz açı 145.544 veya -114.581'den biri mi? Değilse niye değil, durumu kurtarabiliyor musunuz? ;)

Ultra Bonus Soru: 1-4 adımlarını 1000 kere tekrarlayın -- her seferinde aynı özvektöre mi gitti? Eğer öyle ise bu metotla diğerine nasıl ulaşırız?

import numpy as np

A = np.matrix([[1,2],[3,4]])
[ozdegerler, ozvektorler] = np.linalg.eig(A)

#1
k = np.random.rand(2)*2 - 1
k = k.reshape(2,1)
print("k_ilk:",k)

## Kontrol icin rastgele k vektoru yerine
## elimizle tanimli degerler almak istersek:
## k=np.array([[1],[0]])
## k=np.array([[-1],[0]])

#2 
aci = np.rad2deg(np.arctan2(k[1,0],k[0,0]))
print("Açı: {:.3f} derece".format(aci))

print("-"*45)

#3
aci_epsilon = 1E-4
aci_fark = 100
while(aci_fark > aci_epsilon):
    k_ussu = np.dot(A,k)
    aci_ussu = np.rad2deg(np.arctan2(k_ussu[1,0],k_ussu[0,0]))
    aci_fark = np.abs(aci - aci_ussu)
    k = k_ussu
    aci = aci_ussu

#4
k_boy = np.linalg.norm(k)
k_normalize = k/k_boy
print("k_son:",k_normalize)

aci = np.rad2deg(np.arctan2(k_normalize[1,0],k_normalize[0,0]))
print("Açı: {:.3f} derece".format(aci))

# Buradaki sonuc 65.419 veya -114.581 derecelerden biri cikacak.
# Kontrollu olarak elde etmek isterseniz ilk k vektorunu 
# rastlantisal olarak degil de, 2. adim oncesinde:
# k=np.array([[1],[0]])  olarak tanimladiginizda 65.419 dereceyi
# k=np.array([[-1],[0]]) olarak tanimladiginizda -114.581 dereceyi 
# bulursunuz.
#   65.419 derece : [ 0.41597355,  0.90937671]
# -114.581 derece : [-0.41597355, -0.90937671]

k_ilk: [[0.31510986]
 [0.38840985]]
Açı: 50.948 derece
---------------------------------------------
k_son: [[0.41597358]
 [0.9093767 ]]
Açı: 65.419 derece

#5
Yukarıdaki kodu kaç kere çalıştırırsanız çalıştırın sonuç olarak %99.999999... ihtimalle şu iki vektör (ve açıdan) birini bulacaksınız:

$$a_1 = \begin{bmatrix}{0.416\\ 0.909}\end{bmatrix}, \theta_1 = 65.419^o$$$$a_2 = \begin{bmatrix}{-0.416\\ -0.909}\end{bmatrix}, \theta_2 = -114.581^o$$

İkinci vektörü zaten 0. adımdan yaptığımız hesaptan biliyorduk. Buradaki iki açının arasındaki farkın 180 derece olduğu, buna ek olarak da iki vektörün birbirlerinin zıttı oldukları da bize bir şey anlatıyor: bu vektörlerin aynı doğru üzerinde olduklarını. Yani bu iki vektör iki farklı özvektör değil. Bildiğiniz üzere $\vec{u}$ bir özvektörse, herhangi bir skalerle çarpılmış hali de aynı özdeğerli özvektördür, yani matrisimizin etkimesi altında dönmez, boyu da özdeğer kadar ölçeklenir (bu gerçeği "Bariz bilgi" başlığıyla "Lineer Cebir" uygulama dersimizde ispatlamıştık).

import matplotlib.pyplot as plt

plt.arrow(0,0,0.416,0.909,\
          color="blue",\
          head_width=0.2,\
          length_includes_head=True)
plt.arrow(0,0,-0.416,-0.909,\
          color="red",\
          head_width=0.2,\
          length_includes_head=True)
plt.xlim(-1.2,1.2)
plt.ylim(-1.2,1.2)
plt.title("ikimiz bir fidanin +/- sonsuza giden daliyiz...")
plt.show()

Demek ki bu metotla aslında aynı özvektörü bulup duruyormuşuz. 145.544 derecelik açı yapan $$\lambda_1 = -0.372\leftrightarrow u_1 = \begin{bmatrix}-0.825\\0.566\end{bmatrix}$$ $\lambda_1$ özdeğerli $u_1$ özvektörü bu metotla gerçekten de üvey evlat muamelesi görüyor.

Son olarak bizim bulduğumuz ve bulamadığımız iki özvektör yön doğrusunu hatıra olarak yan yana çizelim:

import numpy as np
import matplotlib.pyplot as plt

aci1 = 145.544
aci2 = -114.581

x = np.linspace(-2,2,2)
y1 = x*np.tan(np.deg2rad(aci1))
y2 = x*np.tan(np.deg2rad(aci2))

plt.plot(x,y1,"b-")
plt.plot(x,y2,"r-")
plt.xlim(-4,4)
plt.ylim(-4,4)
plt.show()

Grafikte kırmızı doğru bizim bu metotla bulduğumuz özvektör yönünü (x ekseniyle -114.581 = 65.419 derecelik açı yapan), mavi doğru ise bulamadığımız yönü (x ekseniyle 145.544 = -34.456 derecelik açı yapan) göstermekte.

Ultra bonus kısmı

Bu notkadan sonrası, işlerin biraz ciddiye bindiği kısım, o nedenle sadece ilgisi ve merakı olanların devam etmesini tavsiye ederim, zira kafa karışıklığı yaratabilir.

(2x2)'lik $A$ matrisinin, iki boyutlu bir düzlemdeki vektörleri onların değerlerine (boy & yön) göre çevirip ölçeklendireceğini biliyoruz. Eğer bu vektörler yukarıdaki grafikte işaretlediğimiz doğrulara paralel iseler, sadece boyları ölçekleniyor (kırmızı doğrunun üzerindeyse boyu $\lambda_2 = 5.372$ kadar uzuyor; mavi doğrunun üzerindeyse bu sefer de $\lambda_1 = -0.372$ kadar (yani önce yönü zıt yön oluyor (-), ardından da boyu 0.372 oranında kısalıyor).

Vektörleri -yönlerini ve boylarını değiştirmediğimiz sürece- istediğimiz gibi taşıyabiliriz. O halde gelin başlangıç noktaları orijin olacak şekilde hepsini taşıyalım:

for x in np.arange(-4,5,0.5):
    for y in np.arange(-4,5,0.5):
        plt.arrow(0,0,x,y,\
                  length_includes_head=True,\
                  head_width=.3,\
                  head_length=.5)

plt.xlim(-6,6)
plt.ylim(-6,6)
plt.show()

Böyle de çok karışık oldu. Zaten orijinden başladıklarını bildiğimizden, sadece uçlarını işaretleyelim:

(önce vektörler + uç noktaları yapalım, ardından sadece uç noktalarını gösterelim)

for x in np.arange(-4,5,0.5):
    for y in np.arange(-4,5,0.5):
        plt.arrow(0,0,x,y,\
                  length_includes_head=True,\
                  head_width=.3,\
                  head_length=.5)
        plt.plot(x,y,"or")

plt.xlim(-6,6)
plt.ylim(-6,6)
plt.show()

for x in np.arange(-4,5,0.5):
    for y in np.arange(-4,5,0.5):
        plt.plot(x,y,"or")

plt.xlim(-6,6)
plt.ylim(-6,6)
plt.show()

(Bu kırmızı noktaların her birinin bir vektörü temsil ettiğini aklımızda tutarak), $A$ operatörünü bu vektörlere vurduralım ve çıkan vektörün boyunu ve yönünü kullanarak, vektörümüzün ne yöne nasıl bir dönüşüm geçireceğini oklarla gösterelim:

A = np.array([[1,2],[3,4]])

vec_list = []
for x in np.arange(-4,5,.5):
    for y in np.arange(-4,5,.5):
        plt.plot(x,y,"or")
        v = np.array([x,y])
        vp = np.dot(A,v)
        vp_n = np.linalg.norm(vp)
        vec_list.append([v,vp,vp_n])
        #plt.arrow(x,y,(vp[0]-x,vp[1]-y)/20)

vec_list = np.array(vec_list,dtype=object)
norm_max = max(vec_list[:,2])
#print(norm_max)
#norm_max = np.sqrt(norm_max)

vec_list[:,1] /= norm_max /2

#plt.show()
#vec_list

for row in vec_list[:]:
    #print(row[0]+row[1])
    plt.arrow(row[0][0],row[0][1],row[1][0],row[1][1],\
              length_includes_head=False,\
              head_width=.3,\
              head_length=.5)
plt.show()

Hazır yapıyorken, bizim özvektörlerimizin yönlerini gösteren doğruları da grafiğe ekleyelim:

A = np.array([[1,2],[3,4]])

vec_list = []
for x in np.arange(-4,5,.5):
    for y in np.arange(-4,5,.5):
        plt.plot(x,y,"or")
        v = np.array([x,y])
        vp = np.dot(A,v)
        vp_n = np.linalg.norm(vp)
        vec_list.append([v,vp,vp_n])
        #plt.arrow(x,y,(vp[0]-x,vp[1]-y)/20)

vec_list = np.array(vec_list,dtype=object)
norm_max = max(vec_list[:,2])
#print(norm_max)
#norm_max = np.sqrt(norm_max)

vec_list[:,1] /= norm_max /2

#plt.show()
#vec_list

for row in vec_list[:]:
    #print(row[0]+row[1])
    plt.arrow(row[0][0],row[0][1],row[1][0],row[1][1],\
              length_includes_head=False,\
              head_width=.3, \
              head_length=.5)

aci1 = 145.544
aci2 = -114.581

x = np.linspace(-6,6,2)
y1 = x*np.tan(np.deg2rad(aci1))
y2 = x*np.tan(np.deg2rad(aci2))

plt.plot(x,y1,"b-")
plt.plot(x,y2,"r-")
plt.xlim(-8,8)
plt.ylim(-8,8)
plt.show()

Bu da epey karmaşık görünüyor -- en iyisi sadece kırmızı noktaların olduğu yeri gösteren konum vektörlerinin $A$ operatörünün vurması sonucu gidecekleri yön ve -orantılı- mesafeyi gösteren oklarla özvektör doğrularını gösterelim:

A = np.array([[1,2],[3,4]])

vec_list = []
for x in np.arange(-4,5,.5):
    for y in np.arange(-4,5,.5):
        #plt.plot(x,y,"or")
        v = np.array([x,y])
        vp = np.dot(A,v)
        vp_n = np.linalg.norm(vp)
        vec_list.append([v,vp,vp_n])
        #plt.arrow(x,y,(vp[0]-x,vp[1]-y)/20)

vec_list = np.array(vec_list,dtype=object)
norm_max = max(vec_list[:,2])
#print(norm_max)
#norm_max = np.sqrt(norm_max)

vec_list[:,1] /= norm_max /2

#plt.show()
#vec_list

for row in vec_list[:]:
    #print(row[0]+row[1])
    plt.arrow(row[0][0],row[0][1],row[1][0],row[1][1],\
              length_includes_head=False,\
              head_width=.3,\
              head_length=.5)

aci1 = 145.544
aci2 = -114.581

x = np.linspace(-6,6,2)
y1 = x*np.tan(np.deg2rad(aci1))
y2 = x*np.tan(np.deg2rad(aci2))

plt.plot(x,y1,"b-")
plt.plot(x,y2,"r-")
plt.xlim(-8,8)
plt.ylim(-8,8)
plt.show()

Bir daha hatırlatmak gerekirse, kırmızı doğru bizim bu ödevde kullandığımız metotla bulduğumuz özvektör yönü; mavi doğru ise bu metotla ulaşamadığımız yönü gösteriyor. Gerçekten de tam mavi doğrunun üzerinde olmayan herhangi bir noktadan başladığımızda oklar bizi kırmızı çizgiye (vektör yönlerinden bahsettiğimizden "kırmızı çizgiyle paralel yönlere" diyoruz aslında) doğru sürüklüyor, bu nedenle kendimizi kırmızı çizginin üzerinde buluveriyoruz (ayrıca yukarıdaki grafikten açıkça görüldüğünü umut ettiğim üzere, mavi çizginin yukarısında bir noktadan (örneğin: (0,1)) başladığımızda kırmızı çizginin 65.419 derecelik açı yapan kısmına (yani yukarı doğru giden kısma); mavi çizginin altında bir noktadan (örneğin: (-1,0)) başladığımızda ise kırmızı çizginin -114.581 derecelik açı yapan kısmına (yani aşağıya doğru giden kısma) denk gelmekteyiz -- metodumuzun başlangıç değerine göre 65.419 veya -114.581 derece vermesinin sebebi bu.

Sonuçta tam olarak mavi çizginin üzerinde, yani birinci özvektörün yönünde bir yerde başlamadığımız sürece "dalgalar" bizi sürükleyecek, kendimizi kırmızı çizginin, yani ikinci özvektörün yönünde bulacağız.

Az evvel yapmış olduğumuz "dalgalar" benzetmesi güzel bir analoji olabilir ama daha iyisi var: üzerinde ok figürlerinin olduğu bir masa örtümüz olduğunu düşünün. Bunu mavi çizgi boyunca sıradağların, kırmızı çizgi boyunca da vadilerin/kanyonların/hendeklerin/çukrların olduğu bir yere örttüğümüzü aklımızda canlandıralım. Normal durumda masa örtümüzün üzerindeki oklar birbirlerine paraleller (örtünün ortasının üstündekiler yukarı doğru, altındakiler aşağı doğru çizilmişler). Örtümüzün orta çizgisini mavi çizginin üzerine gelecek şekilde serip, yukarıdan baktığımızda görürüz ki, altındaki yükseklikten ötürü, mavi çizginin etrafındaki oklar birbirinden ayrık dururken, kırmızı çizginin etrafındaki oklar giderek birbirine yaklaşırlar.

Yani mavi çizginin (yüksek bir tepenin) yanından bir top bıraksak, bu top gidip gidip nihayetinde kırmızı çizginin (vadinin) içine girecektir.

Peki durum buysa, mavi çizgiyi hiç mi bulma şansımız yok? Tabii ki var: ters yüz etmek! Yani tepeleri çukur, çukurları tepe olarak tanımlarsak (matematiksel olarak okların yönünü tersine gösterirsek) bu işi yapmış oluruz:

A = np.array([[1,2],[3,4]])

vec_list = []
for x in np.arange(-4,5,.5):
    for y in np.arange(-4,5,.5):
        #plt.plot(x,y,"or")
        v = np.array([x,y])
        vp = -np.dot(A,v)
        vp_n = np.linalg.norm(vp)
        vec_list.append([v,vp,vp_n])
        #plt.arrow(x,y,(vp[0]-x,vp[1]-y)/20)

vec_list = np.array(vec_list,dtype=object)
norm_max = max(vec_list[:,2])
#print(norm_max)
#norm_max = np.sqrt(norm_max)

vec_list[:,1] /= norm_max /2

#plt.show()
#vec_list

for row in vec_list[:]:
    #print(row[0]+row[1])
    plt.arrow(row[0][0],row[0][1],row[1][0],row[1][1],\
              length_includes_head=False,\
              head_width=.3, \
              head_length=.5)

aci1 = 145.544
aci2 = -114.581

x = np.linspace(-6,6,2)
y1 = x*np.tan(np.deg2rad(aci1))
y2 = x*np.tan(np.deg2rad(aci2))

plt.plot(x,y1,"b-")
plt.plot(x,y2,"r-")
plt.xlim(-5,5)
plt.ylim(-5,5)
plt.show()

Artık bütün oklar mavi çizgiye getiriyor! 8)

Bu dediğimizi matematiksel olarak nasıl ifade edebiliriz? Elimizde bir $\vec{a}$ vektörü olsun. Bunu $A$ matrisiyle çarptığımızda kırmızı doğrultuya doğru döndürülüp uzatılacak bu yeni vektörü de $\vec b$ ile gösterelim.

Örneğin, $\vec a = \begin{bmatrix}1.0\\0.5\end{bmatrix}$ için:

a = np.array([[1],[0.5]])
b = np.dot(A,a)

plt.arrow(0,0,a[0,0],a[1,0],\
              color="blue",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)
plt.arrow(0,0,b[0,0],b[1,0],\
              color="red",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)
aci1 = 145.544
aci2 = -114.581

x = np.linspace(-6,6,2)
y1 = x*np.tan(np.deg2rad(aci1))
y2 = x*np.tan(np.deg2rad(aci2))

plt.plot(x,y1,"k--")
plt.plot(x,y2,"k:")

plt.xlim(-1,6)
plt.ylim(-1,6)
plt.show()

Burada siyah "eksenler" birinci ve ikinci özvektör doğrultularını göstermektedir (çizgi-çizgi: birincil (145.544 derecelik) doğrultu; nokta-nokta: ikincil (-114.581 derecelik) doğrultu).

Mavi $\color{blue}{\vec a}$ vektörümüzü kırmızı $\color{red}{\vec b}$ vektörüne götüren değişikliğe de $\color{green}{\vec c}$ vektörü deyip, yeşille gösterelim:

a = np.array([[1],[0.5]])
b = np.dot(A,a)
print(a.shape)
plt.arrow(0,0,a[0,0],a[1,0],\
              color="blue",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)
plt.arrow(0,0,b[0,0],b[1,0],\
              color="red",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)
# plt.arrow komutu parametre olarak başlangıç noktası ile
# x ve y yönlerindeki boyu alır
plt.arrow(a[0,0],a[1,0],b[0,0]-a[0,0],b[1,0]-a[1,0],\
              color="green",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)

x = np.linspace(-6,6,2)
y1 = x*np.tan(np.deg2rad(aci1))
y2 = x*np.tan(np.deg2rad(aci2))

plt.plot(x,y1,"k--")
plt.plot(x,y2,"k:")

plt.xlim(-1,6)
plt.ylim(-1,6)
plt.show()

(2, 1)

Grafikten de açıkça görüldüğü üzere: $$\color{blue}{\vec a} + \color{green}{\vec c} = \color{red}{\vec b}$$ olur.

"Tepeleri çukura, çukurları tepelere çevirmek" demek, mavi okun ucundan bıraktığımız topumuzun yeşil ok doğrultusunda yukarı doğru kayması yerine zıt yönde kayması anlamına geliyor. Yani:

(yukarıdaki denklemden) $$\color{green}{\vec c} = \color{red}{\vec b} - \color{blue}{\vec a}$$

$\color{green}{\vec c}$ vektörünü tersine çevirip $\color{blue}{\vec a}$ vektörüne eklemek demek:

$$\vec a - \vec c = \vec a - (\vec b - \vec a) = 2\vec a - \vec b$$

a = np.array([[1],[0.5]])
b = np.dot(A,a)

plt.arrow(0,0,a[0,0],a[1,0],\
              color="blue",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)
plt.arrow(0,0,b[0,0],b[1,0],\
              color="red",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)
# plt.arrow komutu parametre olarak başlangıç noktası ile
# x ve y yönlerindeki boyu alır
plt.arrow(a[0,0],a[1,0],b[0,0]-a[0,0],b[1,0]-a[1,0],\
              color="green",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)
c = b - a
plt.arrow(a[0,0],a[1,0],-c[0,0],-c[1,0],\
              color="purple",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)


aa = 2*a - b
plt.arrow(0,0,aa[0,0],aa[1,0],\
              color="brown",\
              length_includes_head=True,\
              head_width=.3, \
              head_length=.5)
x = np.linspace(-6,6,2)
y1 = x*np.tan(np.deg2rad(aci1))
y2 = x*np.tan(np.deg2rad(aci2))

plt.plot(x,y1,"k--")
plt.plot(x,y2,"k:")

plt.xlim(-6,6)
plt.ylim(-6,6)
plt.show()

Görüldüğü üzere, bu yaklaşım başlangıç vektörümüzü 2. özvektör doğrultusu yerine 1. özvektör doğrultusuna yöneltiyor.

İşimiz neredeyse bitti ama atlatmamız gereken son bir tehlike var: Yönlendirmeyi tersine aldık ama ölçeklendirme genel olarak hala arttırıcı şekilde davranacak (2. özvektör doğrultusuna yaklaştıkça ölçeklendirme de o özvektörün özdeğerine (5.372 katı) yaklaşıyordu). Bu nedenle değişim miktarının boyunu her seferinde $\vec{b}$'nin boyunu yarılayarak kontrol altında tutuyoruz.

Bu bulgularımızı uyguladığımızda, metodumuz artık diğer özvektöre yakınsayıp onu bulabiliyor!

import numpy as np

A = np.matrix([[1,2],[3,4]])
[ozdegerler, ozvektorler] = np.linalg.eig(A)

#1
k = np.random.rand(2)*2 - 1
k = k.reshape(2,1)
print("k_ilk:",k)

#2 
aci = np.rad2deg(np.arctan2(k[1,0],k[0,0]))
print("Açı: {:.3f} derece".format(aci))

print("-"*45)

#3
aci_epsilon = 1E-4
aci_fark = 100
while(aci_fark > aci_epsilon):
    b = np.dot(A,k)
    k_ussu = 2*k - b/2
    aci_ussu = np.rad2deg(np.arctan2(k_ussu[1,0],k_ussu[0,0]))
    #print("Açı: {:.3f} derece".format(aci_ussu))
    aci_fark = np.abs(aci - aci_ussu)
    k = k_ussu
    aci = aci_ussu

#4
k_boy = np.linalg.norm(k)
k_normalize = k/k_boy
print("k_son:",k_normalize)

aci = np.rad2deg(np.arctan2(k_normalize[1,0],k_normalize[0,0]))
print("Açı: {:.3f} derece".format(aci))

k_ilk: [[-0.90657801]
 [ 0.76704902]]
Açı: 139.766 derece
---------------------------------------------
k_son: [[-0.82456501]
 [ 0.56576722]]
Açı: 145.544 derece