2023年全國(guó)大學(xué)生高教杯數(shù)學(xué)建模C題古塔變形論文

2023高教社杯全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽承諾書我們仔細(xì)閱讀了?全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽章程?和?全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽參賽規(guī)那么?〔以下簡(jiǎn)稱為“競(jìng)賽章程和參賽規(guī)那么〞，可從全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽網(wǎng)站下載〕。我們完全明白，在競(jìng)賽開始后參賽隊(duì)員不能以任何方式〔包括、電子郵件、網(wǎng)上咨詢等〕與隊(duì)外的任何人〔包括指導(dǎo)教師〕研究、討論與賽題有關(guān)的問(wèn)題。我們知道，抄襲別人的成果是違反競(jìng)賽章程和參賽規(guī)那么的，如果引用別人的成果或其他公開的資料〔包括網(wǎng)上查到的資料〕，必須按照規(guī)定的參考文獻(xiàn)的表述方式在正文引用處和參考文獻(xiàn)中明確列出。我們鄭重承諾，嚴(yán)格遵守競(jìng)賽章程和參賽規(guī)那么，以保證競(jìng)賽的公正、公平性。如有違反競(jìng)賽章程和參賽規(guī)那么的行為，我們將受到嚴(yán)肅處理。我們授權(quán)全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽組委會(huì)，可將我們的論文以任何形式進(jìn)行公開展示〔包括進(jìn)行網(wǎng)上公示，在書籍、期刊和其他媒體進(jìn)行正式或非正式發(fā)表等〕。我們參賽選擇的題號(hào)是〔從A/B/C/D中選擇一項(xiàng)填寫〕：C我們的參賽報(bào)名號(hào)為〔如果賽區(qū)設(shè)置報(bào)名號(hào)的話〕：所屬學(xué)?！舱?qǐng)?zhí)顚懲暾娜常洪L(zhǎng)春工業(yè)大學(xué)參賽隊(duì)員(打印并簽名)：1.武太彬2.賈光芒3.牛文正指導(dǎo)教師或指導(dǎo)教師組負(fù)責(zé)人(打印并簽名)：李純潔 〔論文紙質(zhì)版與電子版中的以上信息必須一致，只是電子版中無(wú)需簽名。以上內(nèi)容請(qǐng)仔細(xì)核對(duì)，提交后將不再允許做任何修改。如填寫錯(cuò)誤，論文可能被取消評(píng)獎(jiǎng)資格?！橙掌冢?023年9月_16_日賽區(qū)評(píng)閱編號(hào)〔由賽區(qū)組委會(huì)評(píng)閱前進(jìn)行編號(hào)〕：2023高教社杯全國(guó)大學(xué)生數(shù)學(xué)建模競(jìng)賽編號(hào)專用頁(yè)賽區(qū)評(píng)閱編號(hào)〔由賽區(qū)組委會(huì)評(píng)閱前進(jìn)行編號(hào)〕：賽區(qū)評(píng)閱記錄〔可供賽區(qū)評(píng)閱時(shí)使用〕：評(píng)閱人評(píng)分備注全國(guó)統(tǒng)一編號(hào)〔由賽區(qū)組委會(huì)送交全國(guó)前編號(hào)〕：全國(guó)評(píng)閱編號(hào)〔由全國(guó)組委會(huì)評(píng)閱前進(jìn)行編號(hào)〕：古塔的變形摘要本文對(duì)古塔的變形問(wèn)題建立數(shù)學(xué)模型，它實(shí)質(zhì)上是一個(gè)空間解析幾何問(wèn)題。首先建立空間解析幾何模型，并利用這個(gè)模型對(duì)問(wèn)題1進(jìn)行求解，然后對(duì)模型進(jìn)行數(shù)據(jù)處理和圖形分析，精確地給出了中心位置坐標(biāo)；之后對(duì)于問(wèn)題2，在問(wèn)題1的根底上我們計(jì)算出中心坐標(biāo)的擬合直線，計(jì)算出曲率的值；最后對(duì)于問(wèn)題3，運(yùn)用AR自回歸模型給出古塔的變形趨勢(shì)。在問(wèn)題1中，通過(guò)分析這四年古塔的每一層中的8個(gè)離散點(diǎn)，運(yùn)用MATLAB數(shù)學(xué)軟件建立觀測(cè)數(shù)據(jù)模擬圖，依據(jù)此圖作出每一層相應(yīng)散點(diǎn)的投影平面，經(jīng)過(guò)計(jì)算根本近似正八邊形，得出正八邊形的中心，并運(yùn)用空間直線擬合模型，求出古塔的一條軸心擬合直線；并且運(yùn)用空間平面擬合模型，求出每一層的8個(gè)散點(diǎn)擬合的平面，那么直線與每一平面相交的點(diǎn)，即為每一層的中心坐標(biāo)，以所求第一層中心坐標(biāo)為例，1986年第一層中心坐標(biāo)為，1996年第一層中心坐標(biāo)為，2023年第一層中心坐標(biāo)為，2023年第一層中心坐標(biāo)，其它〔算上塔尖〕14層見正文表5。在問(wèn)題2中，首先運(yùn)用MATLAB軟件畫出每一年的俯視圖，即各年的平面圖，可以看出這四年的古塔是逐漸傾斜和彎曲的，且在第五層開始發(fā)生了扭曲。建立曲率數(shù)學(xué)模型，運(yùn)用軟件求解出相應(yīng)的曲率值。在問(wèn)題3中，首先建立帶季節(jié)項(xiàng)的AR自回歸模型，運(yùn)用問(wèn)題2中所求的曲率作為自變量，代入到模型中，運(yùn)用時(shí)間序列分析的統(tǒng)計(jì)軟件SPSS，得出古塔變形的趨勢(shì)，即隨年代的增長(zhǎng)，古塔的傾斜和彎曲將更加嚴(yán)重，且在第五層由于空間中心的改變，將更加扭曲。在本文最后，對(duì)模型的優(yōu)缺點(diǎn)及改良之處進(jìn)行分析。關(guān)鍵詞：空間解析幾何最小二乘法擬合曲率AR預(yù)測(cè)模型MATLAB7.0軟件1問(wèn)題的重述與分析由于許多外在原因，古塔會(huì)產(chǎn)生各種變形，傾斜，彎曲，扭曲等。為了了解各種變形量，測(cè)繪公司先后于1986年7月，1996年8月，2023年3月和2023年3月對(duì)該塔進(jìn)行了四次觀測(cè)。討論3個(gè)問(wèn)題，問(wèn)題1，給出確定古塔各層中心位置的通用方法，列表給出各次測(cè)量的古塔各層中心坐標(biāo)。本文首先對(duì)古塔的變形情況進(jìn)行分析，可以獲取位置信息，且只有四次的觀測(cè)數(shù)據(jù)信息。建立適當(dāng)?shù)淖鴺?biāo)系，先研究一層塔的八個(gè)點(diǎn)的投影所得的八個(gè)點(diǎn)的坐標(biāo)點(diǎn)，然后再確定各中心的坐標(biāo)，從而找到各層的中心點(diǎn)，可以擬合成一條直線。并且運(yùn)用空間擬合平面模型，求出每一層的8個(gè)散點(diǎn)擬合的平面，那么直線與每一平面相交的點(diǎn)，即為每一層的中心坐標(biāo)。問(wèn)題2，分析該塔傾斜，彎曲，扭曲等變形情況。首先運(yùn)用MATLAB軟件畫出每一年的俯視圖，即各年的平面圖，可以看出四年的古塔是傾斜還是彎曲或是扭曲，建立相應(yīng)的數(shù)學(xué)模型。問(wèn)題3，分析該塔的變形趨勢(shì)。首先建立帶季節(jié)項(xiàng)的AR自回歸模型，運(yùn)用問(wèn)題2中所求的結(jié)果作為自變量，代入到模型中，得出古塔變形的趨勢(shì)。2根本假設(shè)1.觀測(cè)的數(shù)據(jù)準(zhǔn)確；

