Precisione di Math.Sin() e Math.Cos() in C#

Question

Sono terribilmente infastidito dall'inaccuratezza delle funzioni trigonometriche intrinseche nel CLR. E 'ben noto che

Math.Sin(Math.PI)=0.00000000000000012246063538223773 

invece di 0. Qualcosa di simile accade con Math.Cos(Math.PI/2).

Ma quando sto facendo una lunga serie di calcoli che in casi particolari restituiscono

Math.Sin(Math.PI/2+x)-Math.Cos(x) 

e il risultato è zero per x = 0.2, ma non zero per x = 0.1 (provare per credere). Un altro problema è quando l'argomento è un numero elevato, l'inaccuratezza diventa proporzionalmente grande.

Quindi mi chiedo se qualcuno ha codificato una migliore rappresentazione delle funzioni trigonometriche in C# per la condivisione con il mondo. Il CLR chiama alcune librerie matematiche C standard che implementano CORDIC o qualcosa di simile? link: wikipedia CORDIC

Answer 1

È necessario utilizzare una libreria decimale di precisione arbitraria. (.Net 4.0 ha un arbitrary integer class, ma non decimale).

Pochi quelli popolari sono disponibili:

• BigNum
• W3B.Sine

Answer 2

Questo non ha nulla a che fare con l'accuratezza delle funzioni trigonometriche ma più con il sistema di tipo CLS. Secondo la documentazione uno double ha una precisione di 15-16 cifre (che è esattamente ciò che si ottiene) quindi non si può essere più precisi con questo tipo. Quindi se vuoi più precisione dovrai creare un nuovo tipo che sia in grado di memorizzarlo.

Si noti inoltre che non si dovrebbe mai essere la scrittura di un codice come questo:

double d = CalcFromSomewhere(); 
if (d == 0) 
{ 
    DoSomething(); 
} 

si dovrebbe fare invece:

double d = CalcFromSomewhere(); 
double epsilon = 1e-5; // define the precision you are working with 
if (Math.Abs(d) < epsilon) 
{ 
    DoSomething(); 
} 

Answer 3

Questo è il risultato di precisione in virgola mobile. Si ottiene un certo numero di cifre significative possibili e tutto ciò che non può essere rappresentato esattamente viene approssimato. Ad esempio, pi non è un numero razionale, quindi è impossibile ottenere una rappresentazione esatta. Dal momento che non è possibile ottenere un valore esatto di pi, non si otterranno esatti seni e coseni di numeri compreso pi (né si otterranno i valori esatti di seni e coseni per la maggior parte del tempo).

La migliore spiegazione intermedia è "What Every Computer Scientist Should Know About Floating-Point Arithmetic". Se non vuoi approfondire, ricorda che i numeri in virgola mobile sono in genere approssimazioni e che i calcoli in virgola mobile sono come muovere cumuli di sabbia sul terreno: con tutto quello che fai con loro, perdi un po 'di sabbia e raccogliere un po 'di sporcizia.

Se si desidera una rappresentazione esatta, è necessario trovare un sistema di algebra simbolico.

Answer 4

Ho sentito. Sono terribilmente infastidito dall'inesattezza della divisione. L'altro giorno ho fatto:

Console.WriteLine(1.0/3.0); 

e ho avuto ,333333333333333, invece che la risposta corretta, che è 0.3333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333 3333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333 3333333333333333 ...

Forse ora si vede qual è il problema. Math.Pi non è uguale a pi più di 1.0/3.0 è uguale a un terzo. Entrambi differiscono dal vero valore di qualche centinaio di quadrilioni, e quindi tutti i calcoli che esegui con Math.Pi o 1.0/3.0 saranno anch'essi spenti di qualche centinaio di quadrilioni, compreso il seno.

Se non vi piace quell'aritmetica approssimativa è approssimativa a quindi non utilizzare l'aritmetica approssimativa. Usa l'aritmetica esatta. Usavo Waterloo Maple quando avevo bisogno dell'aritmetica esatta; forse dovresti comprarne una copia.

Answer 5

mi respingono l'idea gli errori sono dovuti a round-off. Cosa si può fare è definire sin(x) come segue, con l'espansione di Taylor con 6 termini:

const double π=Math.PI; 
    const double π2=Math.PI/2; 
    const double π4=Math.PI/4; 

    public static double Sin(double x) 
    { 

     if (x==0) { return 0; } 
     if (x<0) { return -Sin(-x); } 
     if (x>π) { return -Sin(x-π); } 
     if (x>π4) { return Cos(π2-x); } 

     double x2=x*x; 

     return x*(x2/6*(x2/20*(x2/42*(x2/72*(x2/110*(x2/156-1)+1)-1)+1)-1)+1); 
    } 

    public static double Cos(double x) 
    { 
     if (x==0) { return 1; } 
     if (x<0) { return Cos(-x); } 
     if (x>π) { return -Cos(x-π); } 
     if (x>π4) { return Sin(π2-x); } 

     double x2=x*x; 

     return x2/2*(x2/12*(x2/30*(x2/56*(x2/90*(x2/132-1)+1)-1)+1)-1)+1; 
    } 

errore tipico è 1e-16 e caso peggiore è 1e-11. È peggio del CLR, ma è controllabile aggiungendo più termini. La buona notizia è che per i casi speciali nell'OP e per Sin(45°) la risposta è esatta.