下面的源代码应该是找到两条曲线的交点.

但是,该函数无法检测交点.

我怎么才能修好它呢?

enter image description here

using MathNet.Numerics.Interpolation;
using System.Collections.Generic;
using System.Linq;

public static Vec2 FindIntersectionOfTwoCurves(List<double> xList1, List<double> yList1,
                                               List<double> xList2, List<double> yList2)
{
    if (xList1 == null || yList1 == null || xList2 == null || yList2 == null ||
        xList1.Count != yList1.Count || xList2.Count != yList2.Count)
        return null;

    IInterpolation interpolation1 = Interpolate.Linear(xList1, yList1);
    IInterpolation interpolation2 = Interpolate.Linear(xList2, yList2);

    double lowerBound = Math.Max(xList1.Min(), xList2.Min());
    double upperBound = Math.Min(xList1.Max(), xList2.Max());

    double step = (upperBound - lowerBound) / 1000; // adjust the step size as needed
    for (double x = lowerBound; x <= upperBound; x += step)
    {
        double y1 = interpolation1.Interpolate(x);
        double y2 = interpolation2.Interpolate(x);

        if (Math.Abs(y1 - y2) < 1e-7)
        {
            return new Vec2(x, y1);
        }
    }

    return null;
}

public class Vec2
{
    public double X { get; set; }
    public double Y { get; set; }

    public Vec2(double x, double y)
    {
        X = x;
        Y = y;
    }
}

推荐答案

由于您正在try 寻找具有非连续一阶导数的两条曲线的交点,因此应使用Bisection method.

一个快速编写的C#实现如下,二进制搜索的实现是从here开始,如果X是线性间隔的,则可以用更简单的线性内插来代替二进制搜索.

namespace ConsoleApp2
{
    public class LinearInterpolator
    {
        private static int BinarySearch<T>(IList<T> list, T value)
        {
            if (list == null)
                throw new ArgumentNullException("list");
            var comp = Comparer<T>.Default;
            int lo = 0, hi = list.Count - 1;
            while (lo < hi)
            {
                int m = (hi + lo) / 2;  // this might overflow; be careful.
                if (comp.Compare(list[m], value) < 0) lo = m + 1;
                else hi = m - 1;
            }
            if (comp.Compare(list[lo], value) < 0) lo++;
            return lo;
        }

        private static int FindFirstIndexGreaterThanOrEqualTo<T>
                                  (List<T> sortedList, T key)
        {
            return BinarySearch(sortedList, key);
        }

        List<double> x_values;
        List<double> y_values;
        public LinearInterpolator(List<double> x, List<double> y)
        {
            // quick argsort
            List<int> indicies = x.AsEnumerable().Select((v,i) => new {obj = v, index = i}).OrderBy(c=> c.obj).Select(c => c.index).ToList();
            x_values = indicies.Select(i => x[i]).ToList();
            y_values = indicies.Select(i => y[i]).ToList();
        }

        public double Interpolate(double x)
        {
            int index = FindFirstIndexGreaterThanOrEqualTo(x_values, x);
            if (index == 0)
            {
                return y_values[0];
            }
            double y1 = y_values[index-1];
            double y2 = y_values[index];
            double x1 = x_values[index-1];
            double x2 = x_values[index];
            return (x - x1) / (x2 - x1) * (y2 - y1) + y1;
        }

    }

    class IntersectionFinder
    {
        public static Nullable<(double,double)> FindIntersectionOfTwoCurves(List<double> xList1, List<double> yList1,
                                                   List<double> xList2, List<double> yList2)
        {
            if (xList1 == null || yList1 == null || xList2 == null || yList2 == null ||
                xList1.Count != yList1.Count || xList2.Count != yList2.Count)
                return null;

            LinearInterpolator interpolator1 = new LinearInterpolator(xList1, yList1);
            LinearInterpolator interpolator2 = new LinearInterpolator(xList2, yList2);

            double lowerBound = Math.Max(xList1.Min(), xList2.Min());
            double upperBound = Math.Min(xList1.Max(), xList2.Max());

            if (lowerBound > upperBound) // X axes don't overlap
            {
                return null;
            }

            double diff_start = interpolator1.Interpolate(lowerBound) - interpolator2.Interpolate(lowerBound);
            double diff_end = interpolator1.Interpolate(upperBound) - interpolator2.Interpolate(upperBound);

            if ((diff_start > 0 && diff_end > 0) || (diff_start < 0 && diff_end < 0)) // intersection doesn't exist
            {
                return null;
            }

            double mid = (lowerBound + upperBound) / 2;
            double diff_mid = interpolator1.Interpolate(mid) - interpolator2.Interpolate(mid);

            int iterations = 0;
            while (Math.Abs(diff_mid) > 1e-7)
            {
                if (diff_start > diff_end) // list1 is higher
                {
                    if (diff_mid > 0) // mid is also higher, intersection in right side
                    {
                        lowerBound = mid;
                        diff_start = diff_mid;
                    }
                    else // mid is lower, intersection in left side
                    {
                        upperBound = mid;
                        diff_end = diff_mid;
                    }
                }
                else // list 2 is higher
                {
                    if (diff_mid < 0) 
                    {
                        lowerBound = mid;
                        diff_start = diff_mid;
                    }
                    else 
                    {
                        upperBound = mid;
                        diff_end = diff_mid;
                    }
                }
                mid = (lowerBound + upperBound) / 2;
                diff_mid = interpolator1.Interpolate(mid) - interpolator2.Interpolate(mid);
                iterations++;
                if (iterations > 10000) // prevent infinite loop if Y is discontinuous
                {
                    return null;
                }
            }

            return new (mid, interpolator1.Interpolate(mid));
        }
    }

    internal class Program
    {
        static void Main(string[] args)
        {
            List<double> xList1 = [ 1, 1.5, 2, 3.5, 4 ];
            List<double> xList2 = [ 1, 2, 3, 4, 5 ];
            List<double> yList1 = [ 0, 2, 4, 6,  8 ];
            List<double> yList2 = [ 8, 6, 4, 2, 0 ];
            Console.WriteLine(IntersectionFinder.FindIntersectionOfTwoCurves(xList1, yList1, xList2, yList2).ToString());
        }
    }
}
(2.5999999940395355, 4.799999992052714)

enter image description here

Csharp相关问答推荐

获取ASP.NET核心身份认证cookie名称

通过EFCore上传大量数据.

try 还原包时出错:Dapper已经为System.Data.SQLClient定义了依赖项''''

迭代C#List并在数据库中 for each 元素执行SELECT语句—更有效的方式?

如何使用while循环实现异常处理

当通过Google的Gmail Api发送邮件时,签名会产生dkim = neutral(正文散列未验证)'

AsNoTrackingWithIdentitySolutions()似乎不起作用?

更新产品但丢失产品ASP.NET Core的形象

如何将端点(或с匹配请求并判断其路径)添加到BCL?

带有可选参数的模拟方法返回意外的不同值,具体取决于可选的默认值

当空判断结果赋给变量时,为什么会出现可能空异常警告的解引用?

C#动态设置ServerReport报表参数

源代码生成器:CS8795分部方法';Class1.GetS2(字符串)';必须有实现部分,因为它有可访问性修饰符?

如何在Akka.NET中重新启动执行元时清除邮箱

如何在使用属性 Select 器时判断是否可以为空

处理方法内部过滤与外部过滤

单位中快照的倾斜方向

使用Try-Catch-Finally为API端点处理代码--有什么缺点?

为什么我不能在固定语句中使用外部函数?

无法通过服务控制台启动.NET Core 6.0服务(错误1053)