function L = constraint_surf_local(Aux,C,C1,C2,N_hat,connect_l)
global Invariant ini_Val Material
%
L=zeros(ini_Val.ncontrol-1,2*ini_Val.ncontrol);
connect=Invariant.connect;
%
A_m = C1(:,1).*(C1(:,1).*Aux.au-C(:,1).*Aux.dau)./(C(:,1).*Aux.au.^2)+C2(:,1)./Aux.au+Aux.Hu.*C1(:,2);
B_m = C1(:,2).*(C1(:,1).*Aux.au-C(:,1).*Aux.dau)./(C(:,1).*Aux.au.^2)+C2(:,2)./Aux.au-Aux.Hu.*C1(:,1);
for iel=1:size(connect,1)
    ig=(1:ini_Val.gaussno);
    igglob=(iel-1)*ini_Val.gaussno+ig;
    % r component
    aux1=A_m(igglob).*Invariant.gweight(igglob).*C(igglob,1);
    aux2=C1(igglob,1)./Aux.au(igglob).*C(igglob,1).*Invariant.gweight(igglob);
    L_loc=N_hat(igglob,:)'*diag(aux1)*Invariant.Nshape(igglob,:) ...
         +N_hat(igglob,:)'*diag(aux2)*Invariant.dNsh(igglob,:);
    L(connect_l(iel,:),2*connect(iel,:)-1) = ...
        L(connect_l(iel,:),2*connect(iel,:)-1) + L_loc;
    % z component
    aux1=B_m(igglob).*Invariant.gweight(igglob).*C(igglob,1);
    aux2=C1(igglob,2)./Aux.au(igglob).*C(igglob,1).*Invariant.gweight(igglob);
    L_loc=N_hat(igglob,:)'*diag(aux1)*Invariant.Nshape(igglob,:) ...
         +N_hat(igglob,:)'*diag(aux2)*Invariant.dNsh(igglob,:);
    L(connect_l(iel,:),2*connect(iel,:)) = ...
        L(connect_l(iel,:),2*connect(iel,:)) + L_loc;
end
L = -L*2*pi;