2024年华中杯数学建模
C题 基于光纤传感器的平面曲线重建算法建模
原题再现
光纤传感技术是伴随着光纤及光通信技术发展起来的一种新型传感器技术。它是以光波为传感信号、光纤为传输载体来感知外界环境中的信号，其基本原理是当外界环境参数发生变化时，会引起光纤传感器中光波参量（如波长、相位、强度等）的变化，即外界信号变化会对光信号产生调制。光纤传感器具有质地轻、体积小、弯曲性能好，抗电磁干扰能力强，灵敏度高，易于安装使用等优点。光纤传感技术最重要的是实时获得结构实时应变信息，再通过解调出来的应变参数来重构得到结构的形变或位移。光纤传感器已在许多领域有实际应用，比如能够对结肠部位进行形状重建等。通过光纤传感器解调系统解调出来的应变信息，间接求出曲率等信息，并基于离散曲率信息对曲线进行重构。
为了便于波长测量，在生产光纤时，已在等间距位置布设好传感器，本次传感器间距为0.6米。在测量时，先在光纤水平状态（即初始状态如图1所示）测量各个传感器位置处信号的波长，然后在受到外力后（如图2所示）测量各个传感器位置处信号的波长。波长λ与曲线曲率k之间的关系近似为k=c(λ-λ0)/λ0,其中λ0是水平光纤在初始状态下测量的波长，λ是光纤在受到外力后测量的波长，c为某个常数，这里假设为4200。本次实验分别测量了两组不同初始状态下受力前后的波长值，具体数据见表1，并请解决如下问题。
问题1.请根据表1给出的波长测量数据，构建数学模型，估算平面光栅各个传感点(FBG1-FBG6)的曲率。进一步，假设初始点坐标为原点，初始的水平光纤方向为x轴，垂直方向为y轴，光纤在平面内受力后在初始位置的切线与水平方向的夹角为45o，请建立模型估算下列表格中横坐标x轴相应位置处的曲率。
问题2.请根据表1波长测量数据和问题1求出的曲率，构建数学模型，分别重构平面曲线，并分析曲线的特点。
问题3.请根据平面曲线方程y=x^3+x(0≤x≤1)，以适当的等间距弧长采样，计算这些采样点的曲率。然后以采样的曲率为基础，构建数学模型，重构平面曲线，并分析重构曲线与原始曲线出现误差的原因。
整体求解过程概述(摘要)
随着光纤通信技术在人们日常生活中的进一步普及，光纤传感技术作为一种新型的传感器技术也正逐渐发展起来。光纤传感器主要依靠光波作为传输信号、光纤作为传输载体来感知外部环境信号。由于光纤传感器具有质地轻、体积小、易弯曲、灵敏度高、抗电磁能力强的优点，其在医学领域上已经得到了较为广泛的运用。因此通过传感器返回数据，获得结构实时应变信息，再通过解调出来的应变参数来重构得到结构的形变或位移对于光纤传感器的广泛运用至关重要。
本文主要研究了依据曲率进行平面曲线重构的问题。根据Frenet框架在二维条件下的简化模型，利用欧拉法、龙格-库塔法等算法实现了在已知离散点曲率的条件下完成了曲线的重构，并分析了模型在多方面上的误差。
针对问题一光纤上各点曲率的估算，本文首先论证了在给定条件下曲率与反射波长间的线性关系，并根据题目所给关系式，通过传感器读数估计曲线上特定传感点曲率依次为2.219、2.217、2.233、2.230、2.236、2.222；2.986、2.978、2.973、2.981、2.984、2.975，然后分析得到了样条插值法的优越性，借助三次样条插值法重新构建了光纤曲线，再通过该重构曲线求解特定横坐标下曲线的估计曲率值，如文中图1所示。
针对问题二根据离散点处的曲率进行曲线的重构，本文首先介绍了在三维曲线重构中广泛应用的Frenet框架，并推导出其退化至二维的一些基础算法。再利用其衍生得到的小步长欧拉法构造重构曲线。借助局部放大分析曲线类圆的特点，并绘制重构曲线的曲率图像验证该曲线类似螺旋线的走势，且保持总长相等的物理性质。得出曲线特点后本文又借助多项式拟合法再次进行曲线拟合，并将两种算法得到的曲线及曲率放在同一坐标系中对比，并无太大差异，从而证实了模型的可靠性。
针对问题三对已知三次多项式的曲率采样与曲线重构，本文首先给出了等间距弧长采样的优点和实现过程。对已知三次曲线进行采样得到图7之后又分别采取前向欧拉法、四阶龙格-库塔法以及有限差分法完成了平面曲线的重构，并从不同类型误差的角度综合分析了这三种不同算法重构曲线的效果。以此为启发本文又探究了模拟退火法在本题目中的应用及其拟合结果，其与目标曲线高度重合的结果验证了该算法的优越性。
模型假设：
(1) 本文假设波长λ和曲线曲率κ在题目考查范围内始终满足线性关系式(1)，且关系式(1)中的c视作常数4200。
(2) 本文假设在对获得的离散的传感点曲率进行曲线拟合时，目标曲线应始终保持较高的平滑度，从而可以使用三次样条差分法进行基本的曲线拟合。
(3) 本文假设光纤上的传感器在水平状态和受到外力作用时都能保持稳定的等间距，且始终能保持返回该位置的波长值。
(4) 本文假设在对取样点曲率进行三次样条插值时，函数可在采样点之间用三次多项式近似表示。
问题分析：
问题一：这一问根据题目给出波长λ和曲线曲率κ之间的近似关系κ = c(λ−λ0)/λ0可将表1中各个传感器测得的波长数据转化为平面光栅各传感点的曲线曲率。针对第二小问，在给定起始点位置、方向及切线方向的情况下，本文将依据第一小问得到的各传感点曲线斜率，通过三次样条曲率插值法拟合曲线，并根据此拟合曲线估算出题目要求特定横坐标x轴相应位置曲率。
问题二：这一问需要综合考虑题目所给的各传感点波长以及第一问求得的传感点曲线曲率，构建数学模型获得光栅的平面重构曲线。本文将使用欧拉法对第一问获得的离散点曲率进行曲线重构，借助Matlab将测试1和测试2的两组拟合曲线展示在图表中，并从平滑度、趋势度、拟合度等多方面分析曲线特点。
问题三：这一问需首先对已知曲线y=x³+x进行适当的等间距弧长采样，获取多个离散曲率数据后，再基于这些数据，通过龙格-库塔法或其他算法重构一条拟合曲线。借助Matlab 工具将拟合曲线与给定曲线展示在一张图表上，观察重构曲线产生的误差并分析其产生的原因.
模型的建立与求解整体论文缩略图
全部论文请见下方“ 只会建模 QQ名片” 点击QQ名片即可
部分程序代码：
%q1.根据传感器读数返回曲率值
%插值法求解特定点的曲率值
c = 4200;%定义常数c%测试1的数据
lambda_0_test1=[1529, 1529,1529,1529,1529, 1529];lambda_test1=[1529.808, 1529.807,1529.813,1529.812,1529.814,1529.809];%测试2的数据
lambda_0_test2=[1540, 1540,1540,1540,1540, 1540];lambda_test2=[1541.095, 1541.092,1541.090,1541.093,1541.094,1541.091];%传感器位置
positions =[0, 0.6,1.2, 1.8,2.4,3.0]; %传感器位置
%计算测试1的曲率
kappa_test1= c * (lambda_test1-lambda_0_test1)./ lambda_0_test1;%创建插值函数（测试1）
curvature_interpolation_test1=fit(positions', kappa_test1','cubicinterp');%计算测试2的曲率
kappa_test2= c * (lambda_test2-lambda_0_test2)./ lambda_0_test2;%创建插值函数（测试2）
curvature_interpolation_test2=fit(positions', kappa_test2','cubicinterp');%指定横坐标位置
x_positions= [0.3, 0.4,0.5,0.6,0.7];%在这些位置估算曲率（测试1和测试2）
estimated_curvatures_test1=curvature_interpolation_test1(x_positions);estimated_curvatures_test2=curvature_interpolation_test2(x_positions);%输出估算的曲率（测试1）
disp('测试1的估算曲率值为：');for idx= 1:length(x_positions)fprintf('在x= %.1f米处的曲率:%.4f\n', x_positions(idx),estimated_curvatures_test1(idx));end%输出估算的曲率（测试2）
disp('测试2的估算曲率值为：');for idx= 1:length(x_positions)fprintf('在x= %.1f米处的曲率:%.4f\n', x_positions(idx),estimated_curvatures_test2(idx));end
%q2.1根据问题一求得结果重构完整曲线
c = 4200;%定义常数c和传感器数据
%测试1的数据
lambda_0_test1=[1529, 1529,1529,1529,1529, 1529];lambda_test1=[1529.808, 1529.807,1529.813,1529.812,1529.814,1529.809];positions =[0, 0.6,1.2, 1.8,2.4,3.0]; %传感器位置
%计算每个传感器点的曲率
kappa_test1= c * (lambda_test1-lambda_0_test1)./ lambda_0_test1;%测试2的数据
lambda_0_test2=[1540, 1540,1540,1540,1540, 1540];lambda_test2=[1541.095, 1541.092,1541.090,1541.093,1541.094,1541.091];%计算每个传感器点的曲率
kappa_test2= c * (lambda_test2-lambda_0_test2)./ lambda_0_test2;%设置步长
ds=0.01;%重构测试1和测试2的曲线
[x_points_test1,y_points_test1]=reconstruct_curve(positions,kappa_test1,ds,pi/4);[x_points_test2,y_points_test2]=reconstruct_curve(positions,kappa_test2,ds,pi/4);%绘制测试1和测试2的曲线
figure;plot(x_points_test1,y_points_test1,'LineWidth',2);holdon;plot(x_points_test2,y_points_test2,'r-','LineWidth',2);title('ReconstructedCurves forTest1andTest 2');xlabel('x');ylabel('y');legend('Test1','Test2');gridon;holdoff;%重构函数定义，增加初始角度参数
function [x_points, y_points]=reconstruct_curve(positions,kappa,ds,initial_angle)%预分配数组
total_steps = sum(round(diff(positions)/ ds));x_points=zeros(1, total_steps+1);y_points=zeros(1, total_steps+1);%初始条件，使用给定的初始角度
theta=initial_angle;%初始角度为pi/4，对应斜率为1x=0;y=0;index=1;%存储第一个点
x_points(index) = x;y_points(index) = y;fori =2:length(positions)steps= round((positions(i)-positions(i-1))/ ds);forj= 1:stepstheta=theta+kappa(i-1)*ds;x=x+cos(theta)*ds;y=y+sin(theta)*ds;index=index+1;x_points(index)= x;y_points(index)= y;endendend
