function A = matrix_AP2(Aux,C,C1,C2)


global Invariant ini_Val Material Constraints dt

if isfield(ini_Val,'iopen')
    iopen = ini_Val.iopen;
else
    iopen = 0;
end
 


A=zeros(ini_Val.ncontrol);
connect = Invariant.connect;
%aux =  (C(:,1)./ Aux.au) .* (Material.k0*Constraints.C0* Aux.Hpu).* Invariant.gweight*dt;

aux_1 = (C1(:,1).*Aux.au- C(:,1).*Aux.dau)./Aux.au.^2.* (Material.k0*Constraints.C0* Aux.Hu).* Invariant.gweight*dt;
aux_2 = (C(:,1)./ Aux.au).* (Material.k0*Constraints.C0* Aux.Hu).* Invariant.gweight*dt;


for iel=1:size(connect,1)
    ig=(1:ini_Val.gaussno);
    igglob=(iel-1)*ini_Val.gaussno+ig;
    A1_loc=Invariant.Nshape(igglob,:)'*diag(aux_1(igglob))*Invariant.dNsh(igglob,:);
    A2_loc=Invariant.dNsh(igglob,:)'*diag(aux_2(igglob))*Invariant.dNsh(igglob,:);
    A3_loc=Invariant.Nshape(igglob,:)'*diag(aux_2(igglob))*Invariant.ddNsh(igglob,:);
    A(connect(iel,:),connect(iel,:)) = ...
        A(connect(iel,:),connect(iel,:)) + A1_loc + A2_loc + A3_loc;
end