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=x3+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