%
% !TEX TS-program = xelatex
% !TEX encoding   = UTF-8 Unicode
% !TEX spellcheck = el-GR
%
% System Modeling Lab02 assignment report
%
% Requires compilation with XeLaTeX
%
% authors:
%   Χρήστος Χουτουρίδης ΑΕΜ 8997
%   cchoutou@ece.auth.gr




% Options:
%
% 1) mainlang=<language>
%    Default:  english
%    Set the default language of the document which affects hyphenations,
%    localization (section, dates, etc...)
%
%    example: \documentclass[mainlang=greek]{AUThReport}
%
% 2) <language>
%    Add hyphenation and typesetting support for other languages
%    Currently supports: english, greek, german, frenc
%
%    example: \documentclass[english, greek]{AUThReport}
%
% 3) short: Requests a shorter title for the document
%    Default: no short
%
%    example: \documentclass[short]{AUThReport}
%
\documentclass[a4paper, 11pt, mainlang=greek, english]{AUThReport/AUThReport}

\CurrentDate{\today}

% Greek report document setup suggestions
%---------------------------------
% Document configuration
\AuthorName{Χρήστος Χουτουρίδης}
\AuthorAEM{8997}
\AuthorMail{cchoutou@ece.auth.gr}

%\CoAuthorName{CoAuthor Name}
%\CoAuthorAEM{AEM}
%\CoAuthorMail{CoAuthor Mail}

% \WorkGroup{Ομάδα Χ}

\DocTitle{Εργασία 2}
\DocSubTitle{Μέθοδοι Πραγματικού Χρόνου, Μέθοδος Κλίσης, Μέθοδος Lyapunov}

\Department{Τμήμα ΗΜΜΥ. Τομέας Ηλεκτρονικής}
\ClassName{Προσομοίωση και Μοντελοποίηση Δυναμικών Συστημάτων}

\InstructorName{Γ. Ροβιθάκης}
\InstructorMail{rovithak@auth.gr}

\CoInstructorName{Λ. Μπίκας}
\CoInstructorMail{lnmpikas@ece.auth.gr}


% Local package requirements
%---------------------------------
%\usepackage{tabularx}
%\usepackage{array}
%\usepackage{commath}

\usepackage{amsmath, amssymb, amsfonts}
\usepackage{graphicx}
\usepackage{subcaption}
\usepackage{float}

\begin{document}

% Request a title page or header
\InsertTitle

\section{Εισαγωγή}

Η παρούσα εργασία ασχολείται με την εκτίμηση αγνώστων παραμέτρων σε δυναμικά συστήματα μέσω μεθόδων πραγματικού χρόνου, δίνοντας έμφαση στη μέθοδο Lyapunov.
Στο πρώτο σκέλος, μελετάται ένα γραμμικό σύστημα μάζας-ελατηρίου-απόσβεσης και αναπτύσσονται δύο εκτιμητές — βασισμένοι στις μεθόδους gradient descent και Lyapunov — με στόχο την ταχεία και αξιόπιστη προσέγγιση των φυσικών παραμέτρων του συστήματος.
Εξετάζεται επίσης η επίδραση θορύβου στη μέτρηση της θέσης.
\par
Στο δεύτερο σκέλος, μελετάται ένα μη γραμμικό μοντέλο ρόλισης αεροσκάφους, όπου οι άγνωστες παράμετροι εμπλέκονται μη γραμμικά στις εξισώσεις κίνησης.
Η εκτίμηση συνδυάζεται με ελεγκτή παρακολούθησης τροχιάς, ο οποίος βασίζεται σε κανονικοποιημένα σφάλματα και λογαριθμική κορεσμένη ανάδραση.
Παρουσιάζεται η εκτίμηση των παραμέτρων, καθώς και η ανακατασκευή της εξόδου με βάση τις εκτιμώμενες τιμές.
Τέλος, μελετάται η ευαισθησία του εκτιμητή παρουσία εξωτερικών διαταραχών.

