C#,数值计算——插值和外推,二维三次样条插值(Spline2D_interp)的计算方法与源程序

news2025/2/26 23:43:25

1 文本格式

using System;

namespace Legalsoft.Truffer
{
    /// <summary>
    /// 二维三次样条插值
    /// Object for two-dimensional cubic spline interpolation on a matrix.Construct
    /// with a vector of x1 values, a vector of x2 values, and a matrix of tabulated
    /// function values yij.Then call interp for interpolated values.
    /// </summary>
    public class Spline2D_interp
    {
        private int[,] bcucof_wt = new int[,] {
            {  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, 0 },
            {  0,  0,  0,  0,  0,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0 },
            { -3,  0,  0,  3,  0,  0,  0,  0, -2,  0,  0, -1,  0,  0,  0,  0 },
            {  2,  0,  0, -2,  0,  0,  0,  0,  1,  0,  0,  1,  0,  0,  0,  0 },
            {  0,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0 },
            {  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  1,  0,  0,  0 },
            {  0,  0,  0,  0, -3,  0,  0,  3,  0,  0,  0,  0, -2,  0,  0, -1 },
            {  0,  0,  0,  0,  2,  0,  0, -2,  0,  0,  0,  0,  1,  0,  0,  1 },
            { -3,  3,  0,  0, -2, -1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0 },
            {  0,  0,  0,  0,  0,  0,  0,  0, -3,  3,  0,  0, -2, -1,  0,  0 },
            {  9, -9,  9, -9,  6,  3, -3, -6,  6, -6, -3,  3,  4,  2,  1,  2 },
            { -6,  6, -6,  6, -4, -2,  2,  4, -3,  3,  3, -3, -2, -1, -1, -2 },
            {  2, -2,  0,  0,  1,  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0 },
            {  0,  0,  0,  0,  0,  0,  0,  0,  2, -2,  0,  0,  1,  1,  0,  0 },
            { -6,  6, -6,  6, -3, -3,  3,  3, -4,  4,  2, -2, -2, -2, -1, -1 },
            {  4, -4,  4, -4,  2,  2, -2, -2,  2, -2, -2,  2,  1,  1,  1,  1 }
        };

        private int m { get; set; }
        private int n { get; set; }
        private double[,] y { get; set; }
        private double[] x1 { get; set; }
        private double[] yv { get; set; }
        private Spline_interp[] srp { get; set; }

        public Spline2D_interp(double[] x1v, double[] x2v, double[,] ym)
        {
            this.m = x1v.Length;
            this.n = x2v.Length;
            this.y = ym;
            this.yv = new double[m];
            this.x1 = x1v;
            this.srp = new Spline_interp[m];
            for (int i = 0; i < m; i++)
            {
                srp[i] = new Spline_interp(x2v, y[i, 0]);
            }
        }

        public double interp(double x1p, double x2p)
        {
            for (int i = 0; i < m; i++)
            {
                yv[i] = (srp[i]).interp(x2p);
            }
            Spline_interp scol = new Spline_interp(x1, yv);
            return scol.interp(x1p);
        }

        /// <summary>
        /// Given arrays y[0..3], y1[0..3], y2[0..3], and y12[0..3], containing the
        /// function, gradients, and cross-derivative at the four grid points of a
        /// rectangular grid cell(numbered counterclockwise from the lower left), and
        /// given d1 and d2, the length of the grid cell in the 1 and 2 directions,
        /// this routine returns the table c[0..3][0..3] that is used by routine bcuint
        /// for bicubic interpolation.
        /// </summary>
        /// <param name="y"></param>
        /// <param name="y1"></param>
        /// <param name="y2"></param>
        /// <param name="y12"></param>
        /// <param name="d1"></param>
        /// <param name="d2"></param>
        /// <param name="c"></param>
        private void bcucof(double[] y, double[] y1, double[] y2, double[] y12, double d1, double d2, double[,] c)
        { 
            //int[,] bcucof_wt = new int[16, 16]{ bcucof_wt_d };
            double d1d2 = d1 * d2;
            double[] cl = new double[16];
            double[] x = new double[16];

            for (int i = 0; i < 4; i++)
            {
                x[i] = y[i];
                x[i + 4] = y1[i] * d1;
                x[i + 8] = y2[i] * d2;
                x[i + 12] = y12[i] * d1d2;
            }
            for (int i = 0; i < 16; i++)
            {
                double xx = 0.0;
                for (int k = 0; k < 16; k++)
                {
                    xx += bcucof_wt[i, k] * x[k];
                }
                cl[i] = xx;
            }
            int l = 0;
            for (int i = 0; i < 4; i++)
            {
                for (int j = 0; j < 4; j++)
                {
                    c[i, j] = cl[l++];
                }
            }
        }

