Converting ```ComplexF64``` array to ```Float64``` array in-place?

General Usage
, , ,
1

I have a calculation that manipulates and produces a ComplexF64 array (whose dimension varies by usage in different places in my codebase). I know the final result must be an array of real numbers, or suppose I only care about the real part of the result. Is there any way to convert the complex array to a Float64 array in-place, without a new memory allocation?

For example, I suppose a strided array could be hacked to point to the same memory and access every other Float64, but I don’t know how to go about this, and how to make sure that doesn’t interfere with garbage collection. Alternatively, is there a way to copy the real parts to the first ~0.5 of the array’s memory, and have that pointed to by a regular Array{Float64}? (And perhaps potentially releasing the other half of memory)?

I can imagine wrapping the complex array in a way that returns real values lazily, but I suppose that would incur some performance penalty when using that array downstream, and would require the implementation of said wrapper.

3

I guess you could reinterpret to Float64, move every odd-indexed element into the previous position, and finally resize!

Oops, not into the previous position, but move back by an appropriate amount (gradually increasing).

2 Likes
4

If you can use StructArrays.jl, the real component of each complex number will be stored contiguously in a single array. Extracting them is then a simple matter of accessing the relevant component in the StructArray.

4 Likes
5

something like this?

cv=rand(Complex{Float64},10)
@view reinterpret(Float64, cv)[1:2:end]

@view reinterpret(Float64, cv)[1:2: div(end,2)]
6

Since you mention multidimensional arrays, here’s what I’ll suggest:

julia> x = rand(ComplexF64,3,3)
3×3 Matrix{ComplexF64}:
 0.188452+0.660617im  0.226408+0.288041im  0.267183+0.144009im
 0.559536+0.254207im  0.443003+0.949891im  0.780902+0.411449im
 0.575998+0.902281im  0.809989+0.872257im  0.426678+0.767691im

julia> selectdim(reinterpret(reshape,Float64,x),1,1)
3×3 view(reinterpret(reshape, Float64, ::Matrix{ComplexF64}), 1, :, :) with eltype Float64:
 0.188452  0.226408  0.267183
 0.559536  0.443003  0.780902
 0.575998  0.809989  0.426678

You could copy the real values into the first half of the array and reinterpret that, but there’s almost no benefit.

However, if you own the code that’s filling the complex array, I’ll suggest you go a step further and just throw out the imaginary values at the time you produce them so that the complex array never exists at all. That’s the best solution.

1 Like
7

@mikmoore That looks like a great solution. I thought you could only reinterpret scalars and not arrays - should have read the documentation more thoroughly! Regarding your point about not generating a complex array to begin with - this is certainly something to look into (I am passing in and out of fourier transforms, so unless I adapt to pass through cosine transforms etc, I don’t see an obvious way around it).

@ToucheSir That also seems like a very elegant solution! I especially like that I will be able to extract from it a plain regular array and not a strided view. Any idea how FFTW.jl might behave if I fft! a complex SructArray? (i.e. maybe the extra performance burden there might make this not worthwhile.)