Birch–Murnaghan equation of state

I adapted a function from https://gist.github.com/tboland1/f02c68f2f2d064fdff16e63ed8ad6e86 to perform the Birch–Murnaghan equation of state.

function [fitresult, gof] = createFit(Volume, Energy)
[xData, yData] = prepareCurveData(Volume,Energy);
ft = fittype( 'E=E0+9/16*(B./14703.6).*V0.*((((V0./V).^(1/3)).^2-1).^3.*BP+(((V0./V).^(1/3)).^2-1).^2*(6-4*((V0./V).^(1/3)).^2))', 'independent', 'V', 'dependent', 'E' );
opts = fitoptions( 'Method', 'NonlinearLeastSquares' );
opts.Display = 'Off';
opts.StartPoint = [min(Volume) max(Volume) min(Energy) max(Energy)];
[fitresult, gof] = fit( xData, yData, ft, opts );
figure( 'Name', ' Volume Birch-Murnaghan' );
h = plot( fitresult, xData, yData );
legend( h, 'Energy vs Volume', 'Optimization using the Birch-Murnaghan Equation of State', 'Location', 'NorthEast' );
xlabel('Volume (Angstroms^3)');
ylabel('Energy (eV)');
grid on
end

What is wrong whith my function? Can anyone check it? It does not work and I don't know why.

I will leave here a dataset which can help to calibrate it.

Thank you

load('example.mat');
Volume=example(:,1);
Energy=example(:,2);

[Alex Sha]

Is this the result you want?

Root of Mean Square Error (RMSE): 0.000735608704630528
Sum of Squared Residual: 6.49344199593844E-6
Correlation Coef. (R): 0.999926872433012
R-Square: 0.999853750213665

Parameter	Best Estimate
----------	-------------
e0	-23213.5418089814
b	95476.615142617
v0	1.00127602872857
bp	4.04105215666881

[Miguel Cardoso]

Yes, I just wanna my values fitted according that equation but I wasn't able to run the function. How did you fix it?