        /// <summary>
        /// Bicubic interpolation within a grid square.Input quantities are
        /// y, y1, y2, y12 (as described in bcucof); x1l and x1u, the lower and upper
        /// coordinates of the grid square in the 1 direction; x2l and x2u likewise for
        /// the 2 direction; and x1, x2, the coordinates of the desired point for the
        /// interpolation.The interpolated function value is returned as ansy, and the
        /// interpolated gradient values as ansy1 and ansy2.This routine calls bcucof.
        /// </summary>
        /// <param name="y"></param>
        /// <param name="y1"></param>
        /// <param name="y2"></param>
        /// <param name="y12"></param>
        /// <param name="x1l"></param>
        /// <param name="x1u"></param>
        /// <param name="x2l"></param>
        /// <param name="x2u"></param>
        /// <param name="x1"></param>
        /// <param name="x2"></param>
        /// <param name="ansy"></param>
        /// <param name="ansy1"></param>
        /// <param name="ansy2"></param>
        /// <exception cref="Exception"></exception>
        public void bcuint(double[] y, double[] y1, double[] y2, double[] y12, double x1l, double x1u, double x2l, double x2u, double x1, double x2, ref double ansy, ref double ansy1, ref double ansy2)
        {
            double d1 = x1u - x1l;
            double d2 = x2u - x2l;
            double[,] c = new double[4, 4];
            bcucof(y, y1, y2, y12, d1, d2, c);
            //if (x1u == x1l || x2u == x2l)
            if (Math.Abs(x1u - x1l) <= float.Epsilon || Math.Abs(x2u - x2l) <= float.Epsilon)
            {
                throw new Exception("Bad input in routine bcuint");
            }
            double t = (x1 - x1l) / d1;
            double u = (x2 - x2l) / d2;
            ansy = ansy2 = ansy1 = 0.0;
            for (int i = 3; i >= 0; i--)
            {
                ansy = t * ansy + ((c[i, 3] * u + c[i, 2]) * u + c[i, 1]) * u + c[i, 0];
                ansy2 = t * ansy2 + (3.0 * c[i, 3] * u + 2.0 * c[i, 2]) * u + c[i, 1];
                ansy1 = u * ansy1 + (3.0 * c[3, i] * t + 2.0 * c[2, i]) * t + c[1, i];
            }
            ansy1 /= d1;
            ansy2 /= d2;
        }

    }
}

2 代码格式

using System;

namespace Legalsoft.Truffer
{
    /// <summary>
    /// 二维三次样条插值
    /// Object for two-dimensional cubic spline interpolation on a matrix.Construct
    /// with a vector of x1 values, a vector of x2 values, and a matrix of tabulated
    /// function values yij.Then call interp for interpolated values.
    /// </summary>
    public class Spline2D_interp
    {
        private int[,] bcucof_wt = new int[,] {
            {  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0, 0 },
            {  0,  0,  0,  0,  0,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0 },
            { -3,  0,  0,  3,  0,  0,  0,  0, -2,  0,  0, -1,  0,  0,  0,  0 },
            {  2,  0,  0, -2,  0,  0,  0,  0,  1,  0,  0,  1,  0,  0,  0,  0 },
            {  0,  0,  0,  0,  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0 },
            {  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0,  1,  0,  0,  0 },
            {  0,  0,  0,  0, -3,  0,  0,  3,  0,  0,  0,  0, -2,  0,  0, -1 },
            {  0,  0,  0,  0,  2,  0,  0, -2,  0,  0,  0,  0,  1,  0,  0,  1 },
            { -3,  3,  0,  0, -2, -1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0 },
            {  0,  0,  0,  0,  0,  0,  0,  0, -3,  3,  0,  0, -2, -1,  0,  0 },
            {  9, -9,  9, -9,  6,  3, -3, -6,  6, -6, -3,  3,  4,  2,  1,  2 },
            { -6,  6, -6,  6, -4, -2,  2,  4, -3,  3,  3, -3, -2, -1, -1, -2 },
            {  2, -2,  0,  0,  1,  1,  0,  0,  0,  0,  0,  0,  0,  0,  0,  0 },
            {  0,  0,  0,  0,  0,  0,  0,  0,  2, -2,  0,  0,  1,  1,  0,  0 },
            { -6,  6, -6,  6, -3, -3,  3,  3, -4,  4,  2, -2, -2, -2, -1, -1 },
            {  4, -4,  4, -4,  2,  2, -2, -2,  2, -2, -2,  2,  1,  1,  1,  1 }
        };

