function Bend_Energy = Bending_Energy(Aux,C,phi,d_phi,Inv)
global Invariant Material Constraints

if(nargin < 5)
    Inv=Invariant;
end

if(nargin==3)
    d_phi=zeros(size(phi));
    Lmd=0;
else
    Lmd = Material.Lmd;
end

if(isfield(Constraints,'f0'))
    f0 = Constraints.f0;
else
    f0 = 0;
end

if(isfield(Constraints,'f0z'))
    f0z = Constraints.f0z;
else
    f0z = 0;
end

if(isfield(Material,'phi_max'))
    tmp = (Material.k0/2* Aux.HmC02) .* C(:,1) .* Aux.au .* Inv.gweight;
else
    tmp = (Material.k0/2* Aux.HmC02) .* C(:,1) .* Aux.au .* Inv.gweight;
end

tmp = tmp + Lmd/2*d_phi.^2./Aux.au.* C(:,1).* Inv.gweight;
Bend_Energy= 2*pi*sum(tmp);