new_contents: class AnnotateRight(object):
      def __init__(self, lines, texts):
          self.lines = lines
          self.texts = texts
  
          self.axes = lines[0].axes
  
          self._compute_coordinates()
          self._avoid_collision()
          self._annotate()
  
      def _data_to_axis(self, line):
          ax = line.axes
          xy = line.get_xydata()
  
          xy_fig = ax.transData.transform(xy)
          xy_ax = ax.transAxes.inverted().transform(xy_fig)
  
          return xy_ax
  
      def _compute_coordinates(self):
          self.coordinates = [self._get_last_xy(l) for l in self.lines]
  
      def _avoid_collision(self):
          rtol = 0.02
  
          new_texts = []
          new_coordinates = []
  
          xy_last = None
          for xy, text in zip(self.coordinates, self.texts):
              if (xy_last is None) or (abs(xy_last[1] - xy[1]) > rtol):
                  new_texts.append(text)
                  new_coordinates.append(xy)
              else:
                  new_texts[-1] = ','.join((new_texts[-1], text))
              xy_last = xy
  
          self.coordinates = new_coordinates
          self.texts = new_texts
  
      def _get_last_xy(self, line):
          xy_last = self._data_to_axis(line)[-1]
          xy_last += np.array([0.01, 0])
          return xy_last
  
      def _annotate(self):
          for xy, text in zip(self.coordinates, self.texts):
              self.axes.annotate(text, xy, xycoords='axes fraction',
                  va='center', ha='left', size='small')
  
  
  class TemporalEvolutionGraph(object):
      def __init__(self, ax=None):
          import matplotlib.pyplot as plt
  
          self.ax = ax or plt.gca()
          self.lines = []
  
      def loglog(self, *args, **kwargs):
          line, = self.ax.loglog(*args, **kwargs)
          self.lines.append(line)
  
      def annotate_lines(self, texts):
          AnnotateRight(self.lines, texts)
  
      def set_title(self, text):
          self.ax.text(0.02, 0.98, text, transform=self.ax.transAxes, va='top')
  
      @property
      def figure(self):
          return self.ax.figure
  
  
  import numpy as np
  import atomic
  
  def time_dependent_z(element_name, times, temperature, density):
      rt = atomic.RateEquations(atomic.element(element_name))
  
      y = rt.solve(times, temperature, density)
  
      ne_tau = np.array([ 1e14, 1e15, 1e16, 1e17, 1e18])
      indices = np.searchsorted(y.times, ne_tau/density)
  
      g = TemporalEvolutionGraph()
  
      for i in indices:
          z_mean = y[i].mean_charge()
          g.loglog(temperature, z_mean, color='black', ls='--')
  
      g.ax.set_xlabel(r'$T_\mathrm{e}\ \mathrm{(eV)}$')
      g.ax.set_ylim(0.4, y.atomic_data.nuclear_charge + 4)
      g.annotate_lines(['$10^{%d}$' % i for i in np.log10(ne_tau)])
  
      z_mean = y.y_coronal.mean_charge()
      g.loglog(temperature, z_mean, color='black')
  
      def Uppercase(text):
          return text[0].upper() + text[1:]
      g.set_title(Uppercase(y.atomic_data.element) + r' time dependent $\left<Z\right>$')
  
  
  if __name__ == '__main__':
      import matplotlib.pyplot as plt
      times = np.logspace(-7, 0, 100)
      temperature = np.logspace(0, 3, 100)
      density = 1e19
  
      time_dependent_z('carbon', times, temperature, density)
  
      plt.draw()
      plt.show()
  
subject: Add time dependent Z graph.