% =========================================================================
%   Praktikum MATLAB/Simulink I
%   Technische Universität Darmstadt
%   IAT, FG rtm (Prof. Dr.-Ing. U. Konigorski)
%
%   Copyright (c) 2017, FG Regelungstechnik und Mechatronik,
%                       Technische Universität Darmstadt
%   All rights reserved.
%
%   Redistribution and use in source and binary forms, with or withou
%   modification, are permitted provided that the following conditions are
%   met:
%
%   1. Redistributions of source code must retain the above copyright
%      notice, this list of conditions and the following disclaimer.
%
%   2. Redistributions in binary form must reproduce the above copyright
%      notice, this list of conditions and the following disclaimer in the
%      documentation and/or other materials provided with the distribution.
%
%
%   THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
%   "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
%   LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
%   A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
%   HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
%   SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
%   LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
%   DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
%   THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
%   (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
%   OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

% =========================================================================

%Versuch 1
clear all;
clc;


%% Teil 1.1.1
fprintf('Aufgabenteil 1.1.1 \n')

% Initialisierung
A = [4 -1; -2 1];
B = [2 1 -1];
C = [-1 2;3 0;0 1];
D = [7;4];

%Berechnen der Produkte

%AB
try 
    y = A * B;
    SchreibeMatrix(y, 'AB')
catch
    disp('AB nicht definiert')
end

%AC
try
    y = A * C;
    SchreibeMatrix(y, 'AC');
catch
    disp('AC nicht definiert')
end

%BC
try
    y = B * C;
    SchreibeMatrix(y, 'BC');
catch
    disp('BC nicht definiert')
end

%BA
try
    y = B * A;
    SchreibeMatrix(y, 'BA');
catch
    disp('BA nicht definiert')
end

%CA
try
    y = C * A;
    SchreibeMatrix(y, 'CA');
catch
    disp('CA nicht definiert')
end

%AC^T
try
    y = A * C.';
    SchreibeMatrix(y, 'AC^T');
catch
    disp('AC^T nicht definiert')
end

%% Teil 1.1.3
fprintf('\nAufgabenteil 1.1.3 \n')

y = (B * C) * (-A * D + 3 * D);
SchreibeMatrix(y, '(BC)*(-AD+3D)');


%% Teil  1.1.4

E = [1 3 0; 0 2 1; 2 4 3];


fprintf('\nAufgabenteil 1.1.4 \n')

% Darstellung der inversen von E
if (det(E)~=0)
    disp('Inverse existiert, da |det(E)|>0');
    SchreibeMatrix(inv(E), 'E^-1');
end

% Berechnung aller Vektoren v mit E*v=[1;2;3]
b = [1;2;3];
v = E \ b; % = inv(E)*b;
SchreibeMatrix(v, 'v');


%% Teil 1.1.5


fprintf('\nAufgabenteil 1.1.5 \n')

A = [1 0 2; 3 2 1; 4 1 3];
b = [1;1;7];

%1. Möglichkeit
disp('Lösung von Ax=b durch linksseitige Matrixdivision')
x = A \ b;
SchreibeMatrix(x, 'x');

%2. Möglichkeit
disp('Lösung  von Ax=b durch Multiplikation mit der Inversen')
x = inv(A) * b;
SchreibeMatrix(x, 'x');

%% Teil 1.1.6
fprintf('\nAufgabenteil 1.1.6 \n')

% Initialisierung
z1 = 1 - 1j;
z2 = -1 + sqrt(3) * 1j;
z3 = 3 / (1 - 1j);


% Real- und Imaginärteil
fprintf('Re(z1)=%g , Im(z1)=%g, Arg(z1)=%g\n', real(z1), imag(z1), angle(z1));
fprintf('Re(z2)=%g , Im(z2)=%g, Arg(z2)=%g\n', real(z2), imag(z2), angle(z2));
fprintf('Re(z3)=%g , Im(z3)=%g, Arg(z3)=%g\n', real(z3), imag(z3), angle(z3));

% Berechnungen
y1 = (z1 + z2)^2 / (2 * (z1 + z2^2) + 4 * sqrt(3) * 1j);
fprintf('(z1+z2)^2/(2*(z1+z2^2)+4*sqrt(3)*i)=%s\n',num2str(y1));
y2 = conj(z1 * z2) * (3 + 1j) * abs(z1 * z2^2);
fprintf('conj(z1*z2)*(3+i)*abs(z1*z2^2)=%s\n',num2str(y2));

%% Teil 1.1.7
fprintf('\nAufgabenteil 1.1.7 \n')

F = [1 3 7; 2 4 1;2 8 0];


% charakteristisches Polynom
p = poly(F);
SchreibeMatrix(p, 'Koeffizienten des char. Polynoms von F (abfallende Potenzen)')

% Nullstellen
SchreibeMatrix(roots(p), 'Nullstellen')

% Eigenwerte und Eigenvektoren
[V, D] = eig(F);
SchreibeMatrix(eig(D), 'Eigenwerte');
SchreibeMatrix(V, 'Eigenvektoren');

%% Teil 1.1.8
fprintf('\nAufgabenteil 1.1.8 \n')

A = zeros(3, 3, 3);
A(:,:,1) = ones(3, 3);
SchreibeMatrix(A(:,:,1), 'A(:,:,1)')
A(:,:,2) = eye(3,3);
SchreibeMatrix(A(:,:,2), 'A(:,:,2)')
A(:,:,3) = triu([3 3 3;3 3 3;3 3 3]);
SchreibeMatrix(A(:,:,3), 'A(:,:,3)')


%% Teil 1.1.9
fprintf('\nAufgabenteil 1.1.8 \n')

x = [6j; 3; -7j];
y = [2; 4; 9j];
z = [-3j; -2j; 1j];

disp('a)')

SchreibeMatrix(x', 'x^H');
SchreibeMatrix(y*z', 'y*z^H');
SchreibeMatrix(y'*z, 'y^H*z');
SchreibeMatrix(y.'*z, 'y^T*z');

%% Teil 1.2.1
fprintf('\nAufgabenteil 1.2.1 \n')

% Graph 4. Ordnung, da 4 Nullstellen
figure();
x = [-6 -4 -2 2 6];
y = [768 0 0 0 0];
t = (-6: 0.01: 6);
y_1 = t.^4 - 2 * t.^3 - 28 * t.^2 + 8 * t + 96;
plot(t, y_1, 'linewidth', 2)
grid on
xlim([-6 6]);
hold on
plot(-6, 768, 'Linestyle', 'none', 'Marker', '.', 'Color', 'black', 'Markersize', 20); 
text(-5.8, 760, '(-6;768)', 'interpreter', 'latex')          

% Stützstellen

p = polyfit(x, y, 4);
SchreibeMatrix(p, 'Koeffizienten des Polynoms:')

%% Teil 1.2.2
fprintf('\nAufgabenteil 1.2.2 \n')

% Koeffizienten
p = [10 -1 3 -7];


% Nullstellen
SchreibeMatrix(roots(p), 'Nullstellen:')

% erste Ableitung
dp = polyder(p);
SchreibeMatrix(dp, '1. Ableitung');

% zweite Ableitung
dp2 = polyder(dp);
SchreibeMatrix(dp2, '2. Ableitung');

x = (-5:5).';
SchreibeMatrix([x, polyval(p, x)], 'Polynom1')
SchreibeMatrix([x, polyval(dp, x)], 'Polynom2')
SchreibeMatrix([x, polyval(dp2, x)], 'Polynom3')

%% Teil 1.2.3
fprintf('\nAufgabenteil 1.2.3 \n')

% 1. Zerlegung
zaehler = [3 9];
nenner = [6 3 4];

[r, p, k] = residue(zaehler, nenner);
fprintf('1. Zerlegung:\n')
SchreibeMatrix(r, 'ri')
SchreibeMatrix(p, 'pi')

% 2. Zerlegung
zaehler = [5 9];
nenner = [4 6 0 3 0];

[r, p, k] = residue(zaehler, nenner);
fprintf('2. Zerlegung:\n')
SchreibeMatrix(r, 'ri')
SchreibeMatrix(p, 'pi')

%% Teil 1.2.4
fprintf('\nAufgabenteil 1.2.4 \n')

% 1. Zerlegung
nenner = [3 9];
zaehler = [6 3 4];


[r, p, k] = residue(zaehler, nenner);
fprintf('1. Zerlegung:\n')
SchreibeMatrix(r, 'ri')
SchreibeMatrix(p, 'pi')
SchreibeMatrix(k, 'kj')

% 2. Zerlegung
nenner = [5 9];
zaehler = [4 6 0 3 0];

[r, p, k] = residue(zaehler, nenner);
fprintf('2. Zerlegung:\n')
SchreibeMatrix(r, 'ri')
SchreibeMatrix(p, 'pi')
SchreibeMatrix(k, 'kj')

%% Teil 1.3.2
fprintf('\nAufgabenteil 1.3.2 \n')

% Initialisierung
k1 = 600;
k2 = 800;
m1 = 2000;
m2 = 2000;
m3 = 5000;

% Matrizen der Bewegungsgleichung Mq''+Kq=0
M = [m1 0 0;0 m2 0; 0 0 m3];
K = [k1 -k1 0; -k1 (k1+k2) -k2; 0 -k2 k2];

% spezielles EW-Problem Ax=lambda*x mit A=M^-1*K und x=q
A = M \ K;

% Eigenwerte
lambda = eig(A);
[x, lambda] = eig(A);
SchreibeMatrix(lambda, 'Eigenwerte des speziellen Eigenwertproblems')
SchreibeMatrix(x, 'Eigenvektoren des speziellen Eigenwertproblems')

% Bestimmung der mechanischen Eigenwerte durch Rücksubstitution -w^2=lambda
w=[sqrt(-lambda), -sqrt(lambda)];
SchreibeMatrix(w, 'mechanischen Eigenwerte')
disp('Die Eigenvektormatrix bleibt gleich');




