測量程序設(shè)計三課件

測量程序設(shè)計曾喆中國石油大學地球科學與技術(shù)學院測量程序設(shè)計曾喆中國石油大學地球科學與技術(shù)學院第三講法方程式的解算及系數(shù)陣的求逆主要內(nèi)容：高斯約化法高斯-約旦法變量循環(huán)重新編號法求逆平方根法求法方程式系數(shù)陣的逆改正的平方根法求法方程式系數(shù)陣的逆補充內(nèi)容第三講法方程式的解算及系數(shù)陣的求逆主要內(nèi)容：平差法方程式條件平差法方程間接平差法方程平差法方程式條件平差法方程1.高斯約化法主要步驟1、消元2、回代求解1.高斯約化法主要步驟實例解方程組第一步,將第一行乘-2加與第二行相加；第一行乘-1加到第三行,得到實例解方程組第二步,將上面第二行乘-2/3加到第三行,得到回代:解第三行得x3,將x3代入第二行得x2,將x2,x3代入第一行得x1,得到解x*=(2,1,-1)T第二步,將上面第二行乘-2/3加到第三行,得到設(shè)方程組為：其初始的矩陣形式為：設(shè)方程組為：其初始的矩陣形式為：如果，我們首先保留矩陣的第一行，并利用它來消去其余三行中的第一列。即

消元過程如果，我們首先保留矩陣的第一行，并利用第i行減去乘以第0行，得到其中，第i行減去乘以第0行，得到同理，如果保留第0行和第1行基礎(chǔ)上消去第1列對角線下方元素。即：同理，如果保留第0行和第1行基礎(chǔ)上消去第1第i行減去乘以第1行，得到其中，第i行減去乘以第1行，得到同理，若，還可以進一步消元同樣地，第3行減去乘以第2行，得到同理，若，還可以進一步消元其中，經(jīng)過這個消元過程，最開始的方程組就變?yōu)橐陨蠈⒃嫉姆匠探M化為上式的過程為消元過程。下面將從方程組的最下一開始上推出所有x。其中，

回代過程如果，根據(jù)下式回代得到所有的x的解回代過程如果，根據(jù)下式回代得到所有的x的小結(jié)從上面可以看出整個高斯消去法解方程的方法主要分為兩步1、消去過程

從k=0開始，在第k步中，定義因子

然后通過下式計算k+1步的方程系數(shù)

通過n-1步，即k=n-2時就可以得到原方程組系數(shù)的上三角陣

小結(jié)從上面可以看出整個高斯消去法解方程的方法主要分為兩步2、回代過程2、回代過程代碼//消元過程doublem;for(intk=0;k<n-1;k++){try{for(inti=k+1;i<n;i++){m=a[i,k]/a[k,k];for(intj=k+1;j<n;j++)a[i,j]=a[i,j]-a[k,j]*m;c[i]=c[i]-c[k]*m;}}catch(DivideByZeroExceptione){Console.WriteLine(e.Message);}}代碼//消元過程/*回代過程*/

x[n-1]=c[n-1]/a[n-1][n-1];for(i=n-2;i>=0;i--){doublesum=0.0;for(j=i+1;j<n;j++)sum=sum+a[i][j]*x[j];x[i]=(c[i]-sum)/a[i][i];}/*回代過程*/2.高斯-約旦法高斯消去法在消元時始終消去對角線下方的元素，而高斯——約當消去法則同時消去對角線上方和下方的元素。2.高斯-約旦法高斯消去法在消元時始終消去對角線下方的元素，第一次消元第二次消元第一次消元第二次消元測量程序設(shè)計三ppt課件故方程組的解為高斯——約當消去法的特點：(1)消元和回代同時進行；(2)乘除法的次數(shù)要比高斯消去法大，所以通常用于同時求解系數(shù)矩陣相同的多個方程組或求逆矩陣。故方程組的解為高斯——約當消去法的特點：(1)消元和回代求逆：則求逆：則不為0則上面方程組第一式除之,得到,交換與并將代入余下各式,得到:反復(fù)運行,可以將x全部放到左邊,y放到右邊,這樣y的系數(shù)陣就是C的逆陣不為0則上面方程組第一式除之,得到,交換計算公式:即計算公式:即3.變量循環(huán)重新編號法求逆方程式逆關(guān)系這里就有3.變量循環(huán)重新編號法求逆方程式逆關(guān)系這里就有通過第一式計算出x0，再代入余下的方程式可得到下式通過第一式計算出x0，再代入余下的方程式可得到下式新系數(shù)通過如下計算可得新系數(shù)通過如下計算可得對上式按下面規(guī)則重新編號經(jīng)過n次這種變換，變量均恢復(fù)原來的狀態(tài)對上式按下面規(guī)則重新編號經(jīng)過n次這種變換，變量均恢復(fù)原來的狀按上法進行了k步，A為正定對稱，就有

