1D heat transfer Finite Element Analysis

I have written the code as follows for the following question. Please help Im confused at the end bit on how to solve the matrix.

CODE:

% control input data

%---------------------------

Ta=72; %Temp of moving air

T1=180;

l=1; % Distance between nodes [m]

units=4; % Length of rod

p=[0:l:units]'; % Nr of nodes

x = [0:l:units]; % Spatial discretisation

nr_nodes = length(p); % Nr of nodes

syms q1 T2 T3 T4 T5

TT=[T1 T2 T3 T4 T5];

for i = 2:nr_nodes % Creates the element matrix

e(i-1,:) = [i-1 i];

end

%%%--- Creating space ---%%%

nr_elements =length(e); % Nr of elements

KK=zeros(nr_nodes,nr_nodes); % Empty stiffness matrix

M=zeros(nr_nodes,nr_nodes); % Empty mass matrix

F=zeros(nr_nodes,1); % Empty forcing vector

Q=zeros(nr_nodes,1); % Initial forcing vector,

%%%--- Geometries ---%%%

K=120; % condv

pc=1; % Thermal conductivity, Btu/(h-ft-F)

hw = 100; % Heat convection water coefficient, Btu

ha = 50; % Heat conv air

P= 1.5708 ; % pi*D Perimeter

A=0.1963; % Element cross section pi*1/4*(D)^2

qs=hw*(T5-40);

%%%--- Assembling matrices ---%%%

for k = 1:nr_elements % Loop over all elements

nodes=e(k,:); % Identify nodes

l=p(k+1)-p(k); % Calculate distance between nodes

Ke = (((A/144)*K)/(l/12))*[1 -1;-1 1]; % Stiffness entry

Me = ((ha*(P/12)*(l/12))/6)*[2 1;1 2]; % Mass entry

Fe = ((ha*(P/12)*Ta*(l/12))/2)*[1;1]; % Forcing entry

KK(nodes,nodes)=KK(nodes,nodes)+Ke; % Position stiffness entry at identified nodes

M(nodes,nodes)=M(nodes,nodes)+Me; % Position mass entry at identified nodes

F(nodes,1)=F(nodes,1)+Fe; % Position forcing entry at identified nodes

kk=KK+M;

end

CODE ENDS

kk*TT'??? syms is not correct i dont know what to do?

This is what I want to achieve:----

Can you take a look how do i establish thiis ans solve for Temperature T1 T2 T3...as below

[Aymon Naqvi]

which book is this?