[element2faces, face2nodes]=getFaces(elements)

function [element2faces, face2nodes]=getFaces(elements)
%function: [element2faces, face2nodes]=getFaces(elements)
%requires: deleterepeatedrows
%generates faces of (tetrahedral 3D) triangulation defined in elements
%elements is matrix, whose rows contain numbers of its element nodes 
%face2edges returns faces numbers of each tetrahedral element
%face2nodes returns three node numbers of each face
%example: [element2faces, face2nodes]=...
%getFaces([1 3 4 5; 7 4 3 5; 5 7 6 4; 6 8 7 4; 2 1 4 5; 2 4 6 5])

%extracts sets of faces
faces1=elements(:,[2 3 4]);
faces2=elements(:,[3 4 1]);
faces3=elements(:,[4 1 2]);
faces4=elements(:,[1 2 3]);

%as sets of their nodes (vertices)
vertices=zeros(size(elements,1)*4,3);
vertices(1:4:end,:)=faces1;
vertices(2:4:end,:)=faces2;
vertices(3:4:end,:)=faces3;
vertices(4:4:end,:)=faces4;

%repeated sets of nodes (joint faces) are eliminated 
[face2nodes,element2faces]=deleterepeatedrows(vertices);
element2faces=reshape(element2faces,size(elements,2),size(elements,1))';