function [xh,eh,EHh] = refine(xH,eH)

% [xh,eh,EHh] = refine(xH,eH)

nH = size(xH,2);
mH = size(eH,1);

nh = 2*nH - 1;

xh = [xH(eH(1,1))];
eh = [];

EHh = zeros(nh,nH);

node = 1;

for i=1:mH
    xH1 = xH(eH(i,1));
    xH2 = xH(eH(i,2));
    
    xh2 = (xH1+xH2)/2;
   
    xh = [ xh , xh2 , xH2 ];
    eh = [ eh ; [node,node+1] ; [node+1,node+2] ];
    
    if (node>1 && i>1)
        EHh(node,i) = 1;   
        EHh(node+1,i) = 1/2;
    end

    if (node+2<nh && i+1<nH)
        EHh(node+1,i+1) = 1/2;
        EHh(node+2,i+1) = 1;
    end
    
    node = node + 2;
end
