实现 Math.Asin 迈克劳林(泰勒)展开式,结果比Math.Asin 慢一倍

项目中需要快速求解Asin(x) 的近似值,原以为用泰勒展开式会快一些,结果比原生的慢一倍。

Math.ASin
        Time Elapsed:   9ms
        Gen 0:          0
        Gen 1:          0
        Gen 2:          0
Maclaurin.ASin
        Time Elapsed:   17ms
        Gen 0:          4
        Gen 1:          0
        Gen 2:          0

各位,谁有能力改进?

附:

http://www.mathportal.org/formulas/pdf/taylor-series-formulas.pdf

http://pages.pacificcoast.net/~cazelais/260/maclaurin.pdf

using System;
using System.Collections.Generic;
using System.Linq;
using System.Runtime.Remoting.Messaging;
using System.Text;
using Diagnostics;

namespace Asin
{
    class Program
    {
        static void Main(string[] args)
        {
            int count = 100000;
            List<double> values = new List<double>(count);

            Random r = new Random();
            for (var i = 0; i <= count; ++i)
            {
                values .Add(r.NextDouble() * 2 - 1);
            }


            CodeTime.Init();
            int? iter = 0;
            CodeTime.Timer("Math.ASin", count, () =>
            {
                var i = iter.Value + 1;
                iter = i;

                Math.Asin(values[i]);

            });

            iter = 0;
            CodeTime.Timer("Maclaurin.ASin", count, () =>
            {
                var i = iter.Value + 1;
                iter = i;

                Maclaurin.Asin(values[i],3);

            });

            while (true)
            {
                iter = 0;
                CodeTime.Timer("Math.ASin", count, () =>
                {
                    var i = iter.Value + 1;
                    iter = i;

                    Math.Asin(values[i]);

                });

                iter = 0;
                CodeTime.Timer("Maclaurin.ASin", count, () =>
                {
                    var i = iter.Value + 1;
                    iter = i;

                    Maclaurin.Asin(values[i], 3);

                });
            }
           

            //var ret = Maclaurin.Asin(0.5, 3);
            //var ret2 = Math.Asin(0.5);
            //Console.WriteLine(ret);
            //Console.WriteLine(ret2);
            Console.ReadLine();
        }
    }

    class Maclaurin
    {
        class ASinImpl
        {
            private List<double> quotieties = new List<double>();
            private IEnumerator<double> computeQuotieties = null;

            public ASinImpl()
            {
                this.computeQuotieties = ComputeQuotiety();
            }

            public double Calc(double v, int precision = 2)
            {
                if (quotieties.Count < precision)
                {
                    for (var i = quotieties.Count; i < precision; ++i)
                    {
                        computeQuotieties.MoveNext();
                        quotieties.Add(computeQuotieties.Current);
                    }
                }

                double ret = 0;
                var values = ComputeValues(v);
                for (int i = 0; i < precision; ++i)
                {
                    values.MoveNext();
                    ret += quotieties[i]*values.Current;
                }

                return ret;
            }

            private IEnumerator<double> ComputeValues(double v)
            {
                double ret = 1;
                double q = v*v;
                for(int i = 0;;++i)
                {
                   
                    if (i == 0)
                    {
                        ret = v;
                        yield return ret;
                    }
                    else
                    {
                        ret *= q;
                        yield return ret;
                    }
                }

                throw new NotImplementedException();
            }

            private IEnumerator<double> ComputeQuotiety()
            {
                for (int i = 0;; i++)
                {
                   
                    double up = Factorial(2*i);

                    double down = Math.Pow(Math.Pow(2, i)*Factorial(i), 2)*(2*i + 1);

                    double quotiety = up/down;

                    yield return quotiety;


                }

                throw new NotImplementedException();
            }

            private long Factorial(long v )
            { 
                if( v < 0)
                    throw new ArgumentOutOfRangeException("v");

                if (v == 0)
                    return 1;

                if (v == 1)
                    return 1;

                long ret = 1;
                for (int i = 2; i <= v; ++i)
                {
                    ret *= i;
                }

                return ret;
            }
        }


        private static ASinImpl asinImpl = new ASinImpl();

        /// <summary>
        /// 
        /// </summary>
        /// <param name="v"></param>
        /// <param name="precision"></param>
        /// <returns></returns>
        public static double Asin(double v, int precision)
        {
            if (v < -1 || v > 1)
            {
                throw new ArgumentOutOfRangeException("v");
            }

            return asinImpl.Calc(v, precision);
        }
    }
}

 

通过一下优化:基本持平

  class ASinImpl
    {
        private readonly int _precision;
        private double[] _quotieties = null;
        private long[] _factorials =null;
        public ASinImpl(int precision = 3)
        {
            _precision = precision;
            _quotieties = new double[precision + 1];
            _factorials = new long[(precision + 2)*2 + 1];
            Factorial();
            ComputeQuotiety();
        }

        public double Calc(double v)
        {
            unchecked
            {
                double retVal = 0;

                double vVal = 1;
                double q = v * v;
                unsafe
                {
                    fixed (double* pq = _quotieties)
                    {
                        for (int i = 0; i < _precision; ++i)
                        {

                            if (i == 0)
                            {
                                vVal = v;
                                //yield return ret;
                                retVal += pq[i] * vVal;
                            }
                            else
                            {
                                vVal *= q;
                                //yield return ret;
                                retVal += pq[i] * vVal;
                            }
                        }

                        return retVal;
                    }

                   
                }
               
            }
            
        }


        private void ComputeQuotiety()
        {
            unchecked
            {
                int precision = _quotieties.Length;
                for (int i = 0; i < precision; i++)
                {

                    double up = _factorials[2*i];

                    double down = Math.Pow(Math.Pow(2, i)*_factorials[i], 2)*(2*i + 1);

                    double quotiety = up/down;

                    _quotieties[i] = quotiety;


                }
            }

        }

        private void Factorial()
        {
            unchecked
            {
                int precision = _factorials.Length ;
                long ret = 1;
                for (long v = 0; v < precision; ++v)
                {
                    if (v == 0)
                    {
                        this._factorials[v] = 1;
                    }
                    else if (v == 1){
                            this._factorials[v] = 1;
                    }
                    else
                    {
                        ret *= v;

                        this._factorials[v] = ret;
                    }

                }
            }
        }
    }

 

posted @ 2014-10-10 10:14  阿牛  阅读(1234)  评论(0编辑  收藏  举报