2.以古塔的地面為水平面；

3.假設(shè)人在觀測(cè)時(shí)是圍繞塔的中心進(jìn)行八個(gè)不同的方面進(jìn)行大量的測(cè)量；4.兩點(diǎn)的距離作到小數(shù)點(diǎn)后百分位即為等長(zhǎng)，后面忽略；5.不考慮異常點(diǎn)；3符號(hào)說(shuō)明每個(gè)點(diǎn)的橫坐標(biāo)；每個(gè)點(diǎn)的縱坐標(biāo)；每個(gè)點(diǎn)的豎坐標(biāo)；直線方程的待估系數(shù)；空間直線的對(duì)應(yīng)的方向向量；：分別是邊和邊由其原始位置旋轉(zhuǎn)的角度；:變化率；4模型的建立與求解4.1問(wèn)題1的求解首先根據(jù)每一年數(shù)據(jù)中的每一層的8個(gè)散點(diǎn)坐標(biāo)運(yùn)用MATLAB軟件，得到如下每年觀測(cè)數(shù)據(jù)的模擬圖。圖11986年觀測(cè)數(shù)據(jù)模擬圖圖21996年觀測(cè)數(shù)據(jù)模擬圖圖32023年觀測(cè)數(shù)據(jù)模擬圖圖42023年觀測(cè)數(shù)據(jù)模擬圖由圖1到圖4的模擬圖可以直觀地看到，四次測(cè)量古塔的顯著變化，相同坐標(biāo)上塔的上部越來(lái)越小，把圖2和圖3疊加到一個(gè)圖。圖5古塔扭曲變化比照?qǐng)D圖5中，紅線代表1996年的塔圖，黃線代表2023年的塔圖。兩個(gè)塔圖疊加在一起比照更加顯明，而在第五層出現(xiàn)扭曲，從第五層往上越來(lái)越小，到2023年成了一個(gè)點(diǎn)的塔尖。經(jīng)過(guò)以上的圖形分析，本文以空間解析幾何為背景，運(yùn)用最小二乘擬合理論給出空間直線與平面的擬合。下面建立空間內(nèi)曲線的最小二乘問(wèn)題，以空間直線為例研究曲線擬合的方法?？臻g直線的標(biāo)準(zhǔn)方程整理得直線射影式方程，其中;;這樣直線可以看作是用這2個(gè)方程表示的平面相交的直線，所以可以分別對(duì)2個(gè)方程進(jìn)行數(shù)據(jù)擬合。設(shè)表示按擬合方程求得的近似值。一般地，它不同于實(shí)測(cè)值兩者之差，同理可得當(dāng)Q取最小值時(shí)的值即為方程的系數(shù)，即滿足以下方程時(shí)Q值最小有令(1)方程組(1)可寫成其中根據(jù)m組數(shù)據(jù)點(diǎn)解方程組就可以求得的值。圖61996年古塔各層中心點(diǎn)擬合直線分析圖圖72023年古塔各層中心點(diǎn)擬合函數(shù)分析圖從圖6，圖7中可以得出：由1996年古塔各層中心點(diǎn)擬合直線與2023年古塔各層中心點(diǎn)擬合直線，方程為，由一次方程的系數(shù)k看出1996年為0.0105,2023年為0.0112，可以得出古塔的中心發(fā)生了一次傾斜，且2023年比1986年傾斜角度變大。運(yùn)用附件中MATLAB的程序經(jīng)過(guò)調(diào)試，得到如下結(jié)果：1986年1996年2023年2023年接下來(lái)運(yùn)用最小二乘法擬合求空間平面方程：設(shè)空間平面方程為，其中為待估參數(shù)。設(shè)古塔的每層8個(gè)點(diǎn)，設(shè)點(diǎn)坐標(biāo)，考慮到數(shù)據(jù)在三個(gè)方向均存在誤差，得矩陣形式：即。通常采用矩陣奇異值分解解算待定參數(shù)的整體最小二乘解。其中那么參數(shù)的整體最小二乘估計(jì)為：殘差矩陣為：。運(yùn)用MATLAB軟件得到古塔各層平面方程，列表為：表11986年的擬合平面方程z=(0.003417)*x+(-0.000831)*y+(0.471956)z=(0.003629)*x+(-0.000818)*y+(5.887189)z=(0.003726)*x+(-0.000883)*y+(11.308475)z=(0.003838)*x+(-0.000935)*y+(15.602646)z=(0.003941)*x+(-0.000874)*y+(20.156847)z=(0.005443)*x+(-0.017498)*y+(33.310426)z=(0.005668)*x+(-0.018744)*y+(37.502066)z=(0.005887)*x+(-0.020009)*y+(41.620232)z=(0.006115)*x+(-0.021455)*y+(45.825845)z=(0.000887)*x+(-0.021682)*y+(52.003938)z=(0.001868)*x+(-0.022189)*y+(56.048295)z=(0.002175)*x+(-0.023884)*y+(61.121714)z=(0.004664)*x+(-0.027610)*y+(66.050359)表21996年的擬合平面方程z=(0.003715)*x+(-0.001488)*y+(0.684400)z=(0.003782)*x+(-0.000493)*y+(5.617203)z=(0.003735)*x+(-0.001266)*y+(11.515906)z=(0.004045)*x+(0.000715)*y+(14.555943)z=(0.003941)*x+(-0.001286)*y+(20.385763)z=(0.005639)*x+(-0.017135)*y+(32.997144)z=(0.005880)*x+(-0.018358)*y+(37.168214)z=(0.005877)*x+(-0.020489)*y+(41.891314)z=(0.006366)*x+(-0.021014)*y+(45.437916)z=(0.000873)*x+(-0.022225)*y+(52.314084)z=(0.001853)*x+(-0.022770)*y+(56.380505)z=(0.002159)*x+(-0.024511)*y+(61.481252)z=(0.005234)*x+(-0.027989)*y+(65.963401)表32023年的擬合平面方程z=(-0.003770)*x+(-0.002373)*y+(5.079985)z=(-0.000194)*x+(-0.003972)*y+(9.661222)z=(-0.003697)*x+(-0.002317)*y+(15.977842)z=(-0.000065)*x+(-0.004433)*y+(19.616420)z=(-0.000946)*x+(-0.004248)*y+(24.611974)z=(-0.020936)*x+(-0.004708)*y+(39.819633)z=(-0.018595)*x+(-0.006811)*y+(43.402884)z=(-0.021604)*x+(-0.006530)*y+(48.330004)z=(-0.022140)*x+(-0.006879)*y+(52.311658)z=(-0.022761)*x+(-0.001523)*y+(52.915196)z=(-0.022023)*x+(-0.002515)*y+(57.363039)z=(-0.023948)*x+(-0.002804)*y+(62.798782)z=(-0.027663)*x+(-0.005448)*y+(70.356472)表42023年的擬合平面方程z=(-0.003740)*x+(-0.002333)*y+(5.040576)z=(-0.002530)*x+(-0.002555)*y+(10.061014)z=(-0.004408)*x+(-0.002633)*y+(16.522721)z=(-0.002736)*x+(-0.003593)*y+(20.518390)z=(-0.001351)*x+(-0.004259)*y+(24.823896)z=(-0.021377)*x+(-0.004712)*y+(40.045928)z=(-0.018205)*x+(-0.007014)*y+(43.305695)z=(-0.021560)*x+(-0.006476)*y+(48.273096)z=(-0.025300)*x+(-0.004743)*y+(52.729697)z=(-0.026175)*x+(0.000449)*y+(53.562983)z=(-0.024393)*x+(-0.002167)*y+(58.395736)z=(-0.018681)*x+(-0.004368)*y+(60.918705)z=(-0.028368)*x+(-0.005161)*y+(70.556460)聯(lián)立直線與平面方程，即得到表5的各年各層的中心坐標(biāo)。表5各年各層的中心坐標(biāo)1986年1996年2023年2023年1層2層3層4層5層6層7層8層9層10層11層12層13層4.2問(wèn)題2的求解問(wèn)題2是為了分析古塔的變形情況，根據(jù)問(wèn)題1的圖1-圖5，運(yùn)用MATLAB軟件，得到如下關(guān)于軸的俯視圖：圖8古塔模擬俯視圖圖91986年觀測(cè)數(shù)據(jù)俯視模擬圖圖101996年觀測(cè)數(shù)據(jù)模擬俯視圖圖112023年觀測(cè)數(shù)據(jù)模擬俯視圖圖8至圖11分別給出了整體和各個(gè)年份的俯視圖，從四個(gè)圖中可以看到并不是均勻的八邊形環(huán)，有些邊環(huán)密，有些邊環(huán)稀，這就直觀地說(shuō)明了古塔發(fā)生了傾斜和彎曲，加上問(wèn)題1的第五層的中心坐標(biāo)的改變，說(shuō)明古塔在第五層上也發(fā)生了扭曲現(xiàn)象。本文在第一問(wèn)最后得出每層各中心坐標(biāo)的擬合直線，運(yùn)用到第二問(wèn)中古塔的傾斜問(wèn)題上，傾斜是指根底兩端點(diǎn)傾斜方向的沉降差與其距離的比值建筑中心線或其墻、柱等，在不同高度的點(diǎn)對(duì)其相應(yīng)底部點(diǎn)的偏移現(xiàn)象。運(yùn)用公式兩條直線夾角公式：彎曲，即不直?？梢苑譃樾巫儚澢翱臻g彎曲。當(dāng)桿件受到與桿軸線垂直的外力或在軸線平面內(nèi)的力偶作用時(shí)，桿的軸線由原來(lái)的直線變成曲線，這種變形叫彎曲變形。曲率處處不為零的空間稱為彎曲空間。物體發(fā)生彎曲時(shí)產(chǎn)生的形變叫做“彎曲形變〞。物體彎曲得越厲害，產(chǎn)生的彈力就越大。例如，將弓拉得越滿，箭就射得越遠(yuǎn)。把一個(gè)物體放在支持物上，物體越重，支持物被壓彎曲得越厲害，支持力就越大。運(yùn)用曲率公式：扭曲物體因外力作用而扭轉(zhuǎn)變形,也用于比喻數(shù)學(xué)的一種改變物體形狀的方法，扭曲變形是在彎曲變形的根底上，旋轉(zhuǎn)，扭成螺旋狀態(tài)，叫扭曲變形。問(wèn)題1中得出扭曲的發(fā)生主要在第五層。首先,計(jì)算方向的扭曲變形。由于邊和的旋轉(zhuǎn)角度一般都很小(對(duì)于一個(gè)實(shí)際問(wèn)題,一般因此可以認(rèn)為約等于,因此:式中:分別是邊和邊由其原始位置旋轉(zhuǎn)的角度。如果，那么說(shuō)明微元的一邊產(chǎn)生均勻豎向位移,那么變形后各點(diǎn)仍保持在一個(gè)平面內(nèi),意味著這個(gè)微元發(fā)生的是剛性的旋轉(zhuǎn)。否那么,當(dāng),即在兩條邊之間存在一個(gè)旋轉(zhuǎn)角的差，就說(shuō)明產(chǎn)生了扭曲變形。那么在變形后的面上,扭曲變形的大小可以采用這兩個(gè)角度之差隨著兩個(gè)對(duì)邊之間距離的變化率來(lái)表示:或者:如果應(yīng)用以上分析過(guò)程考慮y方向的扭曲變形,即首先計(jì)算微元另外兩條對(duì)邊和的旋轉(zhuǎn)角度的差,將會(huì)得到相同的結(jié)果。最終運(yùn)用MATLAB軟件得出4.3問(wèn)題3的求解在問(wèn)題1的根底上，運(yùn)用MATLAB軟件，給出1996年與2023年的中心擬合直線圖來(lái)。時(shí)間序列分析的目標(biāo)就是通過(guò)分析要素〔變量〕隨時(shí)間變化的歷史過(guò)程，揭示其變化開展規(guī)律，并對(duì)未來(lái)狀態(tài)進(jìn)行分析預(yù)測(cè)。如在變形測(cè)量中，可以采用時(shí)間序列分析方法對(duì)觀測(cè)數(shù)據(jù)進(jìn)行分析，以便建立變形體的動(dòng)態(tài)變形預(yù)測(cè)模型，并對(duì)其變形趨勢(shì)進(jìn)行預(yù)測(cè)。而自回歸(AR)模型的參數(shù)估計(jì)是時(shí)間序列分析的根本問(wèn)題，是在模型結(jié)構(gòu)及階次已確定的條件下，對(duì)模型參數(shù)進(jìn)行估計(jì)，使所建立的模型是實(shí)際時(shí)間序列的“最正確〞擬合模型。但在實(shí)際的觀測(cè)中，觀測(cè)值是由一定觀測(cè)手段得到的，不可防止地含有隨機(jī)誤差。理論上講AR(p)參數(shù)的最小二乘估計(jì)也就是線性最小二乘估計(jì)，所以可以用現(xiàn)在廣泛應(yīng)用的整體最小二乘估計(jì)方法進(jìn)行估計(jì)，即不僅考慮自身觀測(cè)值的誤差，同時(shí)考慮與其有關(guān)的自身前一個(gè)或前幾個(gè)時(shí)刻的觀測(cè)值的誤差，從而進(jìn)行參數(shù)估計(jì)，得到的預(yù)測(cè)觀測(cè)數(shù)據(jù)更準(zhǔn)確。但AR模型是典型的預(yù)測(cè)模型，而最小二乘估計(jì)是常用的數(shù)據(jù)分析工具，在預(yù)測(cè)數(shù)據(jù)上還是常應(yīng)用最小二乘估計(jì)比擬準(zhǔn)確，反而最小二乘估計(jì)的數(shù)據(jù)不夠準(zhǔn)確和穩(wěn)定。自回歸模型，子樣觀察值，白噪聲序列表示為，回歸系數(shù)用表示，那么可以得到模型：，模型參數(shù)的最小二乘估計(jì)。設(shè)樣本觀測(cè)值，記,,那么模型可表示為有最小二乘原理可得到模型參數(shù)的估計(jì)為。那么根據(jù)最小二乘估計(jì)值可以得到噪聲的估計(jì)值為噪聲方差為，由此模型可以求出未來(lái)的變形趨勢(shì)。5．模型的評(píng)價(jià)及推廣5.1模型的優(yōu)點(diǎn)(1).在問(wèn)題的求解中，充分運(yùn)用了數(shù)據(jù)和圖形，即表格，使結(jié)果明了清晰。(2).本文采用了多種專業(yè)軟件對(duì)模型進(jìn)行求解，如Matlab，SPSS等提高了模型的精確度。(3).本文所有模型建立均完全基于實(shí)際的統(tǒng)計(jì)數(shù)據(jù)，擬合成近似的線進(jìn)行處理，科學(xué)合理。(4).本文建立多種模型，多種模型的比照更能表達(dá)不同模型的特點(diǎn)及優(yōu)勢(shì)。(5).應(yīng)用相關(guān)幾何知識(shí)，將空間復(fù)雜的點(diǎn)分布簡(jiǎn)化成直線和面，應(yīng)用幾何問(wèn)題求空間點(diǎn)。5.2模型的缺點(diǎn)(1).背景資料的篩選方法有待進(jìn)一步優(yōu)化和改良；(2).模型建立中采用的數(shù)據(jù)不太準(zhǔn)確，利用擬合估計(jì)值，對(duì)模型產(chǎn)生了相對(duì)較大的誤差。(3).建立的模型不夠完善，忽略了許多的因素，過(guò)于簡(jiǎn)單。(4).由于古塔在現(xiàn)實(shí)生活中受到各種影響，理論值難免與實(shí)際情況有所偏差。5.3模型的改良由于都是求出每個(gè)點(diǎn)都是一個(gè)估計(jì)值，從而有較大的誤差，需要在模型的建立中進(jìn)行改良，使估計(jì)的所有的數(shù)據(jù)更加準(zhǔn)確。5.4模型的推廣用中心位置的模型方法，不難預(yù)測(cè)全國(guó)所有古塔或高樓由于外力和內(nèi)力的作用在以后的變化趨勢(shì)。6.參考文獻(xiàn)[1]姚澤清,鄭旭東,全國(guó)大學(xué)生數(shù)學(xué)建模賽題優(yōu)秀論文評(píng)析<2005年—2023年A題>,北京:國(guó)防工業(yè)出版社，2023:198-218.[2]周培德，計(jì)算幾何[M].北京：清華大學(xué)出版社，2023：100-105.[3]襲楊，空間直線擬合的一種方法[J].齊齊哈爾大學(xué)學(xué)報(bào)，2023,25(2):64-68.[4]王能超.數(shù)值分析簡(jiǎn)明教程[M].北京：高等教育出版社，2003：36-37.[5]呂林根，許子道.解析幾何[M].北京：高等教育出版社，2004：156-157.[6]于曉秋，李欣，張宏禮.MATLAB數(shù)值實(shí)驗(yàn)[M].沈陽(yáng)：遼寧科學(xué)技術(shù)出版社，2004：34-40.[7]王正東.數(shù)學(xué)軟件與數(shù)學(xué)實(shí)驗(yàn)[M].北京：科學(xué)出版社，2004：88-89.[8]吳懷宇．時(shí)間序列分析與綜合[M]．武漢:武漢大學(xué)出版社，2004．[9]俞錦成．關(guān)于整體最小二乘的可解性[J]．南京師范大學(xué)學(xué)報(bào)(自然科學(xué)版)，1996，19(1):13-16．7附錄直線擬合：>>x=[566.6648566.7196566.7735566.8161566.8621566.9084566.9467566.9843567.0218567.0569567.1045567.1518];y=[522.7105522.6684522.6273522.5944522.5591522.5244522.5081522.4924522.4764522.4624522.423522.3836];z=[1.78747.320212.755217.078321.720526.235129.836933.350936.854940.172144.440948.7119];F=[z;111111111111];M=F*F';N=F*x';O=F*y';A=(M\N)'B=(M\O)'A=0.0104566.6411B=-0.0066522.7125//1986年>>x=[566.665566.7205566.7751566.8183566.8649566.9118566.9506566.9884567.0265567.062567.1102567.1578];y=[522.7102522.6674522.6256522.5922522.5563522.521522.5042522.4881522.4714522.4572522.4173522.3775];z=[1.7837.314612.750717.075121.71626.229529.832333.345436.848340.167644.435448.7074];F=[z;111111111111];M=F*F';N=F*x';O=F*y';A=(M\N)'B=(M\O)'A=0.0105566.6411B=-0.0067522.7124//1996年>>x=[566.727566.7642566.8004566.8297566.861566.9478566.98567.0313567.0825567.1381567.181567.2238];y=[522.7014522.669522.6387522.6127522.586522.5335522.5115522.4788522.4457522.3926522.3535522.3147];z=[1.76327.290512.726917.05221.703926.204529.81733.336636.822240.144144.424948.6839];F=[z;111111111111];M=F*F';N=F*x';O=F*y';A=(M\N)'B=(M\O)'A=0.0112566.6643B=-0.0084522.7436//2023年>>x=[566.7268566.764566.8001566.8293566.8603566.9471566.9792567.0305567.0816567.137567.1799567.2225];y=[522.7015522.6693522.6384522.6132522.5866522.5342522.5123522.4797522.4466522.3937522.3547522.316];z=[1.76457.30912.732317.069721.709426.21129.824633.339936.843840.161144.432648.6998];F=[z;111111111111];M=F*F';N=F*x';O=F*y';A=(M\N)'B=(M\O)'A=0.0112566.6642B=-0.0084522.7436//2023年平面擬合公式：86.1z=(0.003417)*x+(-0.000831)*y+(0.471956)z=(0.003629)*x+(-0.000818)*y+(5.887189)z=(0.003726)*x+(-0.000883)*y+(11.308475)z=(0.003838)*x+(-0.000935)*y+(15.602646)z=(0.003941)*x+(-0.000874)*y+(20.156847)z=(0.005443)*x+(-0.017498)*y+(33.310426)z=(0.005668)*x+(-0.018744)*y+(37.502066)z=(0.005887)*x+(-0.020009)*y+(41.620232)z=(0.006115)*x+(-0.021455)*y+(45.825845)z=(0.000887)*x+(-0.021682)*y+(52.003938)z=(0.001868)*x+(-0.022189)*y+(56.048295)z=(0.002175)*x+(-0.023884)*y+(61.121714)z=(0.004664)*x+(-0.027610)*y+(66.050359)96.1z=(0.003715)*x+(-0.001488)*y+(0.684400)z=(0.003782)*x+(-0.000493)*y+(5.617203)z=(0.003735)*x+(-0.001266)*y+(11.515906)z=(0.004045)*x+(0.000715)*y+(14.555943)z=(0.003941)*x+(-0.001286)*y+(20.385763)z=(0.005639)*x+(-0.017135)*y+(32.997144)z=(0.005880)*x+(-0.018358)*y+(37.168214)z=(0.005877)*x+(-0.020489)*y+(41.891314)z=(0.006366)*x+(-0.021014)*y+(45.437916)z=(0.000873)*x+(-0.022225)*y+(52.314084)z=(0.001853)*x+(-0.022770)*y+(56.380505)z=(0.002159)*x+(-0.024511)*y+(61.481252)z=(0.005234)*x+(-0.027989)*y+(65.963401)99.1z=(-0.003770)*x+(-0.002373)*y+(5.079985)z=(-0.000194)*x+(-0.003972)*y+(9.661222)z=(-0.003697)*x+(-0.002317)*y+(15.977842)z=(-0.000065)*x+(-0.004433)*y+(19.616420)z=(-0.000946)*x+(-0.004248)*y+(24.611974)z=(-0.020936)*x+(-0.004708)*y+(39.819633)z=(-0.018595)*x+(-0.006811)*y+(43.402884)z=(-0.021604)*x+(-0.006530)*y+(48.330004)z=(-0.022140)*x+(-0.006879)*y+(52.311658)z=(-0.022761)*x+(-0.001523)*y+(52.915196)z=(-0.022023)*x+(-0.002515)*y+(57.363039)z=(-0.023948)*x+(-0.002804)*y+(62.798782)z=(-0.027663)*x+(-0.005448)*y+(70.356472)11.1z=(-0.003740)*x+(-0.002333)*y+(5.040576)z=(-0.002530)*x+(-0.002555)*y+(10.061014)z=(-0.004408)*x+(-0.002633)*y+(16.522721)z=(-0.002736)*x+(-0.003593)*y+(20.518390)z=(-0.001351)*x+(-0.004259)*y+(24.823896)z=(-0.021377)*x+(-0.004712)*y+(40.045928)z=(-0.018205)*x+(-0.007014)*y+(43.305695)z=(-0.021560)*x+(-0.006476)*y+(48.273096)z=(-0.025300)*x+(-0.004743)*y+(52.729697)z=(-0.026175)*x+(0.000449)*y+(53.562983)z=(-0.024393)*x+(-0.002167)*y+(58.395736)z=(-0.018681)*x+(-0.004368)*y+(60.918705)z=(-0.028368)*x+(-0.005161)*y+(70.556460)平面擬合：>>x=xlsread('D:\data.xls','C4:C11');y=xlsread('D:\data.xls','D4:D11');z=xlsread('D:\data.xls','E4:E11');plot3(x,y,z);scatter3(x,y,z,'filled');holdon;X=[ones(8,1)xy];b=regress(z,X);holdon;xfit=min(x):0.1:max(x);yfit=min(y):0.1:max(y);[XFIT,YFIT]=meshgrid(xfit,yfit);ZFIT=b(1)+b(2)*XFIT+b(3)*YFIT;mesh(XFIT,YFIT,ZFIT);title(sprintf('z=(%f)*x+(%f)*y+(%f)',b(3),b(2),b(1)));兩直線夾角：functiona=JiaJiao(x,y)

%求兩條直線夾角

%x,y是三點(diǎn)的橫坐標(biāo)和縱坐標(biāo)

%eg:x=[123];y=[415];

ifx(2)~=x(1)

k1=(y(2)-y(1))/(x(2)-x(1));

end

ifx(3)~=x(2)

k2=(y(3)-y(2))/(x(3)-x(2));

end

ifx(2)==x(1)&x(3)==x(2)

a=0;

elseifx(3)==x(2)

a=pi/2-atan(abs(k1));

elseifx(1)==x(2)

a=pi/2-atan(abs(k2));

elseif1+k1*k2==0

a=pi/2;

else

a=atan(abs((k2-k1)/(1+k2*k1)));%夾角

end

a=a*360/(2*pi);%轉(zhuǎn)化為角度制向量的夾角：A=[129];

B=[010];

acos(dot(A,B)/(norm(A)*norm(B)))*180/pi曲率:

clc;clearall;closeall;

x0=linspace(0,1);

y0=sin(x0).*cos(x0);

h=abs(diff([x0(2),x0(1)]));

%模擬一階導(dǎo)

figure;boxon;holdon;

ythe1=cos(x0).^2-sin(x0).^2;%理論一階導(dǎo)

yapp1=gradient(y0,h);%matlab數(shù)值近似

plot(x0,ythe1,'.');

plot(x0,yapp1,'r');

legend('理論值','模擬值');

title('模擬一階導(dǎo)');

%模擬二階導(dǎo)

figure;boxon;holdon;

ythe2=(-4)*cos(x0).*sin(x0);%理論二階導(dǎo)

yapp2=2*2*del2(y0,h);%matlab數(shù)值近似

plot(x0,ythe2,'.');

plot(x0,yapp2,'r');

legend('理論值','模擬值');

title('模擬二階導(dǎo)');

%模擬曲率

symsxy

y=sin(x)*cos(x);

yd2=diff(y,2);

yd1=diff(y,1);

k=abs(yd2)/(1+yd1^2)^(3/2);

k1=subs(k,x,