Pregunta

El marco .net proporciona en la clase de Matemáticas un método para alimentar el doble. Pero por requerimiento de precisión necesito elevar un decimal a una potencia decimal [Pow (decimal a, decimal b)]. ¿El marco tiene tal función? ¿Alguien sabe de una biblioteca con este tipo de función?

¿Fue útil?

Solución

Para resolver mi problema, encontré algunas series de expansión , y les hice implementar para resolver la ecuación X ^ n = e ^ (n * ln x).

// Adjust this to modify the precision
public const int ITERATIONS = 27;

// power series
public static decimal DecimalExp(decimal power)
{
    int iteration = ITERATIONS;
    decimal result = 1; 
    while (iteration > 0)
    {
        fatorial = Factorial(iteration);
        result += Pow(power, iteration) / fatorial;
        iteration--;
    }
    return result;
}

// natural logarithm series
public static decimal LogN(decimal number)
{
    decimal aux = (number - 1);
    decimal result = 0;
    int iteration = ITERATIONS;
    while (iteration > 0)
    {
        result += Pow(aux, iteration) / iteration;
        iteration--;
    }
    return result;
}

// example
void main(string[] args)
{
    decimal baseValue = 1.75M;
    decimal expValue = 1/252M;
    decimal result = DecimalExp(expValue * LogN(baseValue));
}

Las funciones Pow () y Factorial () son simples porque el poder es siempre un int (dentro de la serie de poder).

Otros consejos

Esto debería ser el más rápido para un exponente entero positivo y una base decimal:

// From http://www.daimi.au.dk/~ivan/FastExpproject.pdf
// Left to Right Binary Exponentiation
public static decimal Pow(decimal x, uint y){
    decimal A = 1m;
    BitArray e = new BitArray(BitConverter.GetBytes(y));
    int t = e.Count;

    for (int i = t-1; i >= 0; --i) {
        A *= A;
        if (e[i] == true) {
            A *= x;
        }
    }
    return A;
}

Aquí hay un programa C # para implementar manualmente Math.Pow () con un mayor grado de precisión que la implementación basada en el doble de .NET. Corte y pegue en linqpad para ejecutarlo inmediatamente, o cambie los .Dump () s a Console.WriteLines.

He incluido una prueba del resultado. La prueba es la siguiente:

  1. Objetivo = .4% anual con capitalización diaria en 10 000
  2. Answer = debería ser 10 040
  3. Cómo = decimal b = 10000; para (int i = 0; i < 365; i ++) {b * = tasa; } donde rate = (1.004) ^ (1/365)

He probado 3 implementaciones de tasa: (1) Cálculo manual (2) Excel (3) Math.Pow

El cálculo manual tiene el mayor grado de precisión. Los resultados son:

Manually calculated rate:   1.0000109371043837652682334292
Excel rate:                 1.000010937104383712500000M [see formula =(1.004)^(1/365)]
Math.Pow rate:              1.00001093710438

Manual - .4%pa on R10,000:  10040.000000000000000000000131 
Excel - .4%pa on R10,000:   10039.999999999806627646709094 
Math.Pow - .4%pa on R10,000:10039.999999986201948942509648

También dejé algunos trabajos adicionales allí, que usé para establecer cuál era el factorial más alto que puede caber en un ulong (= 22).

Código de Linqpad:

/*
a^b = exp(b * ln(a))
    ln(a) = log(1-x) = - x - x^2/2 - x^3/3 - ...   (where |x| < 1)
        x: a = 1-x    =>   x = 1-a = 1 - 1.004 = -.004
    y = b * ln(a)
    exp(y) = 1 + y + y^2/2 + x^3/3! + y^4/4! + y^5/5! + ...
        n! = 1 * 2 * ... * n        
*/

/*
//
// Example: .4%pa on R10,000 with daily compounding
//

Manually calculated rate:   1.0000109371043837652682334292
Excel rate:                 1.000010937104383712500000M =(1.004)^(1/365)
Math.Pow rate:              1.00001093710438

Manual - .4%pa on R10,000:  10040.000000000000000000000131 
Excel - .4%pa on R10,000:   10039.999999999806627646709094 
Math.Pow - .4%pa on R10,000:10039.999999986201948942509648 

*/

static uint _LOOPS = 10;    // Max = 22, no improvement in accuracy after 10 in this example scenario
//  8: 1.0000109371043837652682333497
//  9: 1.0000109371043837652682334295
// 10: 1.0000109371043837652682334292
// ...
// 21: 1.0000109371043837652682334292
// 22: 1.0000109371043837652682334292

// http://www.daimi.au.dk/~ivan/FastExpproject.pdf
// Left to Right Binary Exponentiation
public static decimal Pow(decimal x, uint y)
{
    if (y == 1)
        return x;

    decimal A = 1m;
    BitArray e = new BitArray(BitConverter.GetBytes(y));
    int t = e.Count;

    for (int i = t-1; i >= 0; --i) {
        A *= A;
        if (e[i] == true) {
            A *= x;
        }
    }
    return A;
}

