function [row,col,sk] = stiffcurElast2D (in,xx,iel,conec,locel,t)
%
%  Generacion de la matriz de rigidez
%  Elemento triangulo lineal p/problema
%  de Elasticidad plana
%
%  in:    Numeros de nodo
%  xx:    Tabla de coordenadas
%  iel:   Numeros de elemento
%  conec: Tabla de conectividades
%  locel: Tabla de vectores de localizacion 
%  t: espesor (si t=0: deformacion plana, t>0: tension plana)

nel = length (iel);

inn        = zeros(max(in),1);
inn(in(:)) = 1:length(in);

for i=1:3
    X1(:,i) = [ xx(inn(conec(:,i)),1) ] ;
    Y1(:,i) = [ xx(inn(conec(:,i)),2) ] ;
end

ind  = [ 1, 2, 3, 1, 2];
row  = zeros(36*nel,1);
col  = zeros(36*nel,1);
sk = zeros(36*nel,1);

in1 = 0;
for iel = 1:nel

    dosdelta  = (X1(iel,2)*Y1(iel,3) - X1(iel,3)*Y1(iel,2)) + ...
        (X1(iel,3)*Y1(iel,1) - X1(iel,1)*Y1(iel,3)) + ...
        (X1(iel,1)*Y1(iel,2) - X1(iel,2)*Y1(iel,1)) ;

    for i=1:3
        beta(i)  =   (Y1(iel,ind(i+1) ) - Y1(iel,ind(i+2) ))/dosdelta;
        gamma(i) =  -(X1(iel,ind(i+1) ) - X1(iel,ind(i+2) ))/dosdelta;
    end

    B = [beta(1)       0   beta(2)       0  beta(3)        0
         0       gamma(1)       0  gamma(2)       0   gamma(3)
         gamma(1) beta(1) gamma(2) beta(2) gamma(3)   beta(3)];

    E  = 2100000;
    nu = 0.3;
    if t>0
        C = E/(1-nu^2)*t* ...            %  tension plana
                         [1 nu 0
                          nu 1 0
                          0 0 (1-nu)/2];
    else
        C = E/(1+nu)/(1-2*nu)* ...       %  deformacion plana
                         [1-nu nu 0
                          nu 1-nu 0
                          0 0 (1-2*nu)/2];
    end

    K = B'*C*B*(dosdelta/2);

    for i = 1:6
        for j=1:6
            in1 = in1+1;
            row(in1) = locel(iel,i);
            col(in1) = locel(iel,j);
            sk(in1)  = K(i,j); 
        end
    end
end

