2017-07-06 7 views
1

私は比較的新しいC++ですが、私はいくつかの(不十分な)コーディングと数値的な経験があります。std :: valarray <double>とgslを統合するには?

私はこの質問が毎回投稿されることを知っていますし、配列をどのように統合しますか? MATLABでは、あなたの配列を関数にすることができます(私は忘れてしまったが、前にやったことは分かっている)。そして組み込みのインテグレータにそれを送る。だから私の質問はC++でどうやってやるかだ。

私はこの整数があります。

I = integral(A(z)*sin(qz)*dz) 

qはちょうどダブルconstので、Zは、積分変数ですが、(z)は、同じ持つ配列(私は今からactualfunctionそれを呼ぶことにします)です私のコードではz軸としての点の数。積分境界はz [0]とz [nz-1]です。

この積分は台形則を用いて計算し、z軸の5000点については0.06秒かかる。私の問題は、この計算はおよそ300 * 30 * 20回(私は3つのループがある)であり、この0.06秒はシミュレーションの3時間に非常に速く成長するということです。私のコードの全体的なボトルネックはこの統合です(私は明らかにzを減らすことでスピードアップすることができますが、それはポイントではありません)。

私は、ライブラリ関数は通常、また、シンプソンのルールとして単純なものを使うこともできないことも知っています。なぜなら、被積分関数は非常に振動的であり、複雑な数値的忠実度を自分で実装したくないからです。

GSLは、フォーム内の関数を必要とします:

F F(ダブルのx、void *型のparams)

を=と私はおそらくGSLからQAWO適応統合を使用することができますが、どのように私は形で私の機能を作るのですか私の配列を関数に変換しますか?

私のように何かを考えています:

F(double z, void *params) 
{ 
    std::valarray<double> actualfunction = *(std::valarray<double> *) params; 
    double dz = *(double *) params; // Pretty sure this is wrong 
    unsigned int actual_index = z/dz; // crazy assumption (my z[0] was 0) 
    return actualfunction[actual_index]; 
} 

は、この可能性のようなものですか?私は数値アルゴリズムが実際の関数と同じ空間的な差異を使用するのではないかと疑問に思っています。

gssの方がいいですか?

答えて

0
template<class F> 
struct func_ptr_helper { 
    F f; 
    void* pvoid(){ return std::addressof(f); } 
    template<class R, class...Args> 
    using before_ptr=R(*)(void*,Args...); 
    template<class R, class...Args> 
    using after_ptr=R(*)(Args...,void*); 
    template<class R, class...Args> 
    static before_ptr<R,Args...> before_func() { 
    return [](void* p, Args...args)->R{ 
     return (*static_cast<F*>(p))(std::forward<Args>(args)...); 
    }; 
    } 
    template<class R, class...Args> 
    static after_ptr<R,Args...> after_func() { 
    return [](Args...args, void* p)->R{ 
     return (*static_cast<F*>(p))(std::forward<Args>(args)...); 
    }; 
    } 
}; 
template<class F> 
func_ptr_helper<F> lambda_to_pfunc(F f){ return {std::move(f)}; } 

使用:その後、

auto f = lambda_to_pfunc([&actualfunction, &dz](double z){ 
    unsigned int actual_index = z/dz; // crazy assumption (my z[0] was 0) 
    return actualfunction[actual_index]; 
}); 

void* pvoid - f.pvoid(); 
void(*pfun)(double, void*) = f.after_func(); 

、あなたが通過pfunpvoidを渡すことができます。

すべてのタイプミスを謝罪します。

アイデアは私たちが望むようにラムダを書くことです。それでlambda_to_pfuncはそれを包みますので、void*とそれをCスタイルAPIへの関数ポインタとして渡すことができます。

関連する問題