1、對于任意k，A11和A22均為正定2、A12=-A21T按上法進行了k步，A為對稱正定矩陣，用“變量循環(huán)重新編號法”經(jīng)過n次變換得到A-1，計算如下：A為對稱正定矩陣，用“變量循環(huán)重新編號法”經(jīng)過n次變換得到A主要代碼for(intk=0;k<n;k++){doublea00=a[0];if(a00+1.0==1.0) { delete[]a0; returnfalse; }for(inti=0;i<n;i++){doubleai0=a[i*(i+1)/2];if(i<=n-k-1)a0[i]=-ai0/a00;elsea0[i]=ai0/a00;for(intj=1;j<=i;j++){a[(i-1)*i/2+j-1]=a[i*(i+1)/2+j]+ai0*a0[j];}for(i=1;i<n;i++){a[(n-1)*n/2+i-1]=a0[i];}a[n*(n+1)/2-1]=1.0/a00;}}主要代碼4.平方根法A是正定對稱陣，則有例：根據(jù)矩陣乘法有：4.平方根法A是正定對稱陣，則有由此可得：由此可得：計算公式為計算公式為即即下三角陣LR也為下三角陣,(i>=k>=j)下三角陣這樣，A的逆就為這樣，A的逆就為平方根法步驟1、對A做分解LLT

2、求R=L-13、求A-1=RRT平方根法步驟平方根法就是正定對稱矩陣的Cholesky分解平方根法就是正定對稱矩陣的Cholesky分解代碼doublesum;for(i=0;i<n;i++){for(j=i;j<n;j++){for(sum=a[i][j],k=i-1;k>=0;k--)sum-=a[i][k]*a[j][k];if(i==j){if(sum<=0.0)throw("Choleskyfailed");a[i][i]=sqrt(sum);}elsea[j][i]=sum/a[i][i];}}for(i=0;i<n;i++)for(j=0;j<i;j++)a[j][i]=0.;}代碼doublesum;5.改進的平方根法平方根法求平方是很耗時的，對于正定對稱陣來說有5.改進的平方根法平方根法求平方是很耗時的，對于正定對稱陣來n階對稱正定矩陣A的分解公式：n階對稱正定矩陣A的分解公式：6.補充內(nèi)容數(shù)值計算的誤差C#編程的問題6.補充內(nèi)容數(shù)值計算的誤差數(shù)值計算的誤差問題1、機器的精度(截斷誤差)2、數(shù)值計算的舍入誤差數(shù)值計算的誤差問題1、機器的精度(截斷誤差)截斷誤差浮點數(shù)的表達

用一個尾數(shù)（Mantissa），一個基數(shù)（Base），一個指數(shù)（Exponent）以及一個表示正負的符號來表達實數(shù)。機器精度

浮點數(shù)在32位機子上有兩種精度，float占32位，double占64位。IEEE標準754（1985）

符號位

階碼

尾數(shù)

長度

float

1

8

23

32

double

1

11

52

64

