C++ Thetafunktion Implementierung

Question

Ich versuche, diese Funktion zu implementieren:

aber es funktioniert nicht. Ein minimales, nachprüfbares Beispiel wie folgt aussieht:

#include <iostream> 
#include <cmath> 

int main() 
{ 
    int N {8};  // change this for testing <1..inf> 
    double q {0.1/N}; 
    int countN {static_cast<int>(floor(N/2))}; 
    static const double PI {3.1415926535897932384626433832795}; 

    // Omega[i] = Theta1(u,m)/Theta4(u,m) 
    double Omega[countN]; 
    for (int i=0; i<countN; ++i) 
    { 
    double micro {!(N % 2) * 0.5}; // 0 for odd N, 1/2 for even N 
    double num[countN] {sin(PI/N * (i + 1 - micro))}; 
    double den[countN] {0.5}; 

    for (int m=1; m<4; ++m) 
    { 
     num[i] += pow(-1, m) * pow(q, m*(m + 1)) * sin((2 * m + 1) * PI/N * (i + 1 - micro)); 
     den[i] += pow(-1, m) * pow(q, m*m) * cos(2 * m * PI/N * (i + 1 - micro)); 

    } 
    Omega[i] = fabs(pow(q, 0.25) * num[i]/den[i]); 
    std::cout << " " << Omega[i] << "\n"; 
    } 

    // testing the values, they should be increasing in value 
    for (const auto &elem: Omega) 
    std::cout << elem << " "; 
    std::cout << "\n"; 

    return 0; 
} 

Es gibt eine geringfügige Vereinfachung im Vergleich zum Original ist: I 2 in Zählern und Nennern berücksichtigt und ich nur die q^0.25 außerhalb der Fraktion. Auch countN ist die r aus dem Originaldokument, micro ist nur die 1/2 noch N oder 0 für ungerade N und i 0 für Array-Index aber i+1 für Berechnungen, aber diese sollten nicht Sache insgesamt.

Ich versuchte dies mit wxMaxima:

Theta[1](x,y):=2*y^0.25*sum((-1)^k*y^(k*(k+1))*sin((2*k+1)*x),k,0,n); 
Theta[4](x,y):=1+2*sum((-1)^k*y^(k^2)*cos(2*k*x),k,1,n); 
n:4$ 
N:8$ 
a:0.05$ 
b(i):=%pi/N*(i-(1-mod(N,2))/2)$ 
for N:8 thru 9 do for i:1 thru N/2 do print(["N=",N,"i=",i],Theta[1](b(i),a)/Theta[4](b(i),a)),numer; 

Und die Ergebnisse, in C++:

(q=0.05; N=8) 
Omega[0]=0.2018370065366672 
Omega[1]=0.06058232646142273 
Omega[2]=0.01205653570636574 
Omega[3]=0.02127667733703158 
(q=0.05; N=9) 
Omega[0]=0.348078726440638 
Omega[1]=0.1178366281313341 
Omega[2]=2.559808325080287e-07 
Omega[3]=0.02178788541277828 

und in wxMaxima:

["N=",8,"i=",1]" "0.2018370065366672" " 
["N=",8,"i=",2]" "0.5439269564954693" " 
["N=",8,"i=",3]" "0.7569342043740249" " 
["N=",8,"i=",4]" "0.850913653939989" " 
["N=",9,"i=",1]" "0.348078726440638" " 
["N=",9,"i=",2]" "0.6165773889432575" " 
["N=",9,"i=",3]" "0.7800391631077094" " 
["N=",9,"i=",4]" "0.8532352152763631 

Zu meiner Überraschung, der erste Term gut, für bith N, so kann ich nicht sagen, was in meinem Code nicht richtig ist. Könnte mir jemand helfen, den Fehler zu finden?

Um es klar zu machen: Ich bin ein Anfänger in C++ und ich suche nicht jemanden, der es für mich tut, aber mich von meinem Fehler in der Programmierung (übersetzen Mathe zu C++ - Code) zu lassen.

Answer 1

Sie haben

double den[countN] {0.5}; 

dies initialisiert das erste Element von den-0.5und alle anderen Elemente zu 0.0 (Standard-Initialisierung). Mit anderen Worten, das Obige ist mit so vielen Nullen wie notwendig zum Füllen des Arrays äquivalent. Wahrscheinlich wollten Sie alle Elemente zu 0.5 initialisieren. In Ihrem Fall ist der einfachste Weg, das zu tun, wenn Sie zuerst das Element verwenden - oder, da Sie nur das einzelne Element zugreifen den[i] während der Lebensdauer der den, es eine einfaches double eher als ein Array machen:

for (int i=0; i<countN; ++i) { 
    double micro {N % 2 ? 0.0 : 0.5}; // 0 for odd N, 1/2 for even N 
    double num{sin(PI/N * (i + 1 - micro))}; 
    double den{0.5}; 

    for (int m=1; m<4; ++m) { 
     num += pow(-1, m) * pow(q, m*(m + 1)) * sin((2 * m + 1) * PI/N * (i + 1 - micro)); 
     den += pow(-1, m) * pow(q, m*m) * cos(2 * m * PI/N * (i + 1 - micro)); 
    } 
    Omega[i] = fabs(pow(q, 0.25) * num/den); 
}