Cos(x) - szereg Taylora nie działa dla x>2

0
 int factorial(int n) {
    int i;
    int product = 1;
    for (i = 2; i <= n; ++i) product *= i;
    return product;
}

double cosx_Taylor(double x) {
	const double eps = 0.00001; // dokładność
	
	double sum = 0;
	int i = 0;
	double sign = 1.0;
	double w = 1.0;
	
	while (eps < fabs(w)) {
          w = (pow(x, (double) (2 * i)) * sign) / ((double) (factorial(2 * i)));
          sign = -sign;
          sum += w;
          ++i;
    } 
    return sum;
}

Proszę was o pomoc, moja funkcja nie działa już dla x minimalnie większych od 2, w czym tkwi problem? Myślę, że algorytm jest poprawny.

2

Włącz myślenie.

  1. Kolejny wyraz szeregu JEST mocno skorelowany z poprzednim wyrazem. Nie ma sensu liczyć znów potęgi i silni, wystarczy pomnożyć i podzielić przez coś...
  2. Za duże liczby. Zwykły int nie pomieści zbyt wysokiej silni (long long da radę maksymalnie do 20!)
0
double cosx_Taylor(double x) {
	const double eps = 0.00001; // dokładność
	
	double i = 0;
	double l = 1.0; // licznik
	double m = 1.0; // mianownik
	double sign = -1.0;
	double sum = 1.0;
	double w = 1.0; // bieżący wyraz
	double x2 = x*x;
	
	while (eps < fabs(w)) {
          l *= x2;
          m *= (i+1)*(i+2);
          w = (sign * l)/m;
          sign = -sign;
          sum += w;
          ++i;
    } 
    return sum;
} 

Dzięki Twojej radzie, program działa tyle, że teraz wyświetla błędny wynik, zrezygnowałem z osobnej funkcji liczącej silnie. Nie wiem co tym razem jest źle :<

0
#include <stdio.h>
#include <math.h>


double Cos(double x) {
	const double eps = 0.00001; // dokładność

	double i = 0;
	double l = 1.0; // licznik
	double m = 1.0; // mianownik
	double sign = -1.0;
	double sum = 1.0;
	double w = 1.0; // bieżący wyraz
	double x2 = x*x;

	while (eps < fabs(w)) {
		l *= x2;
		m *= (i+1)*(i+2);
		w = (sign * l)/m;
		sign = -sign;
		sum += w;
		i += 2; // Tu.
	} 
	return sum;
} 



int main(int argc, char* argv[]) {
	double x;
	for (x = 0.0; x < 3.5; x += 0.1)
		printf("Cos(%.1lf) - cos(%.1lf) = %lf\n",
			x, x, Cos(x)-cos(x) );
	
	return 0;
}

Można jeszcze zwiększyć troszkę wydajność:

double Cos(double x) {
	const double eps = 0.00001;
	const double x2 = x*x;

	int i = 1;
	double sum = 1.0;
	double w = 1.0;

	while (eps < fabs(w)) {
		w *= -x2/(i*(i+1));
		sum += w;
		i += 2;
	} 

	return sum;
}
0

Dzięki wielkie :) Jeszcze tylko jedno pytanie, czy to normalne, że dla większych wartości x np. dla 200, wartość funkcji jest kompletnie różna od prawidłowej?

0

Dla większych argumentów nasza suma szybko zmierza do nieskończoności, a potem trudno ją z tej nieskończoności sprowadzić na ziemię. ;) Typ double ma ograniczoną dokładność, co objawia się zaniedbywaniem małych wartości dodanych do dużych. Przeprowadź test i zauważ co się dzieje pod końcem.

#include <stdio.h>
#include <math.h>

double Cos(double x) {
	const double eps = 0.00001;
	const double x2 = x*x;

	int i = 1;
	double sum = 1.0;
	double w = 1.0;


	puts("  i:                  sum += w");
	puts("-----------------------------------------------------");
	
	while (eps < fabs(w)) {
		w *= -x2/(i*(i+1));
		
		printf("%3d: %20g += %g\n", i, sum, w);
		sum += w;
		
		i += 2;
	}
	
	puts("-----------------------------------------------------");
	printf("Cos(%g) = %g\n", x, sum);
	
	return sum;
}


int main(int argc, char* argv[]) {
	Cos(150);
	
	return 0;
}

Oczywiście można tej sytuacji zaradzić. Podam dwa najbanalniejsze sposoby:

Sposób 1 (pamięciożerny):
Każdy wynik cząstkowy zapisujemy do tablicy, którą potem sortujemy i dodajemy wartości w odpowiedniej kolejności (tak aby nie tracić precyzji).

Sposób 2:
Możemy najpierw wybrać odpowiednio mały x, dla którego wartość cosinusa jest taka sama jak dla dużego x. ;)

0
#define PI 3.1415926535897932384626433832795
main(void) {
	double f, i, x=2000.0;

	f = 2 * PI * modf( x/(2*PI), &i );     // |f| < 2*Pi
	printf("%f == %f", cos(x), cos(f));} 
0

po redukcji zakresu można jeszcze skorzystać ze wzoru cos(2x) = 2 cos(x)^2 - 1

double COS(double x){
        if( fabs(x) >= 1) {
                double t=COS(x/2);
                return 2*t*t - 1; }
        else return cos(x); } 

albo iteracyjnie np.:

double Cos(double x){
        int i=0; double t;
        while( fabs(x) > 1 ) {
                x/=2; i++; }
        t=cos(x);
        while(i>0){
                t=2*t*t-1; i--;}
        return t;} 

razem

	double nic;
	x = 2 * PI * modf( x/(2*PI), &nic );     // |x| < 2*Pi
	x /= 8;
	x = cos(x);                              // |x| < 0.785
	x = 2*x*x-1;
	x = 2*x*x-1;
	x = 2*x*x-1;
 
0

Dziękuję wszystkim za odpowiedzi i pozdrawiam

1 użytkowników online, w tym zalogowanych: 0, gości: 1