1

私は最近、MATLAB で有限要素法に取り組んでいます。

MATLAB でコードを最適化しようとしました。

検索中に、mex 関数を使用すると行列の組み立てが高速化できることがわかりました。

「PoissonAssembler.m」を mex 関数に移植しているときに、いくつかの問題が発生しました。

PoissonAssembler.m は基本的にこのような構造になっています。

for every element{

    loc_stiff_assmbl();

    loc_rhs_assmbl(); // this one involves with function evaluation @(x,y)

    for every edges{
        loc_edge_assmbl();

        if boundary{
            loc_rhs_edge_assmbl(); // this one also have function eval
        }
    }
}

matlabでは、私は持っています

f = @(x,y) some_math_function;

この関数は他の数値シミュレーション用に変更されるため、

関数ハンドルを mex ファイルの入力として使用したかった

mexCallMatlab() と feval を使用してこれを行う方法があることがわかりました。

ただし、matlab の呼び出しによるオーバーヘッドが原因で、コードの速度が低下すると思います。

関数ハンドルを変更するたびに mex ファイルをコンパイルせずに回避する方法はありますか?

より正確なコードは次のようになります

for (k=0;k<nE;k++){ // nE : number of elements (about 10^5 to 10^6)
    // other assembling procedure
    blabla

    // function evaluation for rhs assemble
    for (i=0; i<nP; i++){ // nP : number of points in element
        fxy[i] = f(x[i],y[i])};
    }

    // rhs assemble
    for (j=0; j<nP; j++){
        for (i=0; i<nP; i++){ // matrix vector multiplication
            rhs[k*nE + j] += Mass[ i*nP + j]*fxy[i];
        }
    }

    // face loop
    for (f1 = 0; f1 < nF4E; f1++){ // number of face for each elements
        switch(BCType[k1*nF4E + f1]){ // if boundary

            case 1 : // Dirichlet
                // inner product with evaluation of gD = @(x,y) function
                CODE OMITTED
            case 2 : // Neumann
                // inner product with evaluation of gN = @(x,y) function
                CODE OMITTED
        }
    }
}
4

1 に答える 1

1

C (および C++) は、MATLAB 関数ハンドルと同様に、関数の概念を変数およびパラメーターとしてサポートします。サポートされる関数を合理的なセットに制限できる場合は、MEX 関数に渡される追加の引数を使用して関数を選択できます。擬似コードは次のとおりです。

double fn_linear(double x, double y, double *args)
{
    return arg[0] + x*arg[1] + y*arg[2];
}

double fn_quadratic(double x, double y, double *args)
{
    return arg[0] + x*arg[1] + x*x*arg[2] + /// etc
}

typedef double (*EdgeFunction)(double x, double y, double *args);


void computation_function(other parms, EdgeFunction edge_fcn, double *fcn_parms)
{
    for (k=0;k<nE;k++){ // nE : number of elements (about 10^5 to 10^6)
        // other assembling procedure
        blabla

        // function evaluation for rhs assemble
        for (i=0; i<nP; i++){ // nP : number of points in element
            fxy[i] = (*edge_fcn)(x[i], y[i], fcn_parms);
        }
}

void mexFunction()
{
    EdgeFunction edge_fcn;

    if(strcmp(arg3, "linear")==0) {
       edge_fcn = &fn_linear;
       extra_parms = arg4;
    } else {
      ...
    }
}

このように、エッジ関数を選択する文字列と、関数が必要とするいくつかのパラメーターを使用して関数を呼び出します。

std::function自分の状況で C++、特に C++11 を使用できる場合、 andと lambda を使用すると、この種の作業がはるかに簡単にbindなります。可能であれば、それを理解する価値があります。この場合std::map、簡単に検索できるように、文字列名を関数に定義することもできます。

最後に、これらの実装のいずれかを使用すると、既に把握している MATLAB コールバック用に追加の関数を 1 つ定義できることを指摘しておきます。手書きの選択肢では速く、本当に一般的なケースでは遅くなります。両方の長所。

于 2016-07-08T19:57:54.963 に答える