Conversion from float to complex?

New to Julia
1

I currently have code that takes two arrays, makes a grid of them, and then calculates the square root of the squares of each value:

hx = 1
hy = 1
Nx = 16
Ny = 16
x = hx .* (0:Nx-1)
y = hy .* (0:Ny-1)

function meshgrid(x, y)
    X = [i for i in x, j in 1:length(y)]
    Y = [j for i in 1:length(x), j in y]
    return X, Y
end

xx, yy = meshgrid(x, y)

@. R = sqrt((xx - 0.5)^2 + (yy - 0.5)^2)

However, xx and yy are both matrices of float64s, while R is calculated as complex values. Any idea why the change in variable type? All of the imaginary parts are 0, as expected. Would the best way to fix this be to change the variables back to floats (and how?) or to use a different function? It is not clear to me why my code for R is returning a matrix of complex numbers. Thanks!

2

As of right now, your code does not run in an empty julia REPL. This is because your line

@. R = sqrt((xx - 0.5)^2 + (yy - 0.5)^2)

inserts the computed values on the right hand-side into an existing matrix R, which is not defined in your code snippet. In your session, you probably have an existing complex matrix R, which causes the problem. If you change this line to

R = @. sqrt((xx - 0.5)^2 + (yy - 0.5)^2)

The problem should be solved, because this will create a new matrix which is correctly inferred to have real elements.

6 Likes
3

Or use R = @. hypot(xx - 0.5, yy - 0.5), which calls the hypot function, which can avoid spurious underflow/overflow when the arguments are very large/small, compared to sqrt(a^2 + b^2).

3 Likes
4

Thank you so much!

1 Like
5

BTW, the meshgrid function is entirely redundant, you can just directly use x and y with broadcasting (notice the '):

R = @. sqrt((x' - 0.5)^2 + (y - 0.5)^2)

meshgrid-like functions are a Matlabism (I suspect), but isn’t even needed in Matlab (nor in numpy), you can use the same broadcasting mechanism in those, rendering meshgrid functions obsolete.

3 Likes