function [tU,dU_dCP] = incremental_potential(Aux)
global Invariant ini_Val Material phi dt 


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


b = Material.b;
a = Material.a;

%     y_min = value to control the potential
%     Aux.dz2 = (1./(C(:,2)+y_min)).^2;
%     Aux.dz3 = (1./(C(:,2)+y_min)).^3;
%     Aux.dz4 = exp(-Material.a.*(C(:,2)+y_min));

% potential for interaction membrane - support

tN = 200;

tU = 2*pi*sum(((0.75*Aux.dz2-b.*Aux.dz4))./tN.*Aux.ds.*Invariant.gweight);


if nargout == 2
    %Derivatives of E (local)
    tmp = ((0.75*(1./(Aux.MC2+Aux.y_min)).^2-b.*exp(-a.*(Aux.MC2+Aux.y_min))));
    tmp_1 = (-0.75*2.*(1./(Aux.MC2+Aux.y_min)).^3 + a.*b.*exp(-a.*((Aux.MC2+Aux.y_min)))) ;
    dU_dr_l = (tmp.*(Aux.MC1.*Aux.da_dr_l +Aux.Mau.*Invariant.Nshape)).*Aux.Mgweight;
    dU_dz_l = (tmp_1.*Aux.MC1.*Aux.Mau.*Invariant.Nshape...
        + (tmp.*(Aux.MC1.*Aux.da_dz_l))).*Aux.Mgweight;

    ns = length(dU_dr_l)/gaussno;
    vs = [1:1:ns];
    dU_dCPr = reshape(sum(reshape(dU_dr_l,gaussno,ns,p+1),1),ns,p+1);
    dU_dCPz = reshape(sum(reshape(dU_dz_l,gaussno,ns,p+1),1),ns,p+1);
    dU_dCP = zeros(2,ini_Val.ncontrol);

    for i = 0:p
        dU_dCP(1,vs(:)+i) = dU_dCP(1,vs(:)+i) + dU_dCPr(:,i+1)';
        dU_dCP(2,vs(:)+i) = dU_dCP(2,vs(:)+i) + dU_dCPz(:,i+1)';
    end

    dU_dCP=2*pi.*dU_dCP./tN;


end
