% Galerkin method and Exact Solution- Final Question 1a

clear all;
clc;

d = 1e-3; L = 50e-3; T_inf = 30; h = 20; k = 200; T0 = 300;

A = pi*d^2/4; p = pi*d; m = sqrt(h*p/k/A);

syms x T c2 c3
N1 = x*(x-2*L);     N2 = x*(x^2-3*L^2);
T = T0+ c2*N1+c3*N2;

R = k*diff(diff(T))-(p*h/A)*(T-T_inf);

e1 = int(N1*R,0,L);
e2 = int(N2*R,0,L);

[a b] = solve(e1,e2);

c2 =eval(a);     c3 = eval(b);

y = [0 10 20 30 40 50]*10^-3;

for i = 1:6
    x = y(i);
    Temp(i) = eval(T);
    Temp_ex(i) = T_inf+(T0-T_inf)/cosh(m*L)*cosh(m*(L-x));
end

plot(y,Temp,'s--', y, Temp_ex,'o-');
xlabel('x(m)'); ylabel('Temperature (deg C)');
hold on;