d_eq    : 'diff(theta,t) = 1 - theta*(V/L);
  ode2(d_eq, theta, t)$
  subst(L/V[ref]*exp(alpha), theta, %)$
  % / L * V[ref]$
  ic1(%, t = 0, alpha = alpha[0])$
d_alpha : log(expand(%));

r_eq : 'diff(theta,t) = -(V/L)*theta*log(theta*(V/L));
  ode2(r_eq,theta,t)$
  subst(L/V[ref] * exp(alpha),theta,%)$
  ic1(%,t=0,alpha=alpha[0])$
  map(exp,-(V/L)*%)$
  block([logexpand : all], expand(% - log(V/V[ref])))$
r_alpha : expand(logcontract(%));