// http://stackoverflow.com/questions/429165/raising-a-decimal-to-a-power-of-decimal
// natural logarithm series
public static decimal ln(decimal a)
{
    /*
    ln(a) = log(1-x) = - x - x^2/2 - x^3/3 - ...   (where |x| < 1)
        x: a = 1-x    =>   x = 1-a = 1 - 1.004 = -.004
    */
    decimal x = 1 - a;
    if (Math.Abs(x) >= 1)
        throw new Exception("must be 0 < a < 2");

    decimal result = 0;
    uint iteration = _LOOPS;
    while (iteration > 0)
    {
        result -= Pow(x, iteration) / iteration;
        iteration--;
    }
    return result;
}

public static ulong[] Fact = new ulong[] {
    1L,
    1L * 2,
    1L * 2 * 3,
    1L * 2 * 3 * 4,
    1L * 2 * 3 * 4 * 5,
    1L * 2 * 3 * 4 * 5 * 6,
    1L * 2 * 3 * 4 * 5 * 6 * 7,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17 * 18,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17 * 18 * 19,
    1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17 * 18 * 19 * 20,
    14197454024290336768L, //1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17 * 18 * 19 * 20 * 21,        // NOTE: Overflow during compilation
    17196083355034583040L, //1L * 2 * 3 * 4 * 5 * 6 * 7 * 8 * 9 * 10 * 11 * 12 * 13 * 14 * 15 * 16 * 17 * 18 * 19 * 20 * 21 * 22    // NOTE: Overflow during compilation
};

// http://stackoverflow.com/questions/429165/raising-a-decimal-to-a-power-of-decimal
// power series
public static decimal exp(decimal y)
{
    /*
    exp(y) = 1 + y + y^2/2 + x^3/3! + y^4/4! + y^5/5! + ...
    */

    uint iteration = _LOOPS;
    decimal result = 1; 
    while (iteration > 0)
    {
        //uint fatorial = Factorial(iteration);
        ulong fatorial = Fact[iteration-1];
        result += (Pow(y, iteration) / fatorial);
        iteration--;
    }
    return result;
}

void Main()
{   
    decimal a = 1.004M;
    decimal b = 1/365M;

    decimal _ln = ln(a);
    decimal y = b * _ln;
    decimal result = exp(y);
    result.Dump("Manual rate");

    decimal excel = 1.000010937104383712500000M;    // =(1.004)^(1/365)
    excel.Dump("Excel rate");


    decimal m = (decimal)Math.Pow((double)a,(double)b);
    m.Dump("Math.Pow rate");

    //(result - excel).Dump("Diff: Manual - Excel");
    //(m - excel).Dump("Diff: Math.Pow - Excel");

    var f = new DateTime(2013,1,1);
    var t = new DateTime(2014,1,1);
    Test(f, t, 10000, result, "Manual - .4%pa on R10,000");
    Test(f, t, 10000, excel, "Excel - .4%pa on R10,000");
    Test(f, t, 10000, m, "Math.Pow - .4%pa on R10,000");
}

decimal Test(DateTime f, DateTime t, decimal balance, decimal rate, string whichRate)
{
    int numInterveningDays = (t.Date - f.Date).Days;
    var value = balance;
    for (int i = 0; i < numInterveningDays; ++i)
    {
        value *= rate;
    }
    value.Dump(whichRate);
    return value - balance;
}

/*

// Other workings:

//
// Determine maximum Factorial for use in ln(a)
//

ulong max    =  9,223,372,036,854,775,807 * 2   // see http://msdn.microsoft.com/en-us/library/ctetwysk.aspx
Factorial 21 = 14,197,454,024,290,336,768
Factorial 22 = 17,196,083,355,034,583,040
Factorial 23 = 8,128,291,617,894,825,984 (Overflow)

public static uint Factorial_uint(uint i)
{
    // n! = 1 * 2 * ... * n
    uint n = i;
    while (--i > 1)
    {
        n *= i;
    }
    return n;
}

public static ulong Factorial_ulong(uint i)
{
    // n! = 1 * 2 * ... * n
    ulong n = i;
    while (--i > 1)
    {
        n *= i;
    }
    return n;
}

void Main()
{
    // Check max ulong Factorial
    ulong prev = 0;
    for (uint i = 1; i < 24; ++i)
    {
        ulong cur = Factorial_ulong(i);
        cur.Dump(i.ToString());
        if (cur < prev)
        {
            throw new Exception("Overflow");
        }
        prev = cur;
    }
}
*/

Creo que depende mucho del número que planea conectar. Si 'a' y 'b' no son números 'agradables', entonces es probable que obtenga un valor que no termina y que es imposible de almacenar y si C # BigDecimal se comporta como Java BigDecimal, probablemente arroje una excepción en tal caso.

¿Estás seguro de que realmente quieres hacer esto? Una multiplicación decimal es aproximadamente 40 veces más lenta que la de double, por lo que esperaría que un decimal Math.Pow() sea prácticamente inutilizable.

Sin embargo, si espera solo potencias enteras, le sugiero que use el algoritmo de potencia basado en enteros que ya se discutió aquí en SO.

Licenciado bajo: CC-BY-SA con atribución
No afiliado a StackOverflow
scroll top