add support for variable coefficients in celerite models

I'm interested in being able to model different subsets of the data with different kernels, i.e. K = K1 + K2, where K2 could have zeros everywhere except blocks corresponding to the target data subset. As you mentioned this isn't trivial with celerite, but could work in principle, i.e. by generalizing the kernel to be k(tau=|ti-tj|) = a_1(ti) * a_2(tj) * exp(-c*tau) * cos(d*tau) + b_1(ti) * b_2(tj) * exp(-c*tau) * sin(d*tau)