#Worksheet to start analysis of dry ice model in Section 3.5.1.

#The data, in time (seconds)/mass (grams) pairs.

data = [[0., 25.525], [120.0, 24.512], [240.0, 23.524], [360.0, 22.639], [480.0, 21.765], [600.0, 20.89], [720.0,
20.043], [840.0, 19.221], [960.0, 18.431], [1080.0, 17.677], [1200.0, 16.936], [1320.0, 16.22], [1440.0, 15.548],
[1570.0, 14.828], [1680.0, 14.213], [1800.0, 13.553], [1930.0, 12.91], [2040.0, 12.331], [2170.0, 11.689], [2280.0,
11.188], [2410.0, 10.566], [2530.0, 10.043], [2690.0, 9.377], [2780.0, 9.011], [2880.0, 8.616], [3060.0, 7.945],
[3220.0, 7.404], [3380.0, 6.877], [3480.0, 6.593], [3600.0, 6.244]]

#Number of data points
N = len(data);

#A plot:
plt1 = scatter_plot(data);
plt1.axes_labels(['time (seconds)','mass (grams)'])
show(plt1)

#A linear model (bad) of the form u(t) = b - m*t can be fit as follows:
var('t m b');
u(t,m,b) = b - m*t;
SS = function('SS')(m,b); #Sum of squares function, depends on m, b.
SS(m,b) = add((u(data[i][0],m,b)-data[i][1])^2 for i in range(N));

#A plot of SS helps in finding the minimum
plot3d(log(SS), (m,0,0.01), (b,10,30),aspect_ratio=[1000,1,1]).show()

#Then minimize SS with respect to m and b.

#Redo with your own (better/physically motivated) model for u(t)!