[C con Clase] Operaciones trigonométricas y C++

[Davidson, Steven]

Hola Santiago,

2012/10/8 Santiago Tabarez <santiago230792 en gmail.com>:
>
>

[CORTE]

>
> Tu implementacion me gustóo y luce más limpia. Si no te importa la voy a
> usar.
>

No; no me importa. Tampoco es que sea "mi implementación" :)

>

[CORTE]

>
> Claro, me gustaría saber cómo implementar este otro algoritmo. Hasta tal vez
> podría resultar más eficiente.
>

Sí; así es. El algoritmo trivial es de complejidad O( n ), pero el que
te voy a explicar a continuación es de O( lg n ); por cierto, "lg" es
el logaritmo en base 2.

Básicamente, vamos acumulando cuadrados de la base, en lugar de
acumular la base. Esto significa que tratamos el exponente como una
potencia de 2. Por ejemplo, 5^13.

5^13 = 5^( 8 + 4 + 1 ) = 5^( 2^3 + 2^2 + 2^0 )

Como los exponentes son en base 2, podemos usar operaciones a nivel de
bits para averiguar cuáles exponentes de base 2 afectan el resultado
de la potenciación. Esto sería,

5^( 2^3 + 2^2 + 2^0 ) = 5^( 1 * 2^3 + 1 * 2^2 + 0 * 2^1 + 1 * 2^0 ) =
5^2^3 * 5^2^2 * 5^(0*2^1) * 5^2^0;

Factorizando, obtenemos,

5^2*5^2*5^2 * 5^2*5^2 * 1 * 5

De un modo iterativo, podríamos escribir esto:

base = 5*5;  // "nueva" base

acum = 5;  // 5^2^0
res = acum;  // 5^2^0

acum *= base;  // 5^2^1 * 5^2^0

acum *= base;  // 5^2^2 * 5^2^1 * 5^2^0
res *= acum;  // 5^2^2 * 5^2^1 * 5^2^0

acum *= base;  // 5^2^3 * 5^2^2 * 5^2^1 * 5^2^0
res *= acum;  // 5^2^3 * 5^2^1 * 5^2^0

Al final, 'res' guarda el resultado de 1220703125. Como puedes ver, no
siempre acumulamos lo que llevamos hecho en el resultado, porque
algunos de esos exponentes están "desactivados". Esto lo podemos ver
mejor si convertimos los exponentes a binario, que sería la secuencia:
1101b. Cada vez que veamos un 1, acumulamos lo que tenemos al
resultado, y cada vez que veamos un 0, no se acumula al resultado; en
realidad deberíamos acumular un 1, pero para eso no hacemos el
"esfuerzo" y simplemente nos saltamos ese paso.

El algoritmo general es el siguiente:

1.  acum <- base

// Caso particular: exp = 0
2.  Si exp AND 1 = 1, entonces
3.      res <- acum
4.  Si no, entonces
5.      res <- 1
6.  base <- base * base
7.  exp <- exp SHR 1

// Caso general: exp > 0
8.  Mientras que, exp > 0, repetir
9.      acum <- acum * base
10.    Si exp AND 1 = 1, entonces
11.        res <- res * acum
12.    exp <- exp SHR 1
13. Terminar( res )

Aquí, AND es la operación AND a nivel de bits, mientras que SHR es el
desplazamiento de bits a la derecha (SHift Right). Los operadores en
C++ son: & y >>, respectivamente.

>

[CORTE]

> Pero hay un caso que es cero elevado a la cero potencia ( 0º ) que en
> principio sería indeterminado.
>

Efectivamente. Puedes representar tal resultado con NaN (no es un
número) que viene a ser lo mismo que "indeterminado". En C++, puedes
usar la estructura-plantilla, 'numeric_limits', para invocar la
función miembro estática, 'quiet_NaN()'. Esto es,

#include <limits>
...
double num = numeric_limits<double> :: quiet_NaN();

Claro que sería prudente comprobar que la implementación soporta tal
representación, usando el dato miembro estático, 'has_quiet_NaN'.

Por otro lado, podrías forzar tal representación al intentar hacer esto:

return 0.0 / 0.0;

Así dejas al compilador y enlazador que se encarguen de los detalles.

>

[CORTE]

> Actualmente estoy intentando ver cómo implementar el método de
> Newton-Raphson, que lo veo muy complicado. Por ahora voy a dejar las raíces
> de lado un tiempo.
>

En la versión inglesa del enlace que te di habla justamente de este
problema como un ejemplo. Consulta este enlace:
http://en.wikipedia.org/wiki/Newton%27s_method#Square_root_of_a_number

Básicamente, este método se basa en una función, f(x) = 0, al igual
que su derivada. Por lo tanto, convertimos el problema que tenemos:

x = radical( a ),  donde "a" es un número dado

a un problema que el método N-R "entiende"; esto es,

x^2 = a;
x^2 - a = 0;

Por consiguiente,

f(x) = x^2 - a

y su derivada,

f'(x) = 2x

Ahora es cuestión de implementar unos cuantos términos basados en este método:

x0 : número dado, para que sea la aproximación inicial.

x1 = x0 - f(x0) / f'(x0) = x0 - (x0^2-a) / (2x0)
x2 = x1 - f(x1) / f'(x1)
x3 = x2 - f(x2) / f'(x2)
...

Típicamente, aplicaríamos un nivel de tolerancia para detener este
método, para asegurarnos una precisión deseada.

>

[CORTE]

>
> Además creé una función que no esperaba crear, los factoriales. Utilicé
> valores de tipo long double dada la posibilidad de que los resultados puedan
> tomar valores muy grandes (el máximo valor al que devuelve un resultado es
> 170, luego de eso devuelve un "1.#INF"). Si hay algo mal avisame.
>
> long double factorial (long double x)
>     // cálculo de factoriales
> {
>     long double y = 1;
>
>     while ( x > 1 ) {
>         y *= x;
>         x--;
>     }
>     return y;
> }
>

No veo nada raro con este algoritmo. Eso sí, ten cuidado con 'long
double', ya que no todas las implementaciones realmente las
implementan correctamente, especialmente a la hora de mostrar los
resultados, con 'cout >>' o con 'printf()'.

Si requieres la mayor precisión, entonces sugiero crear tu propia
representación de números. Por ejemplo, puedes crear un array de
dígitos para representar una secuencia de dígitos, junto con las
operaciones aritméticas, puedes representar un número. La cantidad de
dígitos corre de tu cuenta; puede ser diez o diez mil dígitos.


Espero que todo esto te ayude.

Steven