呼叫朱莉婭SVD從C++

Question

我嘗試使用茱莉亞語言從Qt的C++爲獲得SVD功能結果：

jl_value_t *array_type = jl_apply_array_type(jl_float64_type, 2); 
jl_array_t *x = jl_alloc_array_2d(array_type, matrixForSvd.rows(), matrixForSvd.cols()); 

jl_value_t * bol = jl_box_bool(true); 

double *p = (double*)jl_array_data(x); 

int ndims = jl_array_ndims(x); 

size_t size0 = jl_array_dim(x,0); 
size_t size1 = jl_array_dim(x,1); 

// Fill array with data 
for(size_t i=0; i<size1; i++) 
    for(size_t j=0; j<size0; j++) 
     p[j + size0*i] = matrixForSvd(j,i); 

jl_function_t *func = jl_get_function(jl_base_module, "svd"); 
jl_tupletype_t *y = (jl_tupletype_t*)jl_call1(func, (jl_value_t*)x); 

但是，當我試圖解析y我得到很多的垃圾（這是確定），只有U和V但不S：

double *res = (double*)jl_array_data(y); 

for(int t = 0; t < 50; t++){ 
    cout<<res[t]<<endl; 
} 

輸出：

6.94448e-310 
1.97626e-323 
2.7725e-318 
9.88131e-324 
9.88131e-324 
0 
-0.404554 //U 
-0.914514 
-0.914514 
0.404554 
0 
6.94448e-310 
6.94448e-310 
1.97626e-323 
2.7725e-318 
9.88131e-324 
9.88131e-324 
0 
-0.576048 //V 
0.817416 
-0.817416 
-0.576048 
0 
6.94448e-310 
6.94448e-310 
1.97626e-323 
2.7725e-318 

那麼，我怎樣才能正確地從茱莉亞功能像SVD獲取元組？

Answer 1

svd.cpp：

#include <julia.h> 

#include <iostream> 
#include <vector> 


void print_vector(jl_array_t *vec) { 
    double *data = (double *) jl_array_data(vec); 
    size_t n = jl_array_dim(vec, 0); 

    for(int i = 0; i < n; ++i) 
     std::cout << data[i] << " "; 
    std::cout << std::endl; 
} 

// from column major vector 
void print_2d_matrix(jl_array_t *mat) { 
    double *data = (double *) jl_array_data(mat); 
    size_t m = jl_array_dim(mat, 0); 
    size_t n = jl_array_dim(mat, 1); 


    for(int i = 0; i < m; ++i) { 
     for(int j = 0; j < n; ++j) 
      std::cout << data[i+j*m] << " "; 
     std::cout << std::endl; 
    } 
} 

int main() { 

    // matrix represented as a vector in column-major order 
    std::vector<double> mat = { 
     1, 0, 0, 0, 0, 0, 0, 2, 0, 3, 0, 0, 0, 0, 0, 0, 2, 0, 0, 0 
    }; 

    // initialize Julia 
    jl_init(NULL); 

    jl_value_t* dims = (jl_value_t *) jl_eval_string("(4, 5)"); 

    // make sure dims isn't cleaned up by the Julia gc till we're done with it. 
    JL_GC_PUSH(&dims); 
     // get the svd function 
     jl_function_t *svd = jl_get_function(jl_base_module, "svd"); 

     // build a wrapper around the std::vector data to pass our matrix 
     // to the svd function 
     jl_value_t* array_type = jl_apply_array_type(jl_float64_type, 2); 
     jl_array_t *jl_mat = jl_ptr_to_array(array_type, mat.data(), dims, 0); 
    JL_GC_POP(); 

    // call svd 
    jl_value_t *ret = jl_call1(svd, (jl_value_t*)jl_mat); 

    // make sure we don't lose our return data 
    JL_GC_PUSH1(&ret); 
     jl_array_t *jl_U = (jl_array_t*)(jl_fieldref(ret, 0)); 
     jl_array_t *jl_S = (jl_array_t*)(jl_fieldref(ret, 1)); 
     jl_array_t *jl_V = (jl_array_t*)(jl_fieldref(ret, 2)); 

     std::cout << "M: " << std::endl; 
     print_2d_matrix(jl_mat); 

     std::cout << "U: " << std::endl; 
     print_2d_matrix(jl_U); 

     std::cout << "S: " << std::endl; 
     print_vector(jl_S); 

     std::cout << "V: " << std::endl; 
     print_2d_matrix(jl_V); 

    JL_GC_POP(); 

    jl_atexit_hook(0); 

    return 0; 
} 

編譯：

g++ -std=c++14 -fPIC -I$HOME/local/include/julia svd.cpp -L$HOME/local/lib/julia -ljulia 

運行：

LD_LIBRARY_PATH=$HOME/local/lib/julia JULIA_HOME=$HOME/local/bin ./a.out 

輸出：

M: 
1 0 0 0 2 
0 0 3 0 0 
0 0 0 0 0 
0 2 0 0 0 
U: 
0 1 0 0 
1 0 0 0 
0 0 0 -1 
0 0 1 0 
S: 
3 2.23607 2 0 
V: 
-0 0.447214 -0 0 
0 0 1 0 
1 0 0 0 
-0 0 -0 1 
0 0.894427 0 0 

朱莉婭輸出：

   _ 
    _  _ _(_)_  | A fresh approach to technical computing 
    (_)  | (_) (_) | Documentation: http://docs.julialang.org 
    _ _ _| |_ __ _ | Type "?help" for help. 
    | | | | | | |/ _` | | 
    | | |_| | | | (_| | | Version 0.4.6 (2016-06-19 17:16 UTC) 
_/ |\__'_|_|_|\__'_| | 
|__/     | x86_64-redhat-linux 

julia> mat = [1 0 0 0 2; 0 0 3 0 0; 0 0 0 0 0; 0 2 0 0 0] 
4x5 Array{Int64,2}: 
1 0 0 0 2 
0 0 3 0 0 
0 0 0 0 0 
0 2 0 0 0 

julia> svd(mat) 
(
4x4 Array{Float64,2}: 
0.0 1.0 0.0 0.0 
1.0 0.0 0.0 0.0 
0.0 0.0 0.0 -1.0 
0.0 0.0 1.0 0.0, 

[3.0,2.23606797749979,2.0,0.0], 
5x4 Array{Float64,2}: 
-0.0 0.447214 -0.0 0.0 
    0.0 0.0  1.0 0.0 
    1.0 0.0  0.0 0.0 
-0.0 0.0  -0.0 1.0 
    0.0 0.894427 0.0 0.0) 

julia> 

Answer 2

https://discourse.julialang.org/t/calling-functions-defined-a-custom-module-when-embeded-in-c/8918/3

HOWTO構建陣列漁政元組沒有jl_eval_string：

jl_value_t *types[] = {(jl_value_t*)jl_long_type, (jl_value_t*)jl_long_type}; 
jl_tupletype_t *tt = jl_apply_tuple_type_v(types, 2); 
typedef struct { 
    ssize_t a; 
    ssize_t b; 
} ntuple2int; 
ntuple2int *tuple = (ntuple2int*)jl_new_struct_uninit(tt); 
JL_GC_PUSH1(&tuple); 
tuple->a = 3; 
tuple->b = 2; 
jl_array_t *ary = jl_ptr_to_array(atype, data, (jl_value_t*)tuple, own_buffer); 
JL_GC_POP();