\subsection{Παραδοτέα}
Τα παραδοτέα της εργασίας αποτελούνται από:
\begin{itemize}
    \item Την παρούσα αναφορά.
    \item Τον κατάλογο \textbf{scripts/}, που περιέχει τον κώδικα της MATLAB.
    \item Το \href{https://git.hoo2.net/hoo2/SystemModling/src/branch/master/Lab02}{σύνδεσμο} με το αποθετήριο που περιέχει όλο το project με τον κώδικα της MATLAB, της αναφοράς και τα παραδοτέα.
\end{itemize}

\section{Θέμα 1}

\subsection{1α - Gradient Estimation}

Το σύστημα που μελετάται είναι το κλασικό σύστημα μάζας-ελατηρίου-απόσβεσης:
\[
m\ddot{x}(t) + b\dot{x}(t) + kx(t) = u(t)
\]
το οποίο αναδιατυπώνεται γραμμικά ως προς τις άγνωστες παραμέτρους:
\[
y(t) = u(t) = \theta^{\top} u(t), \quad \text{με } \theta = [m,\, b,\, k]^{\top},\ u(t) = [\ddot{x}(t),\, \dot{x}(t),\, x(t)]^{\top}
\]
Ο εκτιμητής gradient έχει τη μορφή:
\[
\dot{\hat{\theta}}(t) = \gamma \cdot u(t) \cdot \left(y(t) - \hat{\theta}^{\top}(t) u(t)\right)
\]

Η προσομοίωση του ``πραγματικού'' συστήματος υλοποιήθηκε μέσω της μεθόδου \textbf{Runge–Kutta} 4ης τάξης (RK4), ενώ οι παράγωγοι υπολογίστηκαν με διαφορικό σχήμα.
Το script επιτρέπει επιλογή μεταξύ \textbf{normalized} και \textbf{unnormalized} εκτιμητή και διαφορετικές τιμές $\gamma$.

\subsubsection{Σύγκριση Normalized vs Unnormalized}

Παρακάτω παρουσιάζονται τα αποτελέσματα για \boldmath $\gamma = 0.66$ και $T=30$ \unboldmath \textbf{sec}, για τις δύο περιπτώσεις.
Οι εκτιμήσεις είναι διαφορετικές, κυρίως στις παραμέτρους $m$ και $b$, γεγονός που υποδηλώνει πως η χρήση normalization επηρεάζει όχι μόνο τη δυναμική, αλλά και το σημείο σύγκλισης.

\begin{figure}[H]
    \centering
    \begin{subfigure}[t]{0.48\textwidth}
        \includegraphics[width=\linewidth]{../scripts/output/Prob1a_estimation_constant_gamma0.660_Unnormalized_30s.png}
        \caption{Constant input, unnormalized}
    \end{subfigure}
    \hfill
    \begin{subfigure}[t]{0.48\textwidth}
        \includegraphics[width=\linewidth]{../scripts/output/Prob1a_estimation_constant_gamma0.660_Normalized_30s.png}
        \caption{Constant input, normalized}
    \end{subfigure}

    \vspace{0.5cm}

    \begin{subfigure}[t]{0.48\textwidth}
        \includegraphics[width=\linewidth]{../scripts/output/Prob1a_estimation_sine_gamma0.660_Unnormalized_30s.png}
        \caption{Sine input, unnormalized}
    \end{subfigure}
    \hfill
    \begin{subfigure}[t]{0.48\textwidth}
        \includegraphics[width=\linewidth]{../scripts/output/Prob1a_estimation_sine_gamma0.660_Normalized_30s.png}
        \caption{Sine input, normalized}
    \end{subfigure}
    \caption{Σύγκριση εκτίμησης παραμέτρων με και χωρίς normalization, για $\gamma = 0.66$ και $T=30$ sec.}
    \label{gradient_1a}
\end{figure}

Από το σχήμα \ref{gradient_1a} παρατηρούμε ότι το σύστημα απαιτεί σημαντικό χρονικό διάστημα για να συγκλίνει σε σταθερές τιμές.
Ειδικά για μικρές τιμές του $\gamma$, η σύγκλιση είναι πιο αργή αλλά ομαλότερη, η τιμή του οποίου επιλέχθηκε με δοκιμαστικές εκτελέσεις.
Παρόλα αυτά το σύστημα παρουσίαζε \textbf{μεγάλες ταλαντώσεις}.
Για το λόγο αυτό επιλέχθηκε να γίνει κανονικοποίηση στη μέθοδο κλίσης (normalized gradient update).
Η χρήση κανονικοποίησης βελτιώνει σημαντικά τη σταθερότητα του εκτιμητή και οδηγεί σε εκτιμήσεις με μικρότερη διακύμανση, ιδιαίτερα για τις παραμέτρους $m$ και $b$.
Ωστόσο, ακόμη και με κανονικοποίηση, το σφάλμα εκτίμησης παραμένει, γεγονός που πιθανόν οφείλεται είτε στην αριθμητική προσέγγιση των παραγώγων είτε στο ότι η είσοδος δεν διεγείρει επαρκώς όλες τις δυναμικές του συστήματος.
\par
Η πρώτη εκτέλεση (για $T=30s$ και $\gamma=0.66$) ανέδειξε αυτές τις διαφορές: χωρίς normalization παρατηρήθηκαν εντονότερες αποκλίσεις, ενώ με normalization η σύγκλιση ήταν πιο σταθερή, ειδικά για την παράμετρο $b$.

\subsubsection{Διερεύνηση επίδρασης της παραμέτρου $\gamma$}

Προκειμένου να κατανοήσουμε καλύτερα την επίδραση του $\gamma$, εκτελέσαμε το script με τιμές $\gamma = 0.33,\ 0.50,\ 0.66$ και normalized update.
Ο πίνακας που ακολουθεί παρουσιάζει τις τελικές τιμές εκτίμησης για κάθε περίπτωση:

\begin{table}[H]
    \centering
    \caption{Τελικές εκτιμήσεις παραμέτρων για διαφορετικές τιμές $\gamma$ (με normalized update)}
    \label{tab:gamma_comparison}
    \begin{tabular}{r c c c c}
        \textbf{Σήμα} & \textbf{$\gamma$}  & $\hat{m}$ & $\hat{b}$ & $\hat{k}$ \\
        \hline
        σταθερό      & 0.33 & 0.9963 & 0.6650          & \textbf{0.7143} \\
        σταθερό      & 0.50 & 0.9726 & 0.5558          & \textbf{0.7119} \\
        σταθερό      & 0.66 & 0.9475 & 0.4795          & \textbf{0.7107} \\
        ημιτονοειδές & 0.33 & 1.1492 & \textbf{0.2452} & 0.5622 \\
        ημιτονοειδές & 0.50 & 1.1033 & \textbf{0.2326} & 0.5097 \\
        ημιτονοειδές & 0.66 & 1.0672 & \textbf{0.2350} & 0.4735 \\
        \hline
        Πραγματικές  & --   & 1.315  & 0.225  & 0.725  \\
    \end{tabular}
\end{table}

\subsubsection{Παρατήρηση επιλεκτικής εκτίμησης ανάλογα με το είδος εισόδου}

Από τη σύγκριση των αποτελεσμάτων προκύπτει ότι η εκτίμηση της $k$ προσεγγίζεται επιτυχώς κυρίως όταν το σύστημα διεγείρεται με σταθερό σήμα.
Αυτό είναι λογικό, καθώς σε σταθερή είσοδο κυριαρχεί η στατική απόκριση του συστήματος, στην οποία η $k$ (σταθερά ελατηρίου) έχει άμεση επίδραση.
\par
Αντίθετα, σε ημιτονοειδές σήμα, οι παραγώγοι της ταχύτητας και της επιτάχυνσης ποικίλουν έντονα και επαναληπτικά, ενισχύοντας τη συμβολή των $m$ και $b$ στη δυναμική του σήματος.
Παρατηρούμε ότι σε αυτές τις περιπτώσεις η εκτίμηση της $b$ είναι σαφώς βελτιωμένη.
\par
Η διαφοροποίηση αυτή υποδεικνύει ότι η επιμονή διέγερσης (persistence of excitation) που απαιτείται για την εκτίμηση όλων των παραμέτρων δεν ικανοποιείται εξίσου από κάθε τύπο εισόδου.
Για να επιτευχθεί συνολική παρατηρησιμότητα, ενδεχομένως να απαιτείται πιο σύνθετο ή ακόμα και στοχαστικό σήμα εισόδου.


\subsection{1β – Εκτίμηση παραμέτρων με μέθοδο Lyapunov}

Η μέθοδος Lyapunov προσφέρει μία πιο θεωρητικά τεκμηριωμένη προσέγγιση για την εκτίμηση αγνώστων παραμέτρων, καθώς βασίζεται στον σχεδιασμό κατάλληλης συνάρτησης Lyapunov που εγγυάται τη σύγκλιση των σφαλμάτων εκτίμησης.
Η εκτίμηση των παραμέτρων βασίστηκε στην ελαχιστοποίηση του τετραγωνικού σφάλματος μεταξύ της πραγματικής εισόδου $u(t)$ και της εκτιμώμενης εξόδου του μοντέλου:
\[
\hat{y}(t) = \hat{\theta}^\top \phi(t)
\]
όπου:
\[
y(t) = u(t), \quad \phi(t) = [\ddot{x}(t),\ \dot{x}(t),\ x(t)]^\top
\]

Η συνάρτηση κόστους ορίζεται ως:
\[
K(\hat{\theta}) = \frac{1}{2} (y(t) - \hat{y}(t))^2 = \frac{1}{2} e^2(t)
\]
και η προσαρμογή των παραμέτρων γίνεται με βάση τον κανόνα gradient descent:
\[
\dot{\hat{\theta}} = -\gamma \cdot \nabla K(\hat{\theta}) = \gamma \cdot e(t) \cdot \phi(t)
\]

\boldmath
Αυτή η μορφή δεν απαιτεί ξεχωριστό βοηθητικό δυναμικό μοντέλο ή υπολογισμό σφαλμάτων κατάστασης, αλλά \textbf{στηρίζεται αποκλειστικά σε παρατηρήσιμα σήματα} και προβλέψιμη γραμμική εξίσωση ως προς τις παραμέτρους.
Η μέθοδος αυτή παρουσιάζει καλύτερη σταθερότητα σε σχέση με τη μέθοδο gradient, ειδικά όταν η είσοδος δεν πληροί πλήρως τη συνθήκη επιμονής διέγερσης.
Επιπλέον, η εισαγωγή της παραμέτρου $\gamma$ επιτρέπει τον έλεγχο της ταχύτητας προσαρμογής, χωρίς να απαιτείται normalization.
\unboldmath
\subsubsection{Υλοποίηση Εκτιμητή}

\boldmath
Ο εκτιμητής Lyapunov υλοποιήθηκε με βάση το παραπάνω θεωρητικό μοντέλο.
Επιπλέον, ανακατασκευάστηκε η έξοδος του συστήματος $x(t)$ με βάση τις εκτιμώμενες παραμέτρους, ώστε να συγκριθεί η πραγματική έξοδος με την εκτιμώμενη $\hat{x}(t)$ και να υπολογιστεί το σφάλμα $e_x(t) = x(t) - \hat{x}(t)$.
\unboldmath
\begin{figure}[H]
    \centering
    \begin{subfigure}[t]{0.48\textwidth}
        \includegraphics[width=\linewidth]{../scripts/output/Prob1b_lyapunov_gamma0.660_40s.png}
        \caption{Εκτίμηση παραμέτρων με τη μέθοδο Lyapunov}
    \end{subfigure}
    \hfill
    \begin{subfigure}[t]{0.48\textwidth}
        \includegraphics[width=\linewidth]{../scripts/output/Prob1b_extrastates_gamma0.660_40s.png}
        \caption{Αντίδραση συστήματος και σφάλμα εκτίμησης θέσης}
    \end{subfigure}
    \caption{Συγκριτική παρουσίαση εκτίμησης παραμέτρων και εξόδου (για $\gamma=0.66$, $T=40s$)}
    \label{lyapunov_1b}
\end{figure}

Όπως φαίνεται και από το σχήμα \ref{lyapunov_1b} η μέθοδος Lyapunov οδήγησε σε σταθερή και σχετικά ακριβή εκτίμηση των παραμέτρων του συστήματος.
Η τελική εκτίμηση των παραμέτρων μετά από $40s$ ήταν:
\[
\hat{m} = 0.8428,\quad \hat{b} = 0.2201,\quad \hat{k} = 0.2620
\]
Παρατηρείται ότι οι εκτίμηση για το $b$ συγκλίνει σε τιμές κοντά στις πραγματικές, ενώ η εκτιμήσεις του $m$ και $k$ αποκλίνουν σημαντικά, ομοίως με τη μέθοδο gradient estimation.
\par
Η ανακατασκευή της εξόδου $\hat{x}(t)$ παρόλα αυτά, με βάση το εκτιμώμενο μοντέλο \textbf{παρουσίασε πολύ καλή συμφωνία με την πραγματική} $x(t)$, με το σφάλμα θέσης $e_x(t)$ να παραμένει κάτω από $0.6$ καθ' όλη τη διάρκεια της προσομοίωσης.
Αυτό ενισχύει την αξιοπιστία της εκτίμησης και αποδεικνύει την αποτελεσματικότητα της μεθόδου σε σταθερές συνθήκες.

\subsection{1γ – Επίδραση διαταραχής στη μέτρηση x(t)}

Για τη μελέτη της ευαισθησίας της μεθόδου Lyapunov στην παρουσία θορύβου, προστέθηκε στο σήμα της θέσης $x(t)$ μια ημιτονοειδής διαταραχή της μορφής:
\[
\eta(t) = \eta_0 \cdot \sin(2\pi f_0 t)
\]
με πλάτος $\eta_0 = 0.1$ και συχνότητα $f_0 = 0.5\,\text{Hz}$.

Ο εκτιμητής παραμένει ίδιος με αυτόν του Θέματος 1β.
Η εκτίμηση έγινε δύο φορές: μία με καθαρό σήμα $x(t)$ και μία με διαταραχή $x(t) + \eta(t)$, ενώ οι παράγωγοι $\dot{x}(t)$ και $\ddot{x}(t)$ υπολογίστηκαν από το καθαρό σήμα ώστε να αποφευχθεί υπερβολική αριθμητική αστάθεια.

\begin{figure}[H]
    \centering
    \includegraphics[width=0.9\textwidth]{../scripts/output/Prob1c_disturbance_eta0.10.png}
    \caption{Σύγκριση εκτιμήσεων με και χωρίς διαταραχή ($\eta_0 = 0.1$)}
    \label{fig:disturbance}
\end{figure}

Από το σχήμα \ref{fig:disturbance} παρατηρούμε ότι η προσθήκη περιορισμένης διαταραχής στη μέτρηση της θέσης επηρεάζει ελάχιστα την εκτίμηση.
Οι τελικές τιμές των παραμέτρων \textbf{συγκλίνουν σχεδόν ταυτόσημα με την καθαρή περίπτωση}:
\[
\hat{m} = 0.8432, \quad \hat{b} = 0.2195, \quad \hat{k} = 0.2590
\]

Η μικρή επίδραση οφείλεται στο γεγονός ότι η διαταραχή:
\begin{itemize}
    \item είναι ομαλή και περιοδική (χαμηλή συχνότητα),
    \item φιλτράρεται έμμεσα μέσα από το μοντέλο,
    \item δεν επηρεάζει τις παραγώγους που χρησιμοποιούνται στον εκτιμητή.
\end{itemize}

Συνεπώς, ο εκτιμητής Lyapunov παρουσιάζει \textbf{ικανοποιητική ανοσία} σε ήπια διαταραχή, γεγονός που ενισχύει τη χρησιμότητά του σε ρεαλιστικά σενάρια μέτρησης.

\section{Θέμα 2}

Το σύστημα που εξετάζεται είναι ένα δυναμικό μη γραμμικό σύστημα δεύτερης τάξης που περιγράφει τη ροπή κύλισης (roll motion) ενός αεροσκάφους:
\[
\ddot{r}(t) = -a_1 \dot{r}(t) - a_2 \sin(r(t)) + a_3 \dot{r}^2 \sin(2r(t)) + b u(t) + d(t)
\]
Όπου:
\begin{itemize}
    \item $r(t)$ είναι η γωνία roll,
    \item $u(t)$ είναι η είσοδος (σήμα ελέγχου),
    \item $d(t)$ είναι εξωτερική διαταραχή,
    \item $a_1, a_2, a_3, b$ είναι άγνωστες σταθερές παραμέτρων.
\end{itemize}

Η δομή του συστήματος επιβάλλει την κατασκευή εκτιμητή παραμέτρων \textbf{σε πραγματικό χρόνο}.
Το βασικό βήμα είναι η αναγραφή της εξίσωσης σε μορφή γραμμική ως προς τις άγνωστες παραμέτρους:
\[
\ddot{r}(t) = \theta^\top \phi(t)
\quad \text{με } \theta = [a_1,\, a_2,\, a_3,\, b]^\top,\quad
\phi(t) = [-\dot{r}(t),\, -\sin(r(t)),\, \dot{r}^2(t)\sin(2r(t)),\, u(t)]^\top
\]
Η είσοδος $u(t)$ πρέπει να σχεδιαστεί κατάλληλα ώστε το σύστημα να ακολουθεί τη ζητούμενη τροχιά ακόμη και παρουσία άγνωστων παραμέτρων.
Για τον σκοπό αυτό, προτείνεται η χρήση ενός ελεγκτή ανάδρασης που βασίζεται σε κανονικοποιημένα σφάλματα και σε μια smooth συνάρτηση κορεσμού της μορφής:
\[
T(z) = \ln\left(\frac{1 + z}{1 - z}\right)
\]
Η συνάρτηση $T(z)$ χρησιμοποιείται για να περιορίσει την έξοδο του ελεγκτή όταν τα σφάλματα γίνονται μεγάλα, αποτρέποντας απότομες ή μη ρεαλιστικές εντολές ελέγχου.
Τα σφάλματα κανονικοποιούνται με χρονικά μεταβαλλόμενα ή σταθερά φράγματα, ώστε να εξασφαλίζεται η σταθερότητα και η συμβατότητα με φυσικούς περιορισμούς του συστήματος.

Ο ελεγκτής βασίζεται σε δύο επίπεδα ανάδρασης.
Πρώτα, υπολογίζεται ένα \textbf{setpoint ταχύτητας}:
\[
z_1(t) = \frac{r(t) - r_d(t)}{\phi(t)}, \quad \alpha(t) = -k_1 T(z_1(t))
\]
όπου $r_d(t)$ είναι η επιθυμητή τροχιά και $\phi(t)$ ένα χρονικά φθίνον φράγμα σφάλματος.

Στη συνέχεια, υπολογίζεται η \textbf{τελική είσοδος ελέγχου}:
\[
z_2(t) = \frac{\dot{r}(t) - \alpha(t)}{\rho}, \quad u(t) = -k_2 T(z_2(t))
\]
με σταθερό φράγμα $\rho$ και ενισχυτικό κέρδος $k_2$.
Το σχήμα αυτό επιτρέπει \textbf{υψηλή ακρίβεια παρακολούθησης} (καθώς $\phi_\infty \downarrow$) και ταυτόχρονα περιορίζει τις μεταβατικές ταλαντώσεις μέσω της δομής κορεσμού.

Η μελέτη επεκτείνεται και σε περιπτώσεις με εξωτερικές διαταραχές $d(t)$, όπου αξιολογείται η ευρωστία του εκτιμητή και του ελεγκτή.

\subsection{2α – Εκτίμηση παραμέτρων σε μη γραμμικό σύστημα roll χωρίς διαταραχές}

Εξετάζουμε ένα μη γραμμικό σύστημα δεύτερης τάξης που μοντελοποιεί τη ροπή κύλισης (roll) αεροσκάφους, του οποίου η τροχιά αναφοράς ακολουθεί την προδιαγραφή:
\[
r_d(t) =
\begin{cases}
    0, & 0 \leq t < 10 \\
    \frac{\pi}{10}, & 10 \leq t < 20 \\
    0, & t \geq 20
\end{cases}
\]
Για την παρακολούθηση τροχιάς χρησιμοποιείται ελεγκτής ανάδρασης που περιγράφηκε παραπάνω

\begin{figure}[H]
    \centering
    \includegraphics[width=0.9\textwidth]{../scripts/output/Problem2a_estimation.png}
    \caption{Εκτίμηση παραμέτρων και παρακολούθηση τροχιάς για το Θέμα 2α (χωρίς διαταραχή)}
    \label{fig:2a}
\end{figure}

Από το σχήμα \ref{fig:2a} φαίνεται ότι η \textbf{σύγκλιση} των παραμέτρων είναι \textbf{ομαλή και σχετικά ακριβής}:
\[
\hat{a}_1 = 1.4970,\quad \hat{a}_2 = 1.0015,\quad \hat{a}_3 = 0.7439,\quad \hat{b} = 2.0051
\]
\boldmath
με τις πραγματικές τιμές να είναι $a_1 = 2.0$, $a_2 = 1.0$, $a_3 = 0.5$, $b = 2.0$.
Παρατηρούμε ότι οι παράμετροι $a_2$ και $b$ εκτιμώνται με υψηλή ακρίβεια, ενώ οι $a_1$ και $a_3$ εκτιμούνται με χαμηλή ακρίβεια.
\par
Η παρακολούθηση της επιθυμητής τροχιάς παρόλα αυτά είναι επιτυχής, με το $r(t)$ να ακολουθεί το $r_d(t)$ εντός των ορίων που ορίζει το φράγμα $\phi(t)$ καθ’ όλη τη διάρκεια της προσομοίωσης.
\unboldmath

\subsection{2β – Εκτίμηση παραμέτρων και ανακατασκευή κατάστασης με μετρήσιμα r(t), ṙ(t), u(t)}

Σε αυτό το Θέμα εξετάζεται η περίπτωση κατά την οποία είναι μετρήσιμα τα σήματα $r(t)$, $\dot{r}(t)$ και $u(t)$, και ζητείται η εκτίμηση των αγνώστων παραμέτρων του συστήματος καθώς και η ανακατασκευή της γωνίας roll $r(t)$ μέσω της εκτιμημένης δυναμικής.

Το σύστημα που προσομοιώνεται παραμένει:
\[
\ddot{r}(t) = -a_1 \dot{r}(t) - a_2 \sin(r(t)) + a_3 \dot{r}^2 \sin(2r(t)) + b u(t)
\]

Οι παράμετροι εκτιμώνται χρησιμοποιώντας τη μορφή:
\[
\hat{\theta}(k+1) = \hat{\theta}(k) + \gamma \cdot e(t) \cdot \phi(t)
\quad \text{όπου } \phi(t) = [-\dot{r}, -\sin(r), \dot{r}^2\sin(2r), u(t)]^\top
\]

Η έξοδος του συστήματος $r(t)$ ανακατασκευάζεται από την εκτιμώμενη επιτάχυνση:
\[
\hat{\ddot{r}}(t) = \hat{\theta}^\top \phi(\hat{r}(t), \hat{\dot{r}}(t), u(t)).
\]

\begin{figure}[H]
    \centering
    \includegraphics[width=0.9\textwidth]{../scripts/output/Problem2b_estimation.png}
    \caption{Εκτίμηση παραμέτρων και ανακατασκευή κατάστασης για το Θέμα 2β}
    \label{fig:2b}
\end{figure}

Η εκτίμηση παραμέτρων \textbf{συγκλίνει σε πολύ καλές τιμές}:
\[
\hat{a}_1 = 1.4970,\quad \hat{a}_2 = 1.0015,\quad \hat{a}_3 = 0.7439,\quad \hat{b} = 2.0051
\]
Όπως φαίνεται στο σχήμα \ref{fig:2b}, ενώ η ανακατασκευή της κατάστασης $r(t)$ είναι εύστοχη, με σφάλμα $e_r(t)$ να διατηρείται κάτω από $0.05\,\text{rad}$.
Η ελαφριά χρονική υστέρηση και υποεκτίμηση σε ορισμένα διαστήματα οφείλονται στο ότι η ανακατασκευή γίνεται με βάση εκτιμώμενες παραμέτρους και αρχικές μηδενικές συνθήκες.
Συνολικά, η μέθοδος Lyapunov αποδείχθηκε \textbf{κατάλληλη για ταυτόχρονη εκτίμηση και παρακολούθηση σε πραγματικό χρόνο} με αξιόπιστα αποτελέσματα.

\subsection{2γ – Εκτίμηση παραμέτρων παρουσία εξωτερικής διαταραχής}

Σε αυτό το ερώτημα, επαναλαμβάνουμε την προηγούμενη διαδικασία, αλλά αυτή τη φορά λαμβάνεται υπόψη η παρουσία εξωτερικής διαταραχής:
\[
d(t) = 0.15 \sin(0.5 t)
\]
Η διαταραχή προστίθεται στο δυναμικό μοντέλο του συστήματος, αλλά θεωρείται άγνωστη από τον εκτιμητή.
Έτσι, η μέθοδος Lyapunov επιχειρεί την εκτίμηση των παραμέτρων $a_1$, $a_2$, $a_3$, $b$ \textbf{χωρίς να λαμβάνει υπόψη} τη συνεισφορά του $d(t)$:
\[
\ddot{r}(t) = \theta^\top \phi(t) + d(t)
\]

Η ανακατασκευή της κατάστασης γίνεται όπως και στο προηγούμενο ερώτημα:
\[
\hat{\ddot{r}}(t) = \hat{\theta}^\top \phi(\hat{r}, \hat{\dot{r}}, u)
\]

\begin{figure}[H]
    \centering
    \includegraphics[width=0.9\textwidth]{../scripts/output/Problem2c_estimation.png}
    \caption{Εκτίμηση παραμέτρων και ανακατασκευή $r(t)$ με εξωτερική διαταραχή}
    \label{fig:2c}
\end{figure}

Όπως παρατηρείται στο σχήμα \ref{fig:2c}, η εισαγωγή διαταραχής οδηγεί σε περιοδικό σφάλμα στη σύγκλιση της εκτιμώμενης κατάστασης $\hat{r}(t)$.
Η ταλάντωση στο $e_r(t)$ συνδέεται άμεσα με τη συχνότητα και πλάτος της $d(t)$.

Ωστόσο, οι τελικές εκτιμήσεις παραμέτρων \textbf{παραμένουν σχεδόν αναλλοίωτες}:
\[
\hat{a}_1 = 1.4996,\quad \hat{a}_2 = 0.9356,\quad \hat{a}_3 = 0.7244,\quad \hat{b} = 1.9556
\]
γεγονός που επιβεβαιώνει ότι η μέθοδος Lyapunov είναι σχετικά ανθεκτική σε διαταραχές χαμηλής ενέργειας ή συστηματικά περιοδικού χαρακτήρα.

Το σφάλμα σε $a_2$ και $a_3$ αυξάνεται ελαφρώς σε σχέση με την περίπτωση χωρίς διαταραχή, καθώς η δομή του σήματος $d(t)$ εμφανίζεται στο σφάλμα.
Παρ’ όλα αυτά, το συνολικό σύστημα εκτίμησης παραμένει σταθερό και προβλέψιμο.

\section{Συμπεράσματα}

Στην παρούσα εργασία μελετήθηκαν και συγκρίθηκαν διάφορες μέθοδοι εκτίμησης παραμέτρων δυναμικών συστημάτων.
Στο αρχικό στάδιο της εργασίας εφαρμόστηκε η μέθοδος \textbf{gradient descent} για την εκτίμηση των παραμέτρων του γραμμικού συστήματος.
Αν και η μέθοδος είναι απλή στον υπολογισμό και απαιτεί ελάχιστους πόρους, παρατηρήθηκε ότι η απόδοσή της εξαρτάται έντονα από τη μορφή του σήματος εισόδου και την επιλογή της σταθεράς εκμάθησης $\gamma$.
Πιο συγκεκριμένα, \textbf{η μέθοδος απέτυχε να εκτιμήσει ταυτόχρονα και με ακρίβεια όλες τις παραμέτρους του συστήματος}.

Στη συνέχεια μελετήθηκε η μέθοδος \textbf{Lyapunov} η οποία έδειξε ότι παρουσιάζει καλύτερη σταθερότητα και λιγότερη εξάρτηση από την επιλογή εισόδου και αρχικών συνθηκών.
Στο πλαίσιο γραμμικού συστήματος, η μέθοδος απέδειξε καλή ακρίβεια και σταθερότητα, ακόμη και παρουσία ήπιου θορύβου.
Στην περίπτωση του μη γραμμικού συστήματος roll, οι εκτιμητές πραγματικού χρόνου ενσωματώθηκαν με ελεγκτές παρακολούθησης τροχιάς και αποδείχθηκαν ικανοί να παρακολουθήσουν τις επιθυμητές κινήσεις και να προσαρμοστούν σταδιακά στις άγνωστες παραμέτρους του μοντέλου.
Παράλληλα, αξιολογήθηκε η ευρωστία των εκτιμητών υπό την επίδραση περιοδικής εξωτερικής διαταραχής, με τα αποτελέσματα να δείχνουν μικρές αποκλίσεις στις εκτιμήσεις αλλά διατήρηση της γενικής σταθερότητας και της λειτουργικότητας της μεθόδου.
Συνολικά, η προσέγγιση Lyapunov \textbf{αποδείχθηκε κατάλληλη για διαδικασίες ταυτοποίησης μοντέλων και ελέγχου σε πραγματικό χρόνο}, με δυνατότητες εφαρμογής σε ευρύτερα μη γραμμικά συστήματα.

\end{document}
