function [tD,dD_dCP] = Incremental_cortex_diss(Aux,C,C_1,C1)
global Invariant ini_Val Material dt psi_m

xi = Material.xi_Cor_perp;

p = ini_Val.p;
gaussno = ini_Val.gaussno;

dr = C(:,1)-C_1(:,1);
dz = C(:,2)-C_1(:,2);


%tmp_s = Material.tmp_s;
%h_cut_off = Material.h_cut_off;
%fact = 4;
%tmp_z  = 0.5*(1-tanh(tmp_s*fact*(C(:,2)-h_cut_off*fact)));



%tDL2 = sum((-dr .*C1(:,2) + dz.* C1(:,1)).^2./Aux.au.*C(:,1).* Invariant.gweight .* tmp_z);
tDL2 = sum((-dr .*C1(:,2) + dz.* C1(:,1)).^2./Aux.au.*C(:,1).* Invariant.gweight .* psi_m);

tD = tDL2* pi/dt * xi;

if nargout == 2
    tempA = (- C1(:,2) .* dr + C1(:,1).*dz).^2*ones(1,p+1);
    tempB = 2*(-dr.*C1(:,2) + dz.*C1(:,1)).*(C(:,1)./Aux.au)*ones(1,p+1);
    
    tempC = (1./Aux.au)*ones(1,p+1);
    tempD = (C(:,1)./Aux.au.^2)*ones(1,p+1);
    tempE = C1(:,1)*ones(1,p+1);
    tempF = C1(:,2)*ones(1,p+1);
    
    
    Auxx = Aux.Mgweight;
    Auxx(:,1) = Auxx(:,1) .* psi_m;
    Auxx(:,2) = Auxx(:,2) .* psi_m;
    Auxx(:,3) = Auxx(:,3) .* psi_m;
    Auxx(:,4) = Auxx(:,4) .* psi_m;
    dD_dr_l2 = ((tempC.*Invariant.Nshape - tempD.*Aux.da_dr_l) .* tempA ...
        + tempB.*(-tempF.*Invariant.Nshape + dz*ones(1,p+1).* Invariant.dNsh)).*Auxx;
    dD_dz_l2 = (- tempD.*Aux.da_dz_l.* tempA ...
        + tempB.*(-dr*ones(1,p+1).*Invariant.dNsh + Invariant.Nshape.*tempE)).*Auxx;

    ns = length(dD_dr_l2)/gaussno;
    vs = [1:1:ns];

    
    dD_dCPr2 = reshape(sum(reshape(dD_dr_l2,gaussno,ns,p+1),1),ns,p+1);
    dD_dCPz2 = reshape(sum(reshape(dD_dz_l2,gaussno,ns,p+1),1),ns,p+1);
    dD_dCP2 = zeros(2,ini_Val.ncontrol);
    
    for i = 0:p
        dD_dCP2(1,vs(:)+i) = dD_dCP2(1,vs(:)+i) + dD_dCPr2(:,i+1)';
        dD_dCP2(2,vs(:)+i) = dD_dCP2(2,vs(:)+i) + dD_dCPz2(:,i+1)';
    end
    
    
    dD_dCP = (dD_dCP2 * pi/dt*xi);
end
