Looked more into why we get back Python lists under Blender, but not …

…when

running the same code under plain Python (where it returns numpy arrays).
No clear cause found, see #2 for details

readme.md

@@ -97,14 +97,18 @@ far more efficient compiled code compared to Python's interpreted execution.
 Possible optimizations on the Julia side:
 
 - Performance annotations, such as `@inbounds`, `@fastmath` and `@simd`
+
 - Reduce the number of allocations. E.g. for the case above `@time` reports
   `0.255132 seconds (277.99 k allocations: 44.658 MiB, 16.89% gc time)`. The
   high number of allocations is partly caused by having separate `HalfEdge` instances.
   The Bunny model has 104288 edges, leading to roughly double that number of `HalfEdge`'s
   being allocated. Pre-allocating all needed `HalfEdge` instances in a single 
   array would be possible (`sum(loops)` gives the required number). 
+
 - Look more into type stability, `@code_warntype`, etc
+
 - ~~Using StaticArrays.jl in strategic places~~
+
 - Using a 2D array instead of a 1D array for holding vertices, which would
   make get_vertex and set_vertex simpler, but this might not matter much.
 
@@ -113,15 +117,19 @@ The transfer of mesh data between Blender and Julia is far from optimal:
 - On the Blender side the mesh data is extracted using calls to `foreach_get()`
   into preallocated NumPy arrays. Perhaps there is a way to directly access
   the underlying data from the mesh object?
+  
 - Even though the Julia code returns a set of `Array{Float32}` values these
-  get turned into Python lists when crossing the boundary to Blender. Which are 
-  then turned into NumPy arrays on the Blender side for setting up the result mesh.
-  All in all there's quite a lot of copying going on.
+  get turned into Python lists when crossing the boundary to Blender, for some
+  reason. See #xxx for more details. The lists are then turned into NumPy arrays 
+  on the Blender side for setting up the result mesh. All in all there's quite a 
+  lot of copying going on.
+  
 - Blender (and Python) use 0-based indexing, while Julia uses 1-based indexing.
   We +1/-1 alter the relevant data on the transfer between the two worlds, which
   shouldn't cost that much time, but I don't see an easy way to stick to a single 
   indexing scheme on both sides. Julia has some support for 0-based indexing,
   but it apparently comes with some caveats.
+  
 - PyJulia uses the Julia PyCall package for data conversion between Python and Julia,
   which should be able to use zero-copy transfer between the two, judging by 
   https://github.com/JuliaPy/PyCall.jl#arrays-and-pyarray. But so far that

test/callit.py

@@ -1,13 +1,38 @@
+import time, gc
+import numpy
+
 from julia.api import Julia
 from julia import Main
 
 jl = Julia(compiled_modules=False)
-jl.eval('include("fn.jl")')
+jl.eval("""
+include("fn.jl")
 
-import numpy
+import Base.convert
+
+""")
+
+print('allocating')
+
+print(gc.get_stats())
+
+x = numpy.ones(200*1024*1024, 'float32')
+
+#x = numpy.array([1,2,3,4,5,6], dtype=numpy.float32)
+#print(x)
+
+time.sleep(5)
+
+print('calling into Julia')
+Main.fn(x)
+
+print('Back in Python')
+time.sleep(5)
+
+print('Attempting to clean up')
+del x
+print(gc.get_stats())
 
-x = numpy.array([[1,2,3], [4,5,6]], dtype=np.float32)
+time.sleep(5)
 
-print(x)
-res = Main.fn(x)
-print(x)
+#print(x)

test/fn.jl

@@ -1,7 +1,7 @@
+
 function fn(x::Array{Float32})
     println("array size: $(size(x))");
-    println("max element: $(maximum(x))")
-    println("min element: $(minimum(x))")
-    x[1,1] = 123
-    #return 2x
+
+    sleep(5)
+    x[1:3] .= -1
 end

test/get_array.py

@@ -0,0 +1,18 @@
+import time, gc
+import numpy
+
+from julia.api import Julia
+from julia import Main
+
+jl = Julia(compiled_modules=False)
+jl.eval("""
+include("return_array.jl")
+""")
+
+print('calling into Julia')
+res = Main.fn(123)
+
+print('Back in Python')
+for v in res:
+    print(type(v))
+    print(v)

test/return_array.jl

@@ -0,0 +1,12 @@
+using InteractiveUtils
+
+function fn(n)
+    a = Array{Float32}(undef, 5)
+    a .= n
+    b = Array{UInt32}([n, n, n, n, n])
+    c = Array{UInt32}([n, n, n, n, n])
+    d = Array{UInt32}([n, n, n, n, n])
+    return a, b, c, d
+end
+
+code_warntype(fn, (UInt32,))

test/subdiv.py

@@ -0,0 +1,28 @@
+import numpy
+print(numpy.__file__)
+print(numpy.__version__)
+
+import julia
+print(julia.__file__)
+print(julia.__version__)
+from julia.api import Julia
+
+print('Initializing Julia (this might take a moment the first time)...')
+jl = Julia(compiled_modules=False)
+print('Done!')
+
+# Needs to come after the creation of the Julia instance above
+from julia import Main
+
+jl.eval('include("../catmull-clark.jl")')
+
+vertices = numpy.array([-1.0, -1.0, -1.0, -1.0, -1.0, 1.0, -1.0, 1.0, -1.0, -1.0, 1.0, 1.0, 1.0, -1.0, -1.0, 1.0, -1.0, 1.0, 1.0, 1.0, -1.0, 1.0, 1.0, 1.0], 'float32')
+loop_start = numpy.array([0x00000001, 0x00000005, 0x00000009, 0x0000000d, 0x00000011, 0x00000015], 'uint32')
+loop_total = numpy.array([0x00000004, 0x00000004, 0x00000004, 0x00000004, 0x00000004, 0x00000004], 'uint32')
+loops = numpy.array([0x00000001, 0x00000002, 0x00000004, 0x00000003, 0x00000003, 0x00000004, 0x00000008, 0x00000007, 0x00000007, 0x00000008, 0x00000006, 0x00000005, 0x00000005, 0x00000006, 0x00000002, 0x00000001, 0x00000003, 0x00000007, 0x00000005, 0x00000001, 0x00000008, 0x00000004, 0x00000002, 0x00000006], 'uint32')
+
+res = Main.subdivide(vertices, loop_start, loop_total, loops)
+
+for v in res:
+    print(type(v))
+    print(v)