old_contents: import fenics
  
  
  def apply_backward_euler(Delta_t, u):
      """ Apply the backward Euler (fully implicit, first order) time discretization method. """
      u_t = (u[0] - u[1])/Delta_t
      
      return u_t
      
      
  def apply_bdf2(Delta_t, u):
      """ Apply the Gear/BDF2 (fully implicit, second order) backward difference formula for time discretization. 
      
      Here we use the variable time step size formula given by equation (12) from
      
          @article{eckert2004bdf2,
            title={A BDF2 integration method with step size control for elasto-plasticity},
            author={Eckert, S and Baaser, H and Gross, D and Scherf, O},
            journal={Computational Mechanics},
            volume={34},
            number={5},
            pages={377--386},
            year={2004},
            publisher={Springer}
          }
          
      which we interpreted in the context of finite difference time discretizations in a PDE,
      rather than as an ODE solver, by comparing to the constant time step size BDF2 formula in
      
          @article{belhamadia2012enhanced,
            title={An enhanced mathematical model for phase change problems with natural convection},
            author={Belhamadia, YOUSSEF and Kane, ABDOULAYE S and Fortin, ANDR{\'E}},
            journal={Int. J. Numer. Anal. Model},
            volume={3},
            number={2},
            pages={192--206},
            year={2012}
          }    
          
      Perhaps we could find a reference which derives variable time step size BDF2 
      in the context of finite difference time discretization,
      rather than in the context of ODE solvers.
      """
      tau = Delta_t[0]/Delta_t[1]
      
      u_t = 1./Delta_t[0]*((1. + 2.*tau)/(1. + tau)*u[0] - (1. + tau)*u[1] + tau*tau/(1. + tau)*u[2])
      
      return u_t
      
new_contents: """ **backward_difference_formulas.py** implements some BDF's for time discretization.
  """
  import fenics
  
  
  def apply_backward_euler(Delta_t, u):
      """ Apply the backward Euler (fully implicit, first order) time discretization method. """
      u_t = (u[0] - u[1])/Delta_t
      
      return u_t
      
      
  def apply_bdf2(Delta_t, u):
      """ Apply the Gear/BDF2 (fully implicit, second order) backward difference formula for time discretization. 
      
      Here we use the variable time step size formula given by equation (12) from
      
          @article{eckert2004bdf2,
            title={A BDF2 integration method with step size control for elasto-plasticity},
            author={Eckert, S and Baaser, H and Gross, D and Scherf, O},
            journal={Computational Mechanics},
            volume={34},
            number={5},
            pages={377--386},
            year={2004},
            publisher={Springer}
          }
          
      which we interpreted in the context of finite difference time discretizations in a PDE,
      rather than as an ODE solver, by comparing to the constant time step size BDF2 formula in
      
          @article{belhamadia2012enhanced,
            title={An enhanced mathematical model for phase change problems with natural convection},
            author={Belhamadia, YOUSSEF and Kane, ABDOULAYE S and Fortin, ANDR{\'E}},
            journal={Int. J. Numer. Anal. Model},
            volume={3},
            number={2},
            pages={192--206},
            year={2012}
          }    
          
      Perhaps we could find a reference which derives variable time step size BDF2 
      in the context of finite difference time discretization,
      rather than in the context of ODE solvers.
      """
      tau = Delta_t[0]/Delta_t[1]
      
      u_t = 1./Delta_t[0]*((1. + 2.*tau)/(1. + tau)*u[0] - (1. + tau)*u[1] + tau*tau/(1. + tau)*u[2])
      
      return u_t
      
subject: Add missing head docstring