Why the function nroots returns a Matrix{Array{Int64,0}}

General Usage
3

Note that this is a continuation of this thread: How to convert a Matrix{Array{Int64, 0}} to a Matrix{Int64}

You wrote your nroots! function in a classic Matlab-ish “vectorized” style where it assumes that z is an array, and allocates the result nr as an array too. If z is a scalar, then size(z) = () and nr is a zero-dimensional array (which is not the same thing as a scalar in Julia).

In this case, it doesn’t seem like there is any advantage in a vectorized style — even if z is an array, the different z’s don’t seem to share any computation. I would suggest just writing a function that works for scalar z and then broadcasting it as needed. If I understand your code correctly, this should be:

function nroots(p, z; eps=1e-05, maxiter=50)
   dp= derivative(p)
   ddp=derivative(dp)
   prod=ddp*p 
   nr=0
   for _ in 1:maxiter 
       z, zprev=  z-p(z)/(dp(z) - prod(z)/(2*dp(z))), z
       nr += abs(z-zprev)^2 > eps
   end 
   return nr 
end

(Note that neither this nor your original nroots! function should end with a !, by Julia convention, since neither of them mutate their arguments. Reassigning z is not mutation.)

Then you can broadcast this as nr = nroots.(p, z) over an array z and it will return an array of scalars.

4 Likes