1 文本格式

using System;

namespace Legalsoft.Truffer
{
    /// <summary>
    /// 重心有理插值对象
    /// Barycentric rational interpolation object. 
    /// After constructing the object, 
    /// call interp for interpolated values.
    /// Note that no error estimate dy is calculated.
    /// </summary>
    public class BaryRat_interp : Base_interp
    {
        private double[] w { get; set; }
        private int d { get; set; }

        /// <summary>
        /// Constructor arguments are x and y vectors of length n, 
        /// and order d of desired approximation.
        /// </summary>
        /// <param name="xv"></param>
        /// <param name="yv"></param>
        /// <param name="dd"></param>
        /// <exception cref="Exception"></exception>
        public BaryRat_interp(double[] xv, double[] yv, int dd) : base(xv, yv[0], xv.Length)
        {
            this.w = new double[n];
            this.d = dd;
            if (n <= d)
            {
                throw new Exception("d too large for number of points in BaryRat_interp");
            }
            for (int k = 0; k < n; k++)
            {
                int imin = Math.Max(k - d, 0);
                int imax = k >= n - d ? n - d - 1 : k;
                double temp = (imin & 1) != 0 ? -1.0 : 1.0;
                double sum = 0.0;
                for (int i = imin; i <= imax; i++)
                {
                    int jmax = Math.Min(i + d, n - 1);
                    double term = 1.0;
                    for (int j = i; j <= jmax; j++)
                    {
                        if (j == k)
                        {
                            continue;
                        }
                        term *= (xx[k] - xx[j]);
                    }
                    term = temp / term;
                    temp = -temp;
                    sum += term;
                }
                w[k] = sum;
            }
        }

        /// <summary>
        /// Use equation(3.4.9) to compute the 
        /// barycentric rational interpolant.
        /// Note that jl is not used since 
        /// the approximation is global; 
        /// it is included only
        /// for compatibility with Base_interp.
        /// </summary>
        /// <param name="jl"></param>
        /// <param name="x"></param>
        /// <returns></returns>
        public override double rawinterp(int jl, double x)
        {
            double num = 0;
            double den = 0;
            for (int i = 0; i < n; i++)
            {
                double h = x - xx[i];
                //if (h == 0.0)
                if (Math.Abs(h) <= float.Epsilon)
                {
                    return yy[i];
                }
                else
                {
                    double temp = w[i] / h;
                    num += temp * yy[i];
                    den += temp;
                }
            }
            return num / den;
        }

        /// <summary>
        /// No need to invoke hunt or locate since 
        /// the interpolation is global, so
        /// override interp to simply call rawinterp 
        /// directly with a dummy value of jl.
        /// </summary>
        /// <param name="x"></param>
        /// <returns></returns>
        public new double interp(double x)
        {
            return rawinterp(1, x);
        }
    }
}
 

2 代码格式

using System;namespace Legalsoft.Truffer
{/// <summary>/// 重心有理插值对象/// Barycentric rational interpolation object. /// After constructing the object, /// call interp for interpolated values./// Note that no error estimate dy is calculated./// </summary>public class BaryRat_interp : Base_interp{private double[] w { get; set; }private int d { get; set; }/// <summary>/// Constructor arguments are x and y vectors of length n, /// and order d of desired approximation./// </summary>/// <param name="xv"></param>/// <param name="yv"></param>/// <param name="dd"></param>/// <exception cref="Exception"></exception>public BaryRat_interp(double[] xv, double[] yv, int dd) : base(xv, yv[0], xv.Length){this.w = new double[n];this.d = dd;if (n <= d){throw new Exception("d too large for number of points in BaryRat_interp");}for (int k = 0; k < n; k++){int imin = Math.Max(k - d, 0);int imax = k >= n - d ? n - d - 1 : k;double temp = (imin & 1) != 0 ? -1.0 : 1.0;double sum = 0.0;for (int i = imin; i <= imax; i++){int jmax = Math.Min(i + d, n - 1);double term = 1.0;for (int j = i; j <= jmax; j++){if (j == k){continue;}term *= (xx[k] - xx[j]);}term = temp / term;temp = -temp;sum += term;}w[k] = sum;}}/// <summary>/// Use equation(3.4.9) to compute the /// barycentric rational interpolant./// Note that jl is not used since /// the approximation is global; /// it is included only/// for compatibility with Base_interp./// </summary>/// <param name="jl"></param>/// <param name="x"></param>/// <returns></returns>public override double rawinterp(int jl, double x){double num = 0;double den = 0;for (int i = 0; i < n; i++){double h = x - xx[i];//if (h == 0.0)if (Math.Abs(h) <= float.Epsilon){return yy[i];}else{double temp = w[i] / h;num += temp * yy[i];den += temp;}}return num / den;}/// <summary>/// No need to invoke hunt or locate since /// the interpolation is global, so/// override interp to simply call rawinterp /// directly with a dummy value of jl./// </summary>/// <param name="x"></param>/// <returns></returns>public new double interp(double x){return rawinterp(1, x);}}
}