        private int m { get; set; }
        private int n { get; set; }
        private double[,] y { get; set; }
        private double[] x1 { get; set; }
        private double[] yv { get; set; }
        private Spline_interp[] srp { get; set; }

        public Spline2D_interp(double[] x1v, double[] x2v, double[,] ym)
        {
            this.m = x1v.Length;
            this.n = x2v.Length;
            this.y = ym;
            this.yv = new double[m];
            this.x1 = x1v;
            this.srp = new Spline_interp[m];
            for (int i = 0; i < m; i++)
            {
                srp[i] = new Spline_interp(x2v, y[i, 0]);
            }
        }

        public double interp(double x1p, double x2p)
        {
            for (int i = 0; i < m; i++)
            {
                yv[i] = (srp[i]).interp(x2p);
            }
            Spline_interp scol = new Spline_interp(x1, yv);
            return scol.interp(x1p);
        }

        /// <summary>
        /// Given arrays y[0..3], y1[0..3], y2[0..3], and y12[0..3], containing the
        /// function, gradients, and cross-derivative at the four grid points of a
        /// rectangular grid cell(numbered counterclockwise from the lower left), and
        /// given d1 and d2, the length of the grid cell in the 1 and 2 directions,
        /// this routine returns the table c[0..3][0..3] that is used by routine bcuint
        /// for bicubic interpolation.
        /// </summary>
        /// <param name="y"></param>
        /// <param name="y1"></param>
        /// <param name="y2"></param>
        /// <param name="y12"></param>
        /// <param name="d1"></param>
        /// <param name="d2"></param>
        /// <param name="c"></param>
        private void bcucof(double[] y, double[] y1, double[] y2, double[] y12, double d1, double d2, double[,] c)
        { 
            //int[,] bcucof_wt = new int[16, 16]{ bcucof_wt_d };
            double d1d2 = d1 * d2;
            double[] cl = new double[16];
            double[] x = new double[16];

            for (int i = 0; i < 4; i++)
            {
                x[i] = y[i];
                x[i + 4] = y1[i] * d1;
                x[i + 8] = y2[i] * d2;
                x[i + 12] = y12[i] * d1d2;
            }
            for (int i = 0; i < 16; i++)
            {
                double xx = 0.0;
                for (int k = 0; k < 16; k++)
                {
                    xx += bcucof_wt[i, k] * x[k];
                }
                cl[i] = xx;
            }
            int l = 0;
            for (int i = 0; i < 4; i++)
            {
                for (int j = 0; j < 4; j++)
                {
                    c[i, j] = cl[l++];
                }
            }
        }

        /// <summary>
        /// Bicubic interpolation within a grid square.Input quantities are
        /// y, y1, y2, y12 (as described in bcucof); x1l and x1u, the lower and upper
        /// coordinates of the grid square in the 1 direction; x2l and x2u likewise for
        /// the 2 direction; and x1, x2, the coordinates of the desired point for the
        /// interpolation.The interpolated function value is returned as ansy, and the
        /// interpolated gradient values as ansy1 and ansy2.This routine calls bcucof.
        /// </summary>
        /// <param name="y"></param>
        /// <param name="y1"></param>
        /// <param name="y2"></param>
        /// <param name="y12"></param>
        /// <param name="x1l"></param>
        /// <param name="x1u"></param>
        /// <param name="x2l"></param>
        /// <param name="x2u"></param>
        /// <param name="x1"></param>
        /// <param name="x2"></param>
        /// <param name="ansy"></param>
        /// <param name="ansy1"></param>
        /// <param name="ansy2"></param>
        /// <exception cref="Exception"></exception>
        public void bcuint(double[] y, double[] y1, double[] y2, double[] y12, double x1l, double x1u, double x2l, double x2u, double x1, double x2, ref double ansy, ref double ansy1, ref double ansy2)
        {
            double d1 = x1u - x1l;
            double d2 = x2u - x2l;
            double[,] c = new double[4, 4];
            bcucof(y, y1, y2, y12, d1, d2, c);
            //if (x1u == x1l || x2u == x2l)
            if (Math.Abs(x1u - x1l) <= float.Epsilon || Math.Abs(x2u - x2l) <= float.Epsilon)
            {
                throw new Exception("Bad input in routine bcuint");
            }
            double t = (x1 - x1l) / d1;
            double u = (x2 - x2l) / d2;
            ansy = ansy2 = ansy1 = 0.0;
            for (int i = 3; i >= 0; i--)
            {
                ansy = t * ansy + ((c[i, 3] * u + c[i, 2]) * u + c[i, 1]) * u + c[i, 0];
                ansy2 = t * ansy2 + (3.0 * c[i, 3] * u + 2.0 * c[i, 2]) * u + c[i, 1];
                ansy1 = u * ansy1 + (3.0 * c[3, i] * t + 2.0 * c[2, i]) * t + c[1, i];
            }
            ansy1 /= d1;
            ansy2 /= d2;
        }

    }
}

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.coloradmin.cn/o/1281752.html

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈,一经查实,立即删除!

相关文章

Fisher信息理论与应用

一、概念介绍 Fisher信息量&#xff0c;是一次观测值所能提供的关于未知参数θ的信息量期望值的一种度量。 Fisher信息矩阵&#xff0c;是用利用最大似然函数估计来计算方差矩阵&#xff0c;表示随机变量的一个样本所能提供的关于状态参数在某种意义下的平均信息量。 Fisher…

python 运用pandas 库处理excel 表格数据

文章目录 读取文件查看数据数据选择数据筛选创建新列计算并总结数据分组统计 读取文件 Pandas 是一个强大的数据分析库&#xff0c;它提供了丰富的数据结构和数据分析工具&#xff0c;其中之一是用于读取不同格式文件的 read_* 函数系列。以下是一个简单介绍如何使用 Pandas 读…

qt-C++笔记之组件-分组框QGroupBox

qt-C笔记之组件-分组框QGroupBox code review! 文章目录 qt-C笔记之组件-分组框QGroupBox1.《Qt 6 C开发指南》p752.《Qt 官方文档》3.《Qt 5.12实战》——5.9 分组框控件 1.《Qt 6 C开发指南》p75 2.《Qt 官方文档》 中间段落翻译&#xff1a; 我把示例补充完整&#xff1a; …

CAN 一: CAN基础知识介绍

1、CAN介绍 1.1、什么是CAN? (1)CAN&#xff08;Controller Area Network:控制器局域网&#xff09;&#xff0c;是ISO国际标准化的串行通信协议。为满足汽车产业的“减少线束的数量”、“通过多个LAN&#xff0c;进行大量数据的高速通信”的需求。 (2)CAN总线的发展历史&a…

Apache Flink(六):Apache Flink快速入门 - Flink案例实现

🏡 个人主页:IT贫道_大数据OLAP体系技术栈,Apache Doris,Clickhouse 技术-CSDN博客 🚩 私聊博主:加入大数据技术讨论群聊,获取更多大数据资料。 🔔 博主个人B栈地址:豹哥教你大数据的个人空间-豹哥教你大数据个人主页-哔哩哔哩视频 目录

最新最全的Postman接口测试: postman实现参数化

什么时候会用到参数化 比如&#xff1a;一个模块要用多组不同数据进行测试 验证业务的正确性 Login模块&#xff1a;正确的用户名&#xff0c;密码 成功&#xff1b;错误的用户名&#xff0c;正确的密码 失败 postman实现参数化 在实际的接口测试中&#xff0c;部分参数…

Java 数组另类用法(字符来当数组下标使用)

一、原因 看力扣的时候发现有位大佬使用字符来当数组下标使用。 class Solution {public int lengthOfLongestSubstring(String s) {int result 0;int[] hash new int[130];int i 0;for(int j 0; j < s.length(); j) {while(hash[s.charAt(j)] > 0) {hash[s.charAt…

【java+vue+微信小程序项目】从零开始搭建——健身房管理平台(3)路由导航卫士、主页实现

项目笔记为项目总结笔记,若有错误欢迎指出哟~ 【项目专栏】 【java+vue+微信小程序项目】从零开始搭建——健身房管理平台(1)spring boot项目搭建、vue项目搭建、微信小程序项目搭建 【java+vue+微信小程序项目】从零开始搭建——健身房管理平台(2)后端跨域、登录模块、sp…

【带头学C++】----- 九、类和对象 ---- 9.2 构造函数

目录 9.2 构造函数 9.2.1 构造函数的概述 9.2.2 构造函数定义方法&#xff08;初始化构造函数&#xff09; 9.2.3 提供构造函数的影响 9.2 构造函数 以下是一些C引入构造函数的原因&#xff1a; 初始化对象&#xff1a;构造函数允许在创建对象时立即初始化该对象的成员变量…

数据挖掘实战-基于word2vec的短文本情感分析

&#x1f935;‍♂️ 个人主页&#xff1a;艾派森的个人主页 ✍&#x1f3fb;作者简介&#xff1a;Python学习者 &#x1f40b; 希望大家多多支持&#xff0c;我们一起进步&#xff01;&#x1f604; 如果文章对你有帮助的话&#xff0c; 欢迎评论 &#x1f4ac;点赞&#x1f4…

微机原理——并行接口8255学习1

目录 并行接口特点 可编程并行接口芯片8255 8255端口地址 8255的三种工作方式 8255的两种命令&#xff08;方式命令和C端口命令&#xff09; 由用户扩展的并行接口8255的应用 声光报警器接口设计 步进电机控制接口设计 PA端口实现跑马灯 PB端口实现按键输入 并行接口特…

最大乘积分解(动态规划)

相较于我上一题写的动态规划&#xff0c;这一题比较简单 代码如下&#xff1a; #include<stdio.h>int main(void) {long long n, max[101] {0, 1};scanf("%lld", &n);for(int i 1; i < n; i)max[i] i;for(int i 1; i < n; i)for(int j 1; j &…

RT-Thread 汇编分析启动流程

文章目录 一、汇编指令二、启动文件三、流程图 一、汇编指令 这里介绍即几条最常见实用的汇编指令 LDR R0,[R1]&#xff1a;将R1指定内存地址数据&#xff0c;存储到寄存器R0中。STR R0,[R1,#4]&#xff1a;将寄存器R0中数据存储到寄存器R1加上偏移量4的位置。MOV r0,#0x01&a…

read()之后操作系统都干了什么

首先说明三个参数 file文件 buff从内存中开辟一段缓冲区用来接收读取的数据 size表示这个缓冲区的大小 有关file的参数&#xff1a; 状态&#xff1a;被打开 被关闭权限&#xff1a;可读可写最重要的是inode: 他包含了 文件的元数据(比如文件大小 文件类型 文件在访问前需要加…

Matlab和python详解数独谜题问题

&#x1f517; 运行环境&#xff1a;Matlab、Python &#x1f6a9; 撰写作者&#xff1a;左手の明天 &#x1f947; 精选专栏&#xff1a;《python》 &#x1f525; 推荐专栏&#xff1a;《算法研究》 &#x1f510;#### 防伪水印——左手の明天 ####&#x1f510; &#x1f4…

c语言调用free,提示已触发了一个断点。

在用c语言写数据结构的链表的时候&#xff0c;执行也没有什么大错&#xff0c;逻辑也是对的&#xff0c;但是一道free函数会自动触发一个断点。如图&#xff1a; 这个断点产生的原因是由于分配的内存太小了在使用的时候没有任何问题&#xff0c;但是在执行程序的时候&#xff0…

JAVA毕业设计113—基于Java+Springboot+Vue的体育馆预约系统(源代码+数据库+12000字论文)

基于JavaSpringbootVue的体育馆预约系统(源代码数据库12000字论文)113 一、系统介绍 本项目前后端分离&#xff0c;本系统分为管理员、用户两种角色 用户角色包含以下功能&#xff1a; 注册、登录、场地(查看/预订/收藏/退订)、在线论坛、公告查看、我的预订管理、我的收藏…

[二分查找]LeetCode2040:两个有序数组的第 K 小乘积

本文涉及的基础知识点 二分查找算法合集 题目 给你两个 从小到大排好序 且下标从 0 开始的整数数组 nums1 和 nums2 以及一个整数 k &#xff0c;请你返回第 k &#xff08;从 1 开始编号&#xff09;小的 nums1[i] * nums2[j] 的乘积&#xff0c;其中 0 < i < nums1.…

45 - 多线程性能优化常见问题

1、使用系统命令查看上下文切换 上下文切换常见的监测工具 1.1、Linux 命令行工具之 vmstat 命令 vmstat 是一款指定采样周期和次数的功能性监测工具&#xff0c;我们可以使用它监控进程上下文切换的情况。 vmstat 1 3 命令行代表每秒收集一次性能指标&#xff0c;总共获取 …

代码随想录算法训练营 ---第五十三天

第一题&#xff1a; 简介&#xff1a; 本题和昨天的最大重复子串问题很相似&#xff0c;只不过本题不一定是连续的。 动规五部曲分析如下&#xff1a; 确定dp数组&#xff08;dp table&#xff09;以及下标的含义 dp[i][j]&#xff1a;长度为i-1 的字符串text1与长度为j-1的…