title: Scipy Newton Krylov Expects Square Matrix
body: <p>I am trying to use scipy.optimize.newton_krylov() to solve a least-squares optimization problem, i.e. finding x such that (Ax - b)**2 = 0. My understanding is that A has to be mxn with m>n, b has to be mx1, and x will be nx1. When I try to run the optimization, I get an error:</p>
  
  <pre><code>ValueError: expected square matrix, but got shape=(40, 6)
  </code></pre>
  
  <p>Presumably this error concerns the computation of the Jacobian and not my input matrix A? But if so, how can I change the values I am providing to the functions to resolve this problem? Any advice would be appreciated.</p>
  
  <p>The following code reproduces the error:</p>
  
  <pre><code>import numpy as np
  from scipy.optimize import newton_krylov
  
  A = np.random.uniform(0, 1, (40,6))
  b = np.arange(40)
  x0 = np.ones(6)
  
  def F(x):
      return (A.dot(x) - b)**2
  
  x = newton_krylov(F, np.ones(6))
  </code></pre>