% Uppgift 2A.1
%a)
a = 0.01; b = 0.02; % enheten meter
Ti = 723; % Dirichlet, enheten K
Te = 293; % Neumann, enheten K
K = 315;
N = 400; % antal delintervall
h = (b-a)/N; % nätstorlek
r = a:h:b; % punkter xj, b+h är spökpunkt;
u_exact = 0; % exakt lösning
% Matrisen A
A = zeros(N+2,N+2); %inkludera rad för spökpunkten
for i = 2:N+1
A(i,i-1)=((1/h.^2)-1/(2.*h.*r(i)));
A(i,i) = -2/h.^2;
A(i,i+1)=((1/h.^2)+(1/(2.*h.*r(i))));
end
A(1,1) = 1; % För Dirichlet RV
A(N+2, N) = -1/(2*h);
A(N+2,N+1) = K;
A(N+2,N+2) = 1/(2*h);
f = [0]';
f(1) = Ti;
f(N+2)=K*Te;
ftrans = f';
T = A\ftrans;
plot( r(1:N+1) , T(1:N+1) ); %plotta inte spökpunkten!
z = find(r==0.02);
T(z)
%% b)
% K100 ska finnas för att yttertemp skall ej nå 100 grader celsius (373 K)
% u(0.02) ska vara en funktion av K, då K ändras kommer u(0.02) att ändras
% tills att vi når ett värde som är mindre än 373K.
K1 = 1;
K2 = 10;
K3 = 20;
K4 = 30;
T1 = 717.1204;
T2 = 670.6470;
T3 = 629.6585;
T4 = 596.6963;
T5 = 373;
KI = [K1 K2 K3 K4 ];
TI = [T1 T2 T3 T4 T5];
Points = [1 2 3 4 5];
I = interp1(KI,TI,Points, 'cubic')
% I calculated by hand and the interpoliation is not linear, så i try cubic but do not know. I want to
% find a K where the final temp och T(z) is less than 373 kelvin.
