function[K_TE,Kdel_TE, A_TE,A,K,Kdel,K_TM,Kdel_TM, A_TM,x,y,xnodes,ynodes]=FEM2v4(lamda);

% Generation of the Node Connectivity Table:
% ==========================================

N_elements=18;
n_1st_local_node=[1 2 2 2  3  3 5  5 6  6 7  8  9  9 10 10 11  12];
n_2nd_local_node=[2 6  3 7 4 8 6 10 7 11 8 12 10 14 11 15 12 16];
n_3rd_local_node=[ 5 5 7 6 8 7 10 9 11 10 11 11 14 13 15 14 15 15];
mu=ones(1:N_elements);
eps=ones(1:N_elements);
ko=2*pi/lamda;	
	
	for e=1:N_elements;
	
	n(1,e)=n_1st_local_node(e);
	n(2,e)=n_2nd_local_node(e);
	n(3,e)=n_3rd_local_node(e);
	end;
	
non_cond=[6 7 10 11]; 
% Nodes locations table:
%======================= 
xnodes=[0 .5 1 1.5 0 .5 1 1.5 0 .5 1 1.5 0 .5 1 1.5 ];
ynodes=[0 0  0  0 .25 .25 .25 .25  .5 .5 .5 .5 .75 .75 .75 .75 ];

%Initialization Process:
%=======================

A=zeros(N_elements);
Kdel=zeros(N_elements);
K=zeros(N_elements);
Kedel=zeros(3);

%Loop through all elements:
%==========================

for e=1:N_elements;

%E_z/TM polarization :
%====================
	
	
	pe=1/mu((e));
	qe=eps((e));
	
	
% H_z/TE polarization :
% =========================

	pe=1/eps((e));
	qe=mu((e));
	
	
% Coordinates of the Element nodes:
%==================================

	for i=1:3;
	x(i)=	xnodes(n(i,e));
	y(i)= ynodes(n(i,e));
	end;
	
% Compute the Element matrix entries:
%====================================

Area=.5*abs((x(2)-x(1))*(y(3)-y(1))- (x(3)-x(1))*(y(2)-y(1)));
	
	% A(n(i,e),n(i,j))=0;
	% Kdel(n(i,e),n(i,j))=0;
	% K(n(i,e),n(i,j))=0;
	
	for i=1:3;
	i1=0;
	if i==3;
	i1=3;
	end;
	ip1=(i+1)-i1;
	i2=0;
	if ip1==3;
	i2=3;
	end;	
	ip2=ip1+1-i2;
	bi=y(ip1)-y(ip2);
	ci=x(ip2)-x(ip1);
	
	for j=1:3;
	j1=0;
	if j==3;
	j1=3;
	end;
	jp1=(j+1)-j1;
	j2=0;
	if jp1==3;
	j2=3;
	end;
	jp2=(jp1+1)-j2;
	bj=y(jp1)-y(jp2);
	cj=x(jp2)-x(jp1);
	

Ae(i,j)=-(pe*(bi*bj+ci*cj))/(4*Area);

Kedel(i,j)=-Ae(i,j);

		if  i==j;
			Ke(i,j)=qe*(Area/6);
			Ae(i,j)=Ae(i,j)+(ko^2)*qe*(Area/6);
		else;
			Ke(i,j)=qe*(Area/12);
			Ae(i,j)=Ae(i,j)+(ko^2)*qe*(Area/12);
		end;
		
% Assemble the Element matrices into the Global FEM System:
%==========================================================

	A(n(i,e),n(j,e))=A(n(i,e),n(j,e))+Ae(i,j);
	Kdel(n(i,e),n(j,e))=Kdel(n(i,e),n(j,e))+Kedel(i,j);
	K(n(i,e),n(j,e))=K(n(i,e),n(j,e))+Ke(i,j);
	end;
	end;
	end;

K_TE=K;
Kdel_TE=Kdel;
A_TE=A;

K_TM=K(non_cond,non_cond);
Kdel_TM=Kdel(non_cond,non_cond);
A_TM=A(non_cond,non_cond);
	
eig_squares_TE=eig(Kdel_TE,K_TE);
eig_values_indeces_TE=find(eig_squares_TE >= 0);
eig_values_TE=sqrt(eig_squares_TE(eig_values_indeces_TE));
eig_values_TE=sort(eig_values_TE);
figure(1)
plot(	eig_values_TE,'*');
xlabel('The Modes (TE case)')
ylabel('The eigen_values');
title('Eigen values for  TE modes in a rectangular waveguide');
grid;
legend('a=1.5 cm and b=.75 cm');


figure(2)
eig_squares_TM =eig(Kdel_TM,K_TM);
eig_values_indeces_TM=find(eig_squares_TM >= 0);
eig_values_TM=sqrt(eig_squares_TE(eig_values_indeces_TM));
eig_values_TM=sort(eig_values_TM);
plot(	eig_values_TM,'o');
xlabel('The Modes (TM case)')
ylabel('The eigen_values');
title('Eigen values for  TM modes in a rectangular waveguide');
grid;
legend('a=1.5 cm and b=.75 cm');

e_no=1:16;
table=[e_no' xnodes' ynodes'];

diary data.data

% The General A matrix
A
% The General K matrix
K
% The General del matrix
Kdel
%***************
% The TE A matrix
A_TE
% The TE K matrix
K_TE
% The TE Kdel matrix
Kdel_TE
%******************
% The TM A matrix
A_TM
% The TM K matrix
K_TM
% The TM Kdel matrix
Kdel_TM
%******************
% The next table indiicates the node coordinates.
% First column ====> Node Numbers.
% Second column ===> x coord.
% Third column  ===> y coord.
table
diary