function [tV,dV_dCP] = Volume_functional_support(Aux)
global Invariant ini_Val


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

tV = 2*pi*sum((Aux.dv_supp).*Invariant.gweight);

if nargout==2
    
    %Derivatives of V (local)
    dV_dr_l = (Aux.MC11.*Aux.MC2.*Invariant.Nshape...
        +Aux.MC1.*Aux.MC2.*Invariant.dNsh).*Aux.Mgweight;
    dV_dz_l = Aux.MC1.*Aux.MC11.*Invariant.Nshape.*Aux.Mgweight;
    
    ns = length(dV_dr_l)/gaussno;
    vs = [1:1:ns];
    dV_dCPr = reshape(sum(reshape(dV_dr_l,gaussno,ns,p+1),1),ns,p+1);
    dV_dCPz = reshape(sum(reshape(dV_dz_l,gaussno,ns,p+1),1),ns,p+1);
    dV_dCP = zeros(2,ini_Val.ncontrol);
    
    for i = 0:p
        dV_dCP(1,vs(:)+i) = dV_dCP(1,vs(:)+i) + dV_dCPr(:,i+1)';
        dV_dCP(2,vs(:)+i) = dV_dCP(2,vs(:)+i) + dV_dCPz(:,i+1)';
    end
    dV_dCP=2*pi*dV_dCP;
end