截斷誤差浮點數(shù)的表達例一書中代碼（C++）P17，P18+0.000001例一書中代碼（C++）例二（C#）如下C#代碼：floata=0.65f;floatb=0.6f;floatc=a-b;此時c為多少？0.05？錯誤！此時c為0.0499999523！小數(shù)部分為1/2,1/4,1/8,…,這樣的沒有問題例二（C#）如下C#代碼：截斷誤差原因float

var=5.2f;二進制表示：整數(shù)部分：

101小數(shù)部分：0.2*2=0.4*2=0.8*2=1.6(0.6)*2=1.2(0.2)*2=0.4*2=0.8*2=1.6(0.6)*2=1.2......

0

0

1

1

0

0

1

1

......

截斷誤差原因float

var=5.2f;float是32位，后面尾數(shù)長度最大23位。5.2的二進制表示就為：101.00110011001100110011一共23位。使用科學計數(shù)法表示，結(jié)果為：1.0100110011001100110011*22最高位為1，可以省略，在最后一位補一位：.01001100110011001100110*22

010000001010011001100110011001101位|←8位→||←--------------23位------------→|符號|階碼||尾數(shù)|float是32位，后面尾數(shù)長度最大23位。5.2的二進制表二進制到浮點數(shù)1

01.001100110011001100110

5.

xxxxxxxxxxxxxxxxxxxxxxxxxxxxx0.

0011001100110011001100+0*2-1+0*2-2+1*2-3+1*2-4+0*2-5+0*2-6+1*2-7+......+0*2-215.1999998二進制到浮點數(shù)數(shù)值計算的舍入誤差高斯消去中的主元選擇線性方程組的形態(tài)及解的誤差數(shù)值計算的舍入誤差高斯消去中的主元選擇高斯消去中的主元選擇方程其解為:(0,-1,1)T

高斯消去中的主元選擇方程高斯消去法，第一步消去后，得到

下一步消去乘數(shù)為2.5*103

，

繼續(xù)采用高斯消去得到最后一行：

高斯消去法，第一步消去后，得到假設(shè)5位的精度，做截斷圓整后得到：

這樣就得到：回代得到：假設(shè)5位的精度，做截斷圓整后得到：高斯消去法中，應(yīng)該選擇最大的元素，做行交換后到主元素。避免不必要的圓整誤差帶來的解得精度降低問題。高斯消去法中，應(yīng)該選擇最大的元素，做行交換后到主元素。避免不線性方程組的形態(tài)及解的誤差設(shè)x*是計算得到的解（計算機的精度導(dǎo)致的誤差）A非奇異，理論解一定精度的計算解與理論解的差異線性方程組的形態(tài)及解的誤差設(shè)x*是計算得到的解（計算機的精度例假設(shè)在十進制精度為3位的機器上計算：換元后高斯消去，乘數(shù)為消去得

例假設(shè)在十進制精度為3位的機器上計算：計算解與理論解的差異殘差值均是小于10-3，殘差很小。將主元素換為較大的值可以保證的僅僅是殘差很小。雖然采用了換元的高斯消去，但是由于截斷圓整對于這樣的方程仍然得到的不是正確的解計算解與理論解的差異殘差值均是小于10-3，殘差很小。將主元將精度擴展到6位的十進制，有

解得，

什么原因？將精度擴展到6位的十進制，有范數(shù)對于方程Ax=bM/m比值稱為A的條件數(shù)

范數(shù)由定義可得由上式可以看出，條件數(shù)可以大致度量出b的相對誤差對于解x的放大的倍數(shù)由定義可得由上式可以看出，條件數(shù)可以大致度量出b的相對誤差改變b值，得到改變b值，得到C#編程中的一些關(guān)鍵問題1、委托機制2、界面設(shè)計C#編程中的一些關(guān)鍵問題1、委托機制delegatedoubleProcessDelegate（doubleparam1，doubleparam2）；staticdoubleMultiply(doubleparam1,doubleparam2){returnparam1*param2;}staticdoubleDivide(doubleparam1,doubleparam2){returnparam1/param2;}delegate