Metode Penalti

-

Metode Penalti

--Catatan kuliah MA3071 20210802

Minimumkan $𝑓(𝐱)$ terhadap kendala $𝐱∈𝑆$
Dengan $𝑓$ merupakan fungsi yang kontinu pada $ℝ^𝑛$ dan $𝑆$ adalah himpunan kendala $𝑆=\{𝑥:h_i(𝐱)= 0,i=1,2,...,m\}$

  • Ide dari Metode penalti adalah mengubah masalah optimisasi berkendala menjadi masalah optimisasi tanpa kendala,
  • Kendala tersebut dikalikan dengan suatu konstanta positif (yang membesar) sehingga berfungsi menjadi suatu penalti pada masalah meminimumkan.

Masalah dimodifikasi menjadi :

Minimumkan $𝑓(𝐱)+𝑐𝑃(𝐱)$
Dengan $𝑐$ suatu konstanta positif dan $𝑃(𝐱)$ adalah fungsi pada $ℝ^𝑛$ yang memenuhi

  • $𝑃(𝐱)$ fungsi kontinu
  • $𝑃(𝐱)≥0$ untuk semua $𝐱∈ℝ^𝑛$
    • Biasanya dipilih fungsi penalti $\displaystyle 𝑃(𝐱)=\frac{1}{2} \sum_{i=1}^{m}[ℎ_i(\mathbf{𝑥})]^2$
  • $𝑃(𝐱)=0$ jika dan hanya jika $𝐱∈𝑆$

Contoh

Minimumkan $$f(x_1,x_2)=x_1^2+x_2^2$$ Terhadap kendala $x_1+x_2-1=0$

=======================================================

Dengan metode penalti, masalah menjadi dimodifikasi menjadi,
Minimumkan : $$\theta(\mathbf{x})=x_1^2+x_2^2+\frac{c}{2}\left(x_1+x_2-1\right)^2$$ Dengan $c$ merupakan konstanta positif yang besar $$\nabla \theta(\mathbf{x})=\begin{bmatrix}2x_1+c\left(x_1+x_2-1\right)\\2x_2+c\left(x_1+x_2-1\right)\end{bmatrix}$$ $$\nabla^2 \theta(\mathbf{x})=\begin{bmatrix}2+c&c\\c&2+c\end{bmatrix}$$ Selesaikan masalah optimisasi nonlinear tanpa kendala ini secara numerik

readme first :
-- Operasi aljabar pada python : https://www.ridwanreza.com/2021/07/aljabar-python.html
-- Metode Newton dan Conjugate Gradien : https://www.ridwanreza.com/2021/07/ma3071-algoritma-newton-conjugate.html

In [1]:
import numpy as np

Selesaikan dengan Metode Newton

In [2]:
#Metode Newton

x=np.array([0,0]) #nilai/tebakan awal
print("x0= ",x)
c=1000000 #coba untuk c=1,c=2, c=10, c=100, c=1000000
nablaf=np.array([2*x[0]+c*(x[0]+x[1]-1)  ,  2*x[1]+c*(x[0]+x[1]-1)]) #menghitung nabla f(x) dengan x awal
nabla2f=np.array([[2+c,c],[c,2+c]]) #kebetulan konstan, tidak bergantung pada x
print("nabla theta(x0)= ",nablaf)
print("norm(nabla theta(x0))= ",np.linalg.norm(nablaf))
print("theta(x0)= ",x[0]**2+x[1]**2+c/2*((x[0]+x[1]-1)**2))
error=0.00001

iter=1
maxiter=100

while  np.linalg.norm(nablaf)>error and iter<maxiter:
  x=x-np.linalg.inv(nabla2f)@nablaf #memperbaharui x
  print("\n","x",iter,"= ",x)
  nablaf=np.array([2*x[0]+c*(x[0]+x[1]-1)  ,  2*x[1]+c*(x[0]+x[1]-1)]) #menghitung nabla f(x) dengan x awal
  nabla2f=np.array([[2+c,c],[c,2+c]]) #kebetulan konstan, tidak bergantung pada x
  print("nabla theta(x",iter,")= ",nablaf)
  print("norm(nabla theta(x",iter,")= ",np.linalg.norm(nablaf))
  print("theta(x",iter,")=",x[0]**2+x[1]**2+c/2*((x[0]+x[1]-1)**2))
  iter=iter+1

x0=  [0 0]
nabla theta(x0)=  [-1000000 -1000000]
norm(nabla theta(x0))=  1414213.562373095
theta(x0)=  500000.0

 x 1 =  [0.4999995 0.4999995]
nabla theta(x 1 )=  [-3.62280384e-05 -3.62279529e-05]
norm(nabla theta(x 1 )=  5.12341228091845e-05
theta(x 1 )= 0.49999950000050064

 x 2 =  [0.4999995 0.4999995]
nabla theta(x 2 )=  [-4.98774355e-11 -4.98774355e-11]
norm(nabla theta(x 2 )=  7.053734576184381e-11
theta(x 2 )= 0.4999995000005

Untuk $c=1$ diperoleh nilai $\mathbf{x}=\begin{bmatrix}\frac{1}{4}\\\frac{1}{4}\end{bmatrix}$, titik ini belum memenuhi kendala

Untuk $c=1000000$ diperoleh nilai $\mathbf{x}=\begin{bmatrix}\frac{1}{2}\\\frac{1}{2}\end{bmatrix}$, dengan $f\left(\begin{bmatrix}\frac{1}{2}\\\frac{1}{2}\end{bmatrix}\right)=0.5$

Selesaikan dengan Metode Conjugate Gradien

In [3]:
x=np.array([0,0]) #nilai/tebakan awal
print("x0= ",x)
c=1000000 #coba untuk c=1,c=2, c=10, c=100, c=1000000
nablaf=np.array([2*x[0]+c*(x[0]+x[1]-1)  ,  2*x[1]+c*(x[0]+x[1]-1)]) #menghitung nabla f(x) dengan x awal
print("norm(nabla f(x0))= ",np.linalg.norm(nablaf))
nabla2f=np.array([[2+c,c],[c,2+c]])
d=-nablaf #menghitung d

iter=1
maxiter=100
error=0.000001

while  np.linalg.norm(nablaf)>error and iter<maxiter:
  alpha=-d@nablaf/(d@nabla2f@d)
  x=x+alpha*d
  print("x",iter,"= ", x)
  nablaf=np.array([2*x[0]+c*(x[0]+x[1]-1)  ,  2*x[1]+c*(x[0]+x[1]-1)])
  nabla2f=np.array([[2+c,c],[c,2+c]])
  print("norm(nabla f(x",iter,"))= ",np.linalg.norm(nablaf)) 
  d=-nablaf+( nabla2f@nablaf@d/(d@nabla2f@d)  )*d
  iter=iter+1

x0=  [0 0]
norm(nabla f(x0))=  1414213.562373095
x 1 =  [0.4999995 0.4999995]
norm(nabla f(x 1 ))=  8.647205711577957e-11

Selesaikan dengan Algoritma Umum

Perhatikan apabila masalah tersebut diselesaikan secara numerik degan algoritma optimisasi umum

In [4]:
#Algoritma umum

x=np.array([0,0]) #nilai/tebakan awal
alpha=0.1  #coba kombinasi (0,0,1) , (0,0,0.1), (0,0,0.5), (1,1,1) , (1,1,0.1), (1,1,0.5)
c=1 #coba untuk c=1,c=2, c=10, c=100, c=1000000
print("x0= ",x)
nablaf=np.array([2*x[0]+c*(x[0]+x[1]-1)  ,  2*x[1]+c*(x[0]+x[1]-1)]) #menghitung nabla f(x) dengan x awal

print("nabla f(x0)= ",nablaf)
print("norm(nabla f(x0))= ",np.linalg.norm(nablaf))
error=0.00001

iter=1
maxiter=300

while  np.linalg.norm(nablaf)>error and iter<maxiter:
  x=x-alpha*nablaf #memperbaharui x
  nablaf=np.array([2*x[0]+c*(x[0]+x[1]-1)  ,  2*x[1]+c*(x[0]+x[1]-1)]) #menghitung nabla f(x) yang baru dengan x baru
  if(iter%10==0):
    print("\n","x",iter,"= ",x)
    print("nabla f(x",iter,")= ",nablaf)
    print("norm(nabla f(xx",iter,")= ",np.linalg.norm(nablaf))
  iter=iter+1

print("\n","x",iter-1,"= ",x)
print("nabla f(x",iter-1,")= ",nablaf)
print("norm(nabla f(xx",iter-1,")= ",np.linalg.norm(nablaf))

x0=  [0 0]
nabla f(x0)=  [-1 -1]
norm(nabla f(x0))=  1.4142135623730951

 x 10 =  [0.24848835 0.24848835]
nabla f(x 10 )=  [-0.00604662 -0.00604662]
norm(nabla f(xx 10 )=  0.00855120861640389

 x 20 =  [0.24999086 0.24999086]
nabla f(x 20 )=  [-3.65615844e-05 -3.65615844e-05]
norm(nabla f(xx 20 )=  5.170588852129734e-05

 x 24 =  [0.24999882 0.24999882]
nabla f(x 24 )=  [-4.73838134e-06 -4.73838134e-06]
norm(nabla f(xx 24 )=  6.701083152474188e-06

Apa makna yang dapat diambil dari percobaan di atas?

Apabila ada masukkan silahkan sampaikan

Terima kasih

Comments

Post a Comment

Popular posts from this blog

Aljabar Python

-
Python dapat mengoprasikan Aljabar dengan bantuan library numpy. Sedikit operasi vektor dan matriks yang penting diketahui In [1]: import numpy as np import math Vektor ¶ Penting : Bedakan ketika menulis u=[1,2,3] dan u=np.array([1,2,3]) yang pertama adalah list, yang kedua adalah vektor Tampilan vektor pada python selalu tertulis menyamping , walaupun kita tambahkan transpose, tampilan akan tetap menyamping. Hal ini sangat penting dan harus berhati-hati. u=np.array([1,2,3]) artinya $\mathbb{u}=\begin{bmatrix}1\\2\\3\end{bmatrix}$, tampilan pada python u= [1 2 3] np.transpose(u)=$\mathbb{u}^T=\begin{bmatrix}1&2&3\end{bmatrix}$, tampilan pada python u= [1 2 3] In [2]: u = [ 1 , 2 , 3 ] #ini adalah list, bukan vektor v = [ 4 , 5 , 6 ] #ini adalah list, bukan vektor print ( "List u=" , u ) print ( "List v=" , v ) print ( "Penjumlahan list : u+v
Read more

Fungsi Rekursif

-
Fungsi Rekursif ¶ Fungsi yang memanggil dirinya sendiri untuk melakukan perhitungan Faktorial Misalkan $f(n)$ menyatakan faktorial dari $n$, dengan $n$ adalah bilangan bulat lebih dari 1. Maka $f(n)=n\times (n-1)\times (n-2)\times (n-3)\times ...\times 1$ Bisa juga dituliskan dalam bentuk rekursif, yaitu $f(n)=n\times f(n-1)$ Contoh : $f(5)=5\times 4\times 3\times 2\times 1$ $f(5)=5 \times f(4)$ $f(4)=4 \times f(3)$ $f(3)=3 \times f(2)$ $f(2)=2 \times f(1)$ In [1]: #Menghitung faktorial dari n, TANPA fungsi dan tanpa rekursif print ( "Masukkan bilangan bulat lebih dari 1 yang akan dihitung faktorialnya" ) n = int ( input ()) if ( n == 0 or n == 1 ): hasil = 1 else : hasil = n print ( "n saat ini =" , n , "-hasil saat ini" , hasil ) while n >= 2 : n = n - 1 hasil = hasil * n print ( "n saat ini =" , n , "-hasil saat ini" ,
Read more

Pengolahan Matriks (manual)

-
Pengolahan Matriks dan Vektor ¶ Pada tulisan berikut, akan diberikan contoh algoritma dan program untuk mengolah matriks dan vektor tanpa menggunakan library atau fungsi yang sudah tersedia. Transpose Matriks ¶ Misalkan matriks $A=\begin{bmatrix}a_{1,1}&a_{1,2}&\cdots&a_{1,n} \\a_{2,1}&a_{2,2}&\cdots&a_{2,n} \\\vdots&\ddots&\ddots&\vdots \\a_{m,1}&a_{m,2}&\cdots&a_{m,n}\end{bmatrix}$ Transpose dari matriks $A$ adalah $A^T=\begin{bmatrix}a_{1,1}&a_{2,1}&\cdots&a_{m,1} \\a_{1,2}&a_{2,2}&\cdots&a_{m,2} \\\vdots&\ddots&\ddots&\vdots \\a_{1,n}&b_{2,n}&\cdots&a_{n,m}\end{bmatrix}$, atau $B=\begin{bmatrix}b_{1,1}&b_{1,2}&\cdots&b_{1,m} \\b_{2,1}&b_{2,2}&\cdots&b_{2,m} \\\vdots&\ddots&\ddots&\vdots \\b_{n,1}&b_{n,2}&\cdots&b_{n,m}\end{bmatrix}$ , untuk setiap $i\in\{1,2,...,m\}$ dan $j\in\{1,2,..
Read more

Contoh Fungsi : Standar Deviasi

-
Hitung Standar Deviasi dari Barisan Bilangan Bulat ¶ Diberikan sebuah barisan bilangan bulat dengan n elemen, $A = (a_1,a_2,...,a_n)$. Hitunglah standar Deviasinya. Menghitung Rata-Rata : $$\bar{a}=\frac{1}{n}\sum_{i=1}^n a_i$$ Menghitung Standar Deviasi : $$\sigma=\sqrt{\frac{1}{n-1}\sum_{i=1}^n (a_i-\bar{a})^2}$$ Algoritma Hitung Rata-Rata Deskripsi : Hitung rata-rata Masukkan : list/array dengan n elemen, Keluaran : rata-rata Langkah-langkah : Masukkan list/array bilangan, A deklarasi tempat penyimpanan jumlah, sum=0 untuk setiap i pada (1,2,...,panjang A) sum=sum+A[i] rata-rata=sum/panjang A 3. Hitung Standar Deviasi Deskripsi : Hitung Standar Deviasi Masukkan : list/array bilangan dengan n elemen, Keluaran : nilai standar deviasinya Langkah-langkah : Masukkan n (banyak elemen list/array) Buat list/array berukuran n, A Hitung rata-rata dari A, simpan pada variabel rata Masukkan list/array bilanga
Read more