// Forward substitution    
// Note that we can simplify by precalculating a[i-1]/b[i-1]
  for (int i=1; i < n; i++) {
     b[i] = b[i] - (a[i-1]*c[i-1])/b[i-1];
     f[i] = g[i] - (a[i-1]*f[i-1])/b[i-1];
  }
  x[n-1] = f[n-1] / b[n-1];
  // Backwards substitution                                                           
  for (int i = n-2; i >= 0; i--) {
     f[i] = f[i] - c[i]*f[i+1]/b[i+1];
     x[i] = f[i]/b[i];
  }