makispaiktis

Texnikes_Ergasia2_Thema2_ZitoumenoC V2

Dec 1st, 2020
561
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
  1. clear all
  2. clc
  3.  
  4. elaxistoMeGammaMetavlhto(-1, -1, 0.8, 0);
  5. elaxistoMeGammaMetavlhto(-1, -1, 1, 1);
  6.  
  7.  
  8. function elaxistoMeGammaMetavlhto(x0, y0, s, flag)
  9.    
  10.     % 1. Ορίζω το ε της συνθήκης τερματισμού ίσο με 2/1000
  11.     epsilon = 0.002;
  12.     syms x y
  13.     f = x^3 * exp(-x^2-y^4);
  14.     klisi = gradient(f, [x,y]);
  15.    
  16.     % 2. Ορίζω τις λίστες που θα τοποθετήσω τα xi, yi. Τις ονομάζω xList,
  17.     % yList, θα βάζω επίσης και τις τιμές της f (fList) και του μέτρου της
  18.     % κλίσης της (normKlisisList) σε άλλες 2 λίστες
  19.     k = 1;
  20.     xList = [];             yList = [];
  21.     xList(1) = x0;          yList(1) = y0;
  22.     fList = [];             normKlisisList = [];        gammaList = [];
  23.     fList(1) = subs(f, {x,y}, {xList(length(xList)), yList(length(yList))});
  24.     normKlisisList(1) = norm(subs(klisi, {x,y}, {xList(length(xList)), yList(length(yList))}));
  25.    
  26.     % Armijo Constants
  27.     a = 1/1000;             b = 0.3;        
  28.     gammaList(1) = s;
  29.    
  30.     while normKlisisList(length(normKlisisList)) > epsilon
  31.         k = k + 1;
  32.         dk = -subs(klisi, {x,y}, {xList(length(xList)), yList(length(yList))});     % η κλίση πριν
  33.         % Το dk είναι 2*1 διάνυσμα, το 1ο στοιχείο του αφορά τον υπολογισμό
  34.         % του x και το 2ο τον υπολογισμό του y
  35.         % Πρέπει το γ που θα επιλέξω να κάνει minimize την f(xk + γkdk)
  36.         xPrin = xList(k-1);
  37.         yPrin = yList(k-1);
  38.        
  39.         % ********************************************************************************************
  40.         % ******************************* Κανόνας Armijo *********************************************
  41.         % ********************************************************************************************
  42.         % Δοκιμάζω τιμές ακέραιες 0,1,2,... στο mk και φτιάχνω το βήμα γk
  43.         % και κατ' επέκταση και το xk+1
  44.         for mk = 0:100
  45.             k
  46.             mk
  47.             gk = s * b^mk;
  48.             xNew = xPrin + gk * dk(1);
  49.             yNew = yPrin + gk * dk(2);
  50.             F = matlabFunction(f);
  51.             DF = matlabFunction(klisi);
  52.             % Ήρθε η ώρα να τσεκάρουμε την ανίσωση
  53.             aristeroMelos = F(xPrin, yPrin) - F(xNew, yNew);
  54.             grad = DF(xPrin, yPrin);
  55.             dexioMelos = a * b^mk * s * grad' * grad;
  56.             if aristeroMelos > dexioMelos
  57.                 gammaList(k) = gk;
  58.                 xList(k) = xNew;
  59.                 yList(k) = yNew;
  60.                 fList(k) = F(xNew, yNew);
  61.                 normKlisisList(k) = norm(subs(klisi, {x,y}, {xList(length(xList)), yList(length(yList))}));
  62.                 display('Success for this step k with mk in value of:')
  63.                 mk
  64.                 break;
  65.             end
  66.         end
  67.            
  68.            
  69.     end
  70.    
  71.     xList
  72.     yList
  73.     fList
  74.     normKlisisList
  75.     k
  76.     gammaList
  77.     % Τα τελευταία xk, yk κάθε λίστας τα ονομάζω εν συντομία xx και yy
  78.     display('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~')
  79.     display('~~~~~~~~ About the last found xk (xx) and yk (yy) ~~~~~~~~')
  80.     display('~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~')
  81.     xx = xList(length(xList))
  82.     yy = yList(length(yList))
  83.     F_xx_yy = fList(length(fList))
  84.     NORM_KLISIS = normKlisisList(length(normKlisisList))
  85.     if flag == 0
  86.         plot(fList, 'ro');
  87.         title('Red for s = 0.8, blue for s = 1');
  88.         xlabel('Steps k');
  89.         ylabel('Values of f');
  90.         hold on
  91.     else
  92.         plot(fList, 'bx');
  93.         title('Red for s = 0.8, blue for s = 1');
  94.         xlabel('Steps k');
  95.         ylabel('Values of f');
  96.         hold on
  97.     end
  98.     display('**********************************************************')
  99. end
  100.    
  101.  
RAW Paste Data