2015-09-10 5 views
5

Ich hoffe, diese Frage (und mögliche Antworten darauf) ist allgemein genug im Umfang, so dass es auch für andere nützlich sein wird.Schreiben einer Vorlage für variierende Funktionen mit verschobenen Argumenten

Ich versuche, ein numerisches Problem zu lösen, das in Form

wo ist eine vordefinierte Funktion ein Produkt von double Funktionen beinhaltet, und i passiere einen Integrator.

Die Komplikation ist, dass meine Funktionen nicht statisch sind; Die Integration wird wiederholt durchgeführt, und bei jeder Integration des s hat eine andere Form, beispielsweise in der ersten Schleife könnte es

die zweiten sein, es wird

und so weiter, wo die Mitglieder eines Vektors von Funktionen bezeichnen. Bei jeder Iteration der Integration weiß ich a priori nicht, was die Form von sein wird, nur dass es eine lineare Kombination von s mit konstanten Koeffizienten sein wird.

Da , die ich an den Integrator passieren beinhaltet mit „verschoben“ Argumente, bei jeder Iteration i zugreifen müssen zu jedem der , da ich etwas von der Form

double Y (double x, double y, double z){ 

    return (f(x+y) - f(x)) * (f(y+z) - f(y)) * A(z+y); 

} 

für den Integrator anrufen müssen (oder vielleicht nicht?).

Meine Frage ist, was wäre der beste Weg, meine darzustellen? Ich habe die offensichtlichen Kandidaten (Templates, std :: function) gelesen, aber ich bin mir nicht sicher, wie ich die verschobenen Argumente auf die Funktion "aufprägen" kann.

Ich bin mir sicher, dass dies in einer anderen Inkarnation ein ziemlich typisches Programmierdesignproblem ist. Kann mir jemand Vorschläge machen, wie ich mich ihm nähern soll?

Bearbeiten: Dank Rostislav und vor allem Gombat, die die ursprüngliche Lösung (obwohl Rostislav ist ziemlich elegant) kam. Bevor ich diese Frage schließe, füge ich hier vielleicht eine weitere Frage hinzu, dh betrachte den (allgemeineren) Fall, in dem man nicht so detaillierte Informationen über die Ausgabe f(x) hat (dass es eine lineare Kombination von f_i(x) mit konstanten Koeffizienten ist) Jede Iteration, und stattdessen wissen wir nur, dass f(x) = F(x), F(x) ist eine double Funktion, die wieder aus der f_i(x) zusammengesetzt ist, aber jetzt in beliebiger Weise. Alle anderen Bedingungen sind die gleichen (beginnen mit bekannten f_0(x)), wie sollte ich meine F(x) in diesem Fall darstellen?

Ich dachte, dass da wir die funktionale Abhängigkeit der Lösung F(x) wissen (es hängt nur von einem Argument x, ist ein double und besteht aus bekannten f_i(x) zusammengesetzt), vielleicht ist es möglich, eine „Vorlage“ (schreiben Entschuldigung für den Missbrauch der Nomenklatur, es ist sehr suggestiv) Funktion F(x), die ich an den Integrator übergeben kann, die dann die tatsächliche Form der Funktion "integriert", um integriert zu werden, obwohl dies kompliziert sein könnte.

+0

Wie bestimmen Sie, welche Form f (x) zur Laufzeit haben soll? – Gombat

+0

@Gombat In der 1. Iteration stelle ich einen Vektor von 'f_0 (x)' Funktionen zur Verfügung, der bekannt ist.Dieser Vektor wird mit einer bestimmten Matrix multipliziert, die dazu dient, sie zu mischen, und daher habe ich zu Beginn der 2. Iteration einen Vektor von Linearkombinationen (mit konstanten Vorfaktoren) der 'f_0 (x)' s, die ich füttere elementweise zur 'Y (x, y, z)' Funktion (dh in diesem Schritt sollte jede der 'f (x)' Funktionen in 'Y (x, y, z)' bereits eine lineare Kombination sein der 'f_0 (x)' s). –

+0

Was ist mit dem Erstellen einer Klasse, die ein Mitglied hat, das die Informationen speichert, die benötigt werden, um 'f_i (x)' s zu kombinieren und 'f (x)' eine Methode der Klasse zu sein. 'f (x)' hätte Zugang zu dem Mitglied, das die benötigten Informationen hält. – Gombat

Antwort

1

auf die Kommentare Basierend hoffe ich folgendes helfen:

class MyF 
{ 
public: 
    // typedef as abbreviation 
    typedef std::function<double(const double)> func; 

    MyF(const std::vector<func>& functions) 
    : m_functions(functions) 
#ifdef _DEBUG 
    , m_lcUpdated(false) 
#endif 
    { } 

    // f combines the inner functions depending on m_linearCombination 
    double operator()(const double x) 
    { 
    assert(m_lcUpdated); 
    assert(m_linearCombination.size() > 0); 
    assert(m_linearCombination.size() <= m_functions.size()); 

    double ret(0.0); 
    // if m_linearCombination has a value only for the first function 
    // no need to run the loop for all functions, so we iterate only 
    // over the linearCombination entries. 
    for (size_t i = 0; i < m_linearCombination.size(); i++) 
    { 
     // if the factor is very small, no need for m_function call 
     if (std::abs(m_linearCombination) < DBL_EPSILON) continue; 

     ret += m_linearCombination[i] * m_functions[i](x); 
    } 

    m_lcUpdated = false; 
    return ret; 
    } 

    void setLinearCombination(const std::vector<double>& lc) 
    { 
     m_linearCombination = lc; 
#ifdef _DEBUG 
     m_lcUpdated = true; 
#endif 
    } 


private: 
    std::vector<double> m_linearCombination; 
    std::vector<func> m_functions; 
#ifdef _DEBUG 
    bool m_lcUpdated; 
#endif 
}; 

Der Konstruktor von MyF bekommt die Funktionen f_i(x) in der ersten Iteration (oder vor Sie auch MyF vor, ohne Argumente bauen könnte und ein schreiben? Methoden einstellen, um die vector von std::function s) einzustellen. Das private Mitglied m_lcUpdated überprüft, ob Sie die lineare Kombination vor dem Aufruf der Funktion mit operator() aktualisiert haben.

std::vector<MyF::func> functions; 

// ...fill the functions 

// construct your function class 
MyF F(functions); 
// and the variable to fill in each iteration on how to linear combine the inner functions 
std::vector<double> linearCombination; 

// Provide the first linear combination (fill it before of course) 
F.setLinearCombination(linearCombination); 

// In each iteration, call the method combining inner functions 
F(x); // it works like this, since it is the operator() which is overloaded for const double as input 
// And provide the linear combination for the next iteration 
F.setLinearCombination(linearCombination); 

Es wäre auch möglich, die linearCombinations als std::map<size_t,double> zu haben, wo zuerst der Index ist, so dass Sie sich über die Elemente nicht wiederholen, deren Wert möglicherweise Null, wenn es möglich ist, dass eine lineare Kombination z.B f_0(x) + f_2(x).

Bearbeiten: Ich implementiert @ Gombat-Lösung, wie es intuitiver war, aber auf jeden Fall auch Rostislavs Antwort unten sehen.

3

Da Sie Y integrieren, ist es wahrscheinlich viele Male bei der Integration aufgerufen werden, so ist es besser, im Grunde nicht die Anrufe aufgrund Typs Löschung zu std::function Inline kann mit std::function wie den Compiler zu vermeiden, dass std::function intern verwendet.

Stattdessen können Sie die Lambdas direkt verwenden. Wenn Sie Kontrolle über die Implementierung der Integrationsprozedur haben (oder die Art der Funktion bereits als Template-Parameter integriert haben), können Sie alle std::function Verwendungen vermeiden. Andernfalls wäre es sinnvoll, std::function nur für Ihre Kommunikation mit dem Integrationsverfahren zu verwenden, aber vermeiden Sie sie wo möglich. Unten ist der Beispielcode es (live Beispiel here) Umsetzung:

#include <iostream> 
#include <string> 

#include <tuple> 
#include <functional> 

template<typename... Fn> 
class LinearCombination { 
public: 
    template<typename... Un> 
    LinearCombination(Un&&... fs) 
     : functions(std::forward<Un>(fs)...) 
    { 
     coefs.fill(0.0); 
    } 

    double operator()(double x) const 
    { 
     return evaluateImpl(x, std::integral_constant<size_t, sizeof...(Fn)-1>()); 
    } 

    void setCoef(size_t i, double c) 
    { 
     coefs[i] = c; 
    } 

private: 
    template<size_t I> 
    double evaluateImpl(double x, std::integral_constant<size_t, I> index) const 
    { 
     return evaluateOne(x, index) + evaluateImpl<I - 1>(x, std::integral_constant<size_t, I - 1>()); 
    } 

    template<size_t I> 
    double evaluateImpl(double x, std::integral_constant<size_t, 0> index) const 
    { 
     return evaluateOne(x, index); 
    } 

    template<size_t I> 
    double evaluateOne(double x, std::integral_constant<size_t, I>) const 
    { 
     auto coef = coefs[I]; 
     return coef == 0.0 ? 0.0 : coef * std::get<I>(functions)(x); 
    } 


    std::tuple<Fn...> functions; 
    std::array<double, sizeof...(Fn)> coefs; 
}; 

// This helper function is there just to avoid writing something like... 
// LinearCombination<decltype(f1), decltype(f2)> when you use lambdas f1 and f2 
template<typename... Fn> 
auto make_linear_combination(Fn&&... fn) 
{ 
    return LinearCombination<Fn...>{std::forward<Fn>(fn)...}; 
} 

double A(double) 
{ 
    return 1.0; 
} 

/// Integration 1.0 
double integrate3D(double(*f)(double, double, double)) 
{ 
    return f(1, 2, 3); 
} 

struct YHelper { 
    static double Y(double x, double y, double z) 
    { 
     return (f(x + y) - f(x)) * (f(y + z) - f(y)) * A(z + y); 
    } 

    static std::function<double(double)> f; 
}; 

std::function<double(double)> YHelper::f; 


/// Integration 2.0 
template<typename Integrable> 
double integrate3D_2(Integrable&& f) 
{ 
    return f(1, 2, 3); 
} 

int main() 
{ 
    auto f1 = [](double x) { return x; }; 
    auto f2 = [](double x) { return 2 * x; }; 

    auto lc = make_linear_combination(std::move(f1), std::move(f2)); 
    lc.setCoef(0, 1.0); 
    lc.setCoef(1, -1.0); 

    std::cout << lc(2.0) << "\n"; 

    YHelper::f = std::ref(lc); 
    std::cout << integrate3D(&YHelper::Y) << "\n"; 

    auto Y = [&lc](double x, double y, double z) { return (lc(x + y) - lc(x)) * (lc(y + z) - lc(y)) * A(z + y); }; 
    std::cout << integrate3D_2(Y) << "\n"; 
} 

Hinweis, dass unter Integration 1.0 ist eine Beispielimplementierung für den Fall, dass Sie keine Kontrolle über die Unterzeichnung des Integrationsverfahrens haben. Der Code unter Integration 2.0 ist nicht nur sauberer, sondern wird auch besser funktionieren.

PS: Natürlich vorsichtig sein, wenn coef zu 0.0 in evaluateOne vergleichen - davon ausgehen, es wird nur funktionieren, wenn Sie die Koeffizienten direkt auf wörtliche 0.0 und keinem Ergebnis einer Berechnung festgelegt. Andernfalls verwenden Sie abs(coef) < epsilon mit einem epsilon Wert, der zu Ihrer Anwendung passt.

Verwandte Themen