Ccall fails on **arrays of ComplexF64 (simple MWE)

General Usage
1

I’m trying to modify 2D arrays using ccall, and this fails when the arrays contain complex double. The same code works when restricted to double.

MWE

Consider a C function which copies a 2D real array into a complex array.

#include <complex.h>
typedef double complex dcmplx;  // like ComplexF64?

int sht(double** origin, dcmplx** target, 
        int array_size, int num_components) {

    // copy origin elements to target for each component
    for (int i=0;i<(int)(num_components);++i) {
        for (int j=0;j<(int)(array_size);++j) {
            (target[i])[j] = (origin[i])[j];
        }
    }
    return 1;
}

I put this in a sht.c file and compile it into sht.so. Then I try to call it from Julia, following instructions from this possibly-out-of-date mailing list conversation.

#  set up 
origin = collect(1.0:3.0)
target = zeros(ComplexF64, 3)

pointer_origin = pointer([origin],1)
pointer_target = pointer([target],1)

# perform the call
ccall(
    (:sht, "./sht.so"),
    Cint,
    (Ptr{Ptr{Cdouble}}, Ptr{Ptr{ComplexF64}}, Cint, Cint),
    pointer_origin, pointer_target, 3, 1)

println(target)  # want this to be [1.0, 2.0, 3.0]

This prints out the super weird

Array{Complex{Float64},1}[[]]

This same example appears to work when every dcmplx and ComplexF64 is replaced with double and Cdouble.

Why is this happening? How can I achieve this functionality of being able to modify complex arrays in Julia using a C library?

Other Details

  • Julia 1.3.1
  • This MWE is based off of a problem I’ve encountered trying to assist in wrapping a C library that performs fast spherical harmonic transforms for Julia, hence the requirement for 2D arrays of complex numbers.
  • The 2D arrays, i.e. double**, are very natural in this context, since depending on spin, a field on a sphere may be defined by more than one component.
  • I’ve put these two snippets into a repo if you want to try running them, with a simple makefile. There’s also a main.c in there to confirm that the C function works.
  • Sorry if this question has been answered before – I think there’s a fair number of past discussions, which helped me get this working for doubles. I’m not sure this has been discussed for complex doubles, though I don’t know why they’re different.
2

… sign…

I’ve said it many times and I’ll say it again. This is one of the MOST misunderstood concept in C and C++. It’s common enough that I’m usually surprised when finding someone that knows it on the internet… It’s not your fault since basically everyone, me included, misunderstand it at some point thanks to C and C++'s infamous array to pointer decay rule.

The thing that surprise a lot of people is that, C and C++ do have propery array type and propery 2D array type. An array in C is completely different from a pointer and a multi-dimensional array in C and C++ is NOT, nor is it equivalent to, a cascade of pointers, no matter what the indexing syntax make you think.

A n-dimensional array, in any language that properly implement them, including C, C++ and julia, is a single pointer to a contagious piece of memory, there’s no double or multiple indirection anywhere.

Now on to the question of what you actually need. What does your C function want, does it really take a double**? Or does it take a 2D array of double.

If it takes a 2D array, you should just pass a pointer instead of a pointer to pointer. That’s the most natural way to do it.

If it takes a pointer to pointer, i.e. not a 2D array but an C array of C pointers, like your sht does, you need to construct just that, an array of pointers. It’s annoy to do and you need to manually preserve the objects you’ve got pointer from. The object you will pass directly as ccall arguments would be [pointer(origin)] and [pointer(target)] and you’ll most likely need to use GC.@preserve origin target ccall(...) to make sure the two objects are valid during the ccall…

4 Likes