%%
%----------------------------------------------%
% Exercises on the Systems of linear equations %
%----------------------------------------------%
%% 
%-----------------------------%
% Ex. 1.0.1
clear all
close all
clc

A=[3,1; 4,2];
b=[3,4]';

xsol=A\b
%which is equivalent to
mldivide(A,b)

inv(A)*b

%the inverse of the matrix A is
%defined as inv(A)=(1/det(A))*Ap
Ap=[A(2,2),-A(2,1); -A(1,2), A(1,1)]';

inv(A)
Ap/det(A)

%verify that x is the solution of the equation Ax=b
%by evaluating
A*xsol

%Plot the lines 3*x(1)+x(2)=3 and 4*x(1)+2*x(2)=4
x=linspace(0,2,100);
hold on
plot(x,3-3*x,'r','LineWidth',2)
plot(x,(4-4*x)/2,'g','LineWidth',2)
Legend('3*x_1+x_2=3','4*x_1+2*x_2=4')
plot(xsol(1),xsol(2),'ks','LineWidth',2,'MarkerFaceColor','g',...
                'MarkerSize',10)
xlabel('x_1'),ylabel('x_2')
hold off
%we can check the intersection is in xsol

print -f -depsc 'intersection.eps'

%%
%-----------------------------%
% Ex. 1.0.2
clear all

A=[1,5,4; 7,9,3; 4,5,2];
b=[3,4,9]';

xsol=A\b
inv(A)*b

%compare the two in the long format

format long
xsol
inv(A)*b

%the two operations A\b (mldivide) 
%and inv(A)*b are not equivalent
help mldivide
doc mldivide
help inv
% you will find explained why a way A\b is a better way
% of solving Ax=b than x=inv(A)*b

%%
%-----------------------------%
% Ex. 1.0.3
clear all
close all
clc
% we want to solve the system of linear
% equations xA=b

A=[3,1; 4,2];
b=[3,4];

xsol=b/A
% which is equivalent to
mrdivide(b,A)

b*inv(A)

% verify that x is the solution of the equation Ax=b
% by evaluating
xsol*A

%Plot the lines 3*x(1)+x(2)=3 and 4*x(1)+2*x(2)=4
x=linspace(-7,-3,100);
hold on
plot(x,(3-3*x)/4,'r')
plot(x,(4-x)/2,'g')
plot(xsol(1),xsol(2),'ks','LineWidth',2,'MarkerFaceColor','g',...
                'MarkerSize',10)
xlabel('x_1'),ylabel('x_2')
hold off

% show that xsol'=(b/A)'=A'\b'
xsol'
%(b/A)'
A'\b'
%%
%-----------------------------%
% Ex. 1.0.4
clear all
close all
clc
A=[1 -2; 2 -4]; b=[5 6]';
det(A)
%two parallel lines
pinv(A)*b

x=linspace(-7,-3,100);
hold on
plot(x,(x-5)/2,'r')
plot(x,(x-3)/2,'g')
xlabel('x_1'),ylabel('x_2')
hold off

%%
%-----------------------------%
% Ex. 1.0.5
clear all
close all
clc
% Overdetermined system %
A=[1 2; 3 -4; 3 2]; b=[5 6 8]';
%Plot
x=linspace(-4,8,100);
hold on
plot(x,(5-x)/2,'r')
plot(x,(3*x-6)/4,'g')
plot(x,(8-3*x)/2,'b')
xlabel('x_1'),ylabel('x_2')
hold off
%%
%-----------------------------%
% Example
clear all

% Undetermined system %

format long
A=[3 5 4; 2 1 7]; b=[7 6]';
xsol=A\b

A*xsol
norm(A*xsol-b)

%%
%-----------------------------%
% Ex. 1.0.6
clear all
close all
clc

A=[2 1; 5 2]
b=[3; 4]

% solve the system of linear equations first
c=inv(A)*b
c=A\b
c=mldivide(A,b)

y=c(1)
z=c(2)

% now we can solve z=3 +8x -x^2
x_1=4+sqrt(16+3-z)
x_2=4-sqrt(16+3-z)

% Let's check graphically that we can the same answer

%f1 = [3 - (3+8x-x^2)]/2
%f1 = [4 - 2(3+8x-x^2)]/5

a=0; b=9; N=50;
x=linspace(a, b, N);
f1=(3-(3+8*x-x.^2))/2;
f2=(4-2*(3+8*x-x.^2))/5;

hold on
plot(x,f1,'r','LineWidth',2)
plot(x,f2,'g','LineWidth',2)
Legend('f_1','f_2')
plot(x_1,y,'ks','LineWidth',2,'MarkerFaceColor','g',...
                'MarkerSize',10)
plot(x_2,y,'ks','LineWidth',2,'MarkerFaceColor','g',...
                'MarkerSize',10)
xlabel('x_1'),ylabel('x_2')
hold off

print -f -depsc 'non-linear.eps'

%%
%-----------------------------%
% Ex. 2.0.7
clear all
close all
clc
% electrical circuit
V=10; R1=2; R2=4; R3=3;

A=[-1 1 1; R1/V R2/V 0; R1/V 0 R3/V];

%check that det(A) is different from zero
det(A)

I=inv(A)*[0; 1; 1]
I1=I(1); I2=I(2); I3=I(3);

% clearly is equivalent to
A2=[1 -1 -1; R1/V R2/V 0; R1/V 0 R3/V];
inv(A2)*[0; 1; 1]

% variable resistence
N=200;
R3=linspace(0,100,N);

for i=1:N
    A=[-1 1 1; R1/V R2/V 0; R1/V 0 R3(i)/V];

    if det(A)~=0
        I=inv(A)*[0; 1; 1];
        I1(i)=I(1); I2(i)=I(2); I3(i)=I(3);
    else
        print('det(A) is zero for R3 equal to')
        R3(i)
        break
    end
end

hold on
plot(R3,I1,'r')
plot(R3,I2,'b')
plot(R3,I3,'g')
legend('I1 [A]','I2 [A]','I3 [A]')
xlabel('R3 [\Omega]')
hold off

