filename = 'ex.txt';
delimiterIn = ' ';
headerlinesIn = 5;
matrix = importdata(filename,delimiterIn,headerlinesIn);
A = matrix.data;
A = A(:,1:3);
temp_ids = sum(A(:) == 0.2);
density_ids = length(A)/temp_ids;
density = A((1:density_ids),1);
pressure = A(:,3);
pressure = reshape(pressure,[density_ids,temp_ids]);
[~,dPdD] = gradient(pressure, mean(diff(density)));
[~,ddPddD_grad] = gradient(dPdD, mean(diff(density)));
ddPddD_del2 = 4 *del2(pressure, mean(diff(density)));
subplot(1,2,1)
plot(density, ddPddD_grad)
grid on;
xlim([0.4 0.8])
xlabel('\rho (g/cm^3)');
ylabel('\partial^2p/\partial\rho^2')
subplot(1,2,2)
plot(density, ddPddD_del2)
grid on;
xlim([0.4 0.8])
xlabel('\rho (g/cm^3)');
ylabel('\partial^2p/\partial\rho^2')
temperature = 100:5:300;
density_spacing = density(2,1) - density(1,1);