从C调用MATLAB复合矩阵函数

Question

我已经使用MATLAB编译器SDK创建了MATLAB函数的DLL MakeBus。在C中使用。该函数将1个整数和2个双矩阵作为输入，并给出一个复杂的双矩阵作为输出。它在MATLAB处理过程中会产生稀疏的矩阵。 这是MATLAB源代码

function Sbus = makeSbus(baseMVA, bus, gen)

on = find(gen(:, 8) > 0); %% which generators are on?
gbus = gen(on, 1);               %% what buses are they at?
    


nb = size(bus, 1);
ngon = size(on, 1);
Cg = sparse(gbus, (1:ngon)', ones(ngon, 1), nb, ngon);  


%% connection matrix.element i, j is 1 if gen on(j) at bus i is ON
    Sbus =  ( Cg * (gen(on,2 ) + 1j * gen(on, 3))   - (bus(:, 3) + 
             1j * bus(:, 4)) ) / baseMVA;

SSBUS的输出是9*1维度的复杂双矩阵。

在此处使用案例9的数据，我从MATLAB获得以下答案。

https://matpower.org/docs.org/docs/ref/ref/matpower5.0/case9。 html

case1='case9.m';
A=loadcase(case1);
Sbus =makeSbus(A.baseMVA,A.bus, A.gen);

sbus = {0.7230 + 0.27030i， 1.63+ 0.0654i， 0.85-0.1095i， 0.0 + 0.0i， -0.90- 0.30i， 0.0+ 0.0i， -1.0-0.35i， 0.0 + 0.0i， -1.25- 0.50i}

我希望从C中调用此代码。使用MATLAB编译器SDK创建必要的库文件。这是我的代码

mxArray* Bus, * Gen,*var;
mxArray* sbus = NULL;

Bus = mxCreateDoubleMatrix(9, 17, mxREAL);
Gen = mxCreateDoubleMatrix(3, 25, mxREAL);
sbus = mxCreateDoubleMatrix(9, 1, mxCOMPLEX);
var = mxCreateDoubleMatrix(1, 1, mxREAL);
int data1[]={100};
double B[9][13], G[3][21];
// reads B and G matrix from a csv file


memcpy(mxGetPr(Bus), B, 9 * 13 * sizeof(double));
memcpy(mxGetPr(Gen), G, 3 * 21 * sizeof(double));
memcpy(mxGetPr(var), data1, sizeof(int));
    

if (!makeSbusInitialize()) 
   {
        fprintf(stderr, "Could not initialize the library.\n");
        return -2;
    }
    else
    {
        /* Call the library function */

        mlfMakeSbus(1, &sbus, var, Bus, Gen);
      size_t i = 0, j = 0; /* loop index variables */
      size_t r = 0, c = 0; /* variables to store the row and column 
      length of the matrix */
      double* data; /* variable to point to the double data stored 
      within the mxArray */

      /* Get the size of the matrix */
      r = mxGetM(sbus);
      c = mxGetN(sbus);
      /* Get a pointer to the double data in mxArray */
      data = mxGetPr(sbus);
       for (i = 0; i < c; i++)
       {
          for (j = 0; j < r; j++)
          {
            printf("%4.2f\t", data[j * c + i]);
          }
          printf("\n");
       }
      printf("\n");

      makeSbusTerminate();
   }

在我编译时从c调用函数，它显示“使用稀疏的错误 索引必须为正。

矩阵中的 我在这里采取的方法是正确的吗？

Answer 1

这是我在您的最新编辑中看到的问题：

稀疏数据的访问不正确（假设SBUS稀疏）：
稀疏矩阵仅物理存储非零元素。例如，如果矩阵大小为10x10，但只有7个元素是非零的，则物理上只有7个元素存储在内存中，而不是100个元素。如果您尝试使用两个嵌套的前面的元素访问这些元素，例如您正在通过代码数据[J * c + i]进行i和j范围基于矩阵大小的代码数据，则您将访问无效的内存和Crash Matlab。您必须具有特殊代码才能访问稀疏矩阵的元素，该矩阵考虑了存储在变量中的非零元素的索引。请参阅MXGetir和MXGetJC函数以访问这些索引。例如，此代码片段只需使用1个基于1的基于1的索引来打印数据：

mwIndex *Ir, *Jc, nrow, k=0; /* k will count how far we are into the sparse data */
Ir = mxGetIr(sbus); /* pointer to row indexes */
Jc = mxGetJc(sbus); /* pointer to cumulative number of elements up to this column */
for( j=0; j<n; j++ ) { /* for each column */
    nrow = Jc[j+1] - Jc[j];  /* number of non-zero row elements stored for this column */
    while(nrow--) { /* for each stored element in this column */
        printf("sbus(%d,%d) = %4.2f\n", Ir[k]+1,j+1,data[k]); /* 1-based index output */
        k++; /* increment counter */
    }
}

在上面，请注意，即使您在MATLAB级别将显示的索引显示为0时，稀疏索引数据也存储为基于0 1基。上面的代码在第一列，然后是第二列等中打印出元素。即，元素以列顺序存储在内存中，因此，当我们按内存顺序浏览数据时，它以列顺序自然而然地打印出来。

警告：上面的代码假定稀疏矩阵是真实的。如果稀疏矩阵很复杂，则需要考虑假想数据。访问真实＆amp;的方法假想数据取决于您是使用R2018A+交错复杂存储模型还是R2017B-独立的复杂存储器模型。您正在使用哪个版本的MATLAB以及使用哪种内存模型？这将影响如何编写元素访问的代码。对于R2018A+内存模型，它只是每个点而不是一个元素。例如，类似的东西：

    printf("sbus(%d,%d) = %4.2f + %4.2f i\n", Ir[k]+1,j+1,data[2*k],data[2*k+1]); /* 1-based index output */

矩阵大小的不匹配：
MxArrays Bus和Gen的创建为9x17和3x25，但您的C变量B和G被声明为[9] [13]和[3] [21]。他们不匹配。因此，在您的数据复制过程中，有没有设置的元素。

内存中元素顺序的不匹配：
MXARRAY矩阵元素在内存列中排序，而本机C矩阵元素在内存行中排序。从本质上讲，它们是在内存中彼此转置的。您不能仅仅是元素，因为订购是不同的。您需要以不同的方式读取数据，或者在执行数据副本时转换矩阵。

INT的不正确副本加倍：
data1是type int，但是mxarray var是type double。您无法通过memcpy将int位模式（data1）复制到双重（MxArray var的数据区域）中...结果将是垃圾。改用分配。例如，类似的事情会起作用：

*mxGetPr(var) = data1[0];

或者您可以直接从data1 [0]这样创建var：

var = mxCreateDoubleScalar(data1[0]);