function mat = jac(xx)
global nb nv barras restric massa x0
global extforces totforces tensoes grav
%%% calcula forcas internas no estado xx
dintforces=zeros(nv*2,4);
dbarforces=zeros(2*2,4);
mat=zeros(2*nv,2*nv);
noloc=zeros(2,1);
dxyloc=zeros(2,2,4);
dll=zeros(4);
dxyloc(1,1,1)=1;dxyloc(1,2,2)=1;dxyloc(2,1,3)=1;dxyloc(2,2,4)=1;
%% compute internal forces
for i=1:nb
  noloc(1)=barras(i,1);noloc(2)=barras(i,2);
  xyloc(1,1) = xx(2*noloc(1)-1);
  xyloc(1,2) = xx(2*noloc(1));
  xyloc(2,1) = xx(2*noloc(2)-1);
  xyloc(2,2) = xx(2*noloc(2));
  l0=barras(i,3);
  for k=1:4
    ll=0;
    dll(k)=0;
    for j=1:2
     dd(j)=xyloc(1,j)-xyloc(2,j);
     ddd(j,k)=dxyloc(1,j,k)-dxyloc(2,j,k);
     ll=ll+dd(j)*dd(j);
     dll(k)=dll(k)+2*dd(j)*ddd(j,k);
    end
    ll=sqrt(ll);
    dll(k)=0.5/ll*dll(k);
    kkk=barras(i,4);
    for j=1:2
     auxcos(j)=dd(j)/ll;
     dauxcos(j,k)=ddd(j,k)/ll-dd(j)/ll/ll*dll(k);
    end
    tens=kkk*(ll-l0);
    dtens(k)=kkk*dll(k);
    dbarforces(1,k)=-dtens(k)*auxcos(1)-tens*dauxcos(1,k);
    dbarforces(2,k)=-dtens(k)*auxcos(2)-tens*dauxcos(2,k);
    dbarforces(3,k)=-dbarforces(1,k);
    dbarforces(4,k)=-dbarforces(2,k);
  end
  for ino=1:2
   for ico=1:2
       iunk=(noloc(ino)-1)*2+ico;
       iloc=(ino-1)*2+ico;
     for jno=1:2
      for jco=1:2
         junk=(noloc(jno)-1)*2+jco;
         jloc=(jno-1)*2+jco;
         mat(iunk,junk)=mat(iunk,junk)+dbarforces(iloc,jloc);
      end
     end
    end
   end
end

%%%
%% Restricoes
for i=1:2*nv
    if (restric(i) != 0) mat(i,:)=0; mat(i,i)=-100000000; end
end

% Aqui mat tem o jacobiano
end