% 2D heat flow in L-shape domain
clear
format compact

% data
k=0.7;     % conductivity coefficient

% discretization
Nodes = 15;  % number of nodes
Nel = 18;    % number of elemenets
Ndofs = Nodes;    % one dof at each node

Coord(1,1:2)=[8,0];  % coordinates of nodes
Coord(2,1:2)=[6,2];
Coord(3,1:2)=[8,4];
Coord(4,1:2)=[5,4];
Coord(5,1:2)=[4,3];
Coord(6,1:2)=[4,0];
Coord(7,1:2)=[0,0];
Coord(8,1:2)=[2,2];
Coord(9,1:2)=[3,4];
Coord(10,1:2)=[4,4];
Coord(11,1:2)=[4,5];
Coord(12,1:2)=[0,4];
Coord(13,1:2)=[2,6];
Coord(14,1:2)=[4,8];
Coord(15,1:2)=[0,8];

Inc = [1 3 2;     % nodes of elements (topology)
    2 3 4;
    1 2 6;
    2 4 5;
    2 5 8;
    6 2 8;
    6 8 7;
    8 5 9;
    5 4 10;
    5 10 9;
    9 10 11;
    8 9 13;
    8 13 12;
    8 12 7;
    12 13 15;
    15 13 14;
    14 13 11;
    11 13 9];
Edof=[[1:Nel]',Inc];

temp=40;       % Dirichlet b.c.
bc=[1 temp; 6 temp*1.5; 7 temp*2; 12 temp*1.5; 15 temp];

% plot discretization
figure(1)   % elements
[Ex,Ey]=coordxtr(Edof,Coord,[1:Nodes]',3);
eldraw2(Ex,Ey,[1,1,1],Edof), axis([-1,9,-1,9]), axis equal;

figure(2)  % nodes
eldraw2(Ex,Ey,[3,4,0]), hold on, axis([-1,9,-1,9]), axis equal;
for in=1:Nodes
    h=text(Coord(in,1),Coord(in,2),int2str(in));
    set(h,'fontsize',14,'fontweight','bold');
end

% assembling
KG=zeros(Ndofs,Ndofs); % global matrix
f=zeros(Ndofs,1);      % global vector
f(3)=-15;
f(4)=-20;
f(10)=-10;
f(11)=-20;
f(14)=-15;
for iel=1:Nel
    Kel=flw2te(Ex(iel,:),Ey(iel,:),[1],[k 0;0 k]);
    KG=assem(Edof(iel,:),KG,Kel);
end

% solution
u=solveq(KG,f,bc);

% plot solution
figure(3)
Ed = extract(Edof,u);
fill(Ex',Ey',Ed')
axis([-1,9,-1,9]), axis equal, colorbar;

Tmin=min(u), Tmax=max(u)

Tmin =
    3.4999
Tmax =
    80

Published with MATLAB® 7.6