title: Can I speed up this basic linear algebra code?
body: <p>I was wondering whether it is possible to optimise the following using <code>Numpy</code> or mathematical trickery.</p>
  
  <pre><code>def f1(g, b, dt, t1, t2):
    p = np.copy(g)
    for i in range(dt):
      p += t1*np.tanh(np.dot(p, b)) + t2*p
    return p
  </code></pre>
  
  <p>where <code>g</code> is a vector of length <code>n</code>, <code>b</code> is an <code>n</code>x<code>n</code> matrix, <code>dt</code> is the number of iterations, and <code>t1</code> and <code>t2</code>are scalars.</p>
  
  <p>I have quickly ran out of ideas on how to optimise this further, because <code>p</code> is used within the loop, in all three terms of the equation: when added to itself; in the dot product; and in a scalar multiplication.</p>
  
  <p>But maybe there is a different way to represent this function or there are other tricks to improve its efficiency. If possible, I would prefer not to use <code>Cython</code> etc., but I'd be willing to use it if the speed improvements are significant. Thanks in advance, and apologies if the question is out of scope somehow.</p>
  
  <h2>Update:</h2>
  
  <p>The answers provided so far are more focused on what the values of the input/output could be to avoid unnecessary operations. I have now updated the MWE with proper initialisation values for the variables (I didn't expect the optimisation ideas to come from that side -- apologies). <code>g</code> will be in the range <code>[-1, 1]</code> and <code>b</code> will be in the range <code>[-infinity, infinity]</code>. Approximating the output is not an option because the returned vectors are later given to an evaluation function -- approximation may return the same vector for fairly similar input, so it is not an option.</p>
  
  <hr>
  
  <h2>MWE:</h2>
  
  <pre><code>import numpy as np
  import timeit
  
  iterations = 10000
  
  setup = """
  import numpy as np
  n  = 100
  g  = np.random.uniform(-1, 1, (n,)) # Updated.
  b  = np.random.uniform(-1, 1, (n,n)) # Updated.
  dt = 10
  t1 = 1
  t2 = 1/2
  
  def f1(g, b, dt, t1, t2):
    p = np.copy(g)
    for i in range(dt):
      p += t1*np.tanh(np.dot(p, b)) + t2*p
    return p
  """
  
  functions = [
    """
      p = f1(g, b, dt, t1, t2)
    """
  ]
  
  if __name__ == '__main__':
    for function in functions:
      print(function)
      print('Time = {}'.format(timeit.timeit(function, setup=setup,
                                             number=iterations)))
  </code></pre>