当前位置:文档之家› 最小二乘法拟合曲线的VB程序源代码

最小二乘法拟合曲线的VB程序源代码

'按教材定义LSS方法求系数的函数(P123),输入G,y()矩阵,不求误差(另外定义函数)
Private Function answerX(G() As Single, y() As Single) As Single()
Dim i As Integer, j As Integer, k As Integer
Dim xx() As Single, w() As Single
Dim sigma As Single, t As Single, bat As Single
Dim m As Integer, n As Integer
Dim all As Single
m = UBound(G, 1)
n = UBound(G, 2)
Dim GG() As Single
ReDim GG(m, n + 1)
ReDim xx(n)
ReDim w(m)

'将y放入G最后一列,定义为GG矩阵
For i = 1 To m
For j = 1 To n
GG(i, j) = G(i, j)
Next j
GG(i, n + 1) = y(i)
Next i

'开始LSS方法
For k = 1 To n

'形成矩阵Qk
all = 0
For i = k To m
all = all + GG(i, k) ^ 2
Next i
all = Sqr(all)
sigma = -Sgn(GG(k, k)) * all
w(k) = GG(k, k) - sigma
For j = k + 1 To m
w(j) = GG(j, k)
Next j
bat = sigma * w(k)

'变换GGk-1到GGk
GG(k, k) = sigma
For j = k + 1 To n + 1
all = 0
For i = k To m
all = all + w(i) * GG(i, j)
Next i
t = all / bat
For i = k To m
GG(i, j) = GG(i, j) + t * w(i)
Next i
Next j

'解三角方程Rx=h1
xx(n) = GG(n, n + 1) / GG(n, n)
For i = n - 1 To 1 Step -1
all = 0
For j = i + 1 To n
all = all + GG(i, j) * xx(j)
Next j
xx(i) = (GG(i, n + 1) - all) / GG(i, i)
Next i
answerX = xx

' Dim shuchu As String
' Open fullpath & "xishua.txt" For Output As 1
' shuchu = ""
' For i = 1 To n
' shuchu = shuchu & xx(i) & Space(5)
' Next i
' Print #1, shuchu
' Close 1

Next k
End Function


'定义n次多项式最小二乘法拟合时求系数矩阵函数
'系数形势为p(x)=a(1)+a(2)x+a(3)x^2+...+a(n)x^(n-1)+a(n+1)x^n
'因此,g1(x)=x^0,g2(x)=x,g3(x)=x^2...gn(x)=x^(n-1),g(n+1)=x^n
'需给出x(1-m)和y(1-m)矩阵
Private Function answerA(x() As Single, y() As Single, n As Integer) As Single()
Dim m As Integer
Dim i As Integer, j As Integer
Dim aa() As Single, G() As Single
m = UBound(x)
ReDim aa(n + 1)
ReDim G(m, n + 1)

'求出G矩阵
For j = 1 To n + 1
For i = 1 To m
G(i, j) = x(i) ^ (j - 1)
Next i
Next j

' Dim shuchu As String
' Open fullpath & "MATG.txt" For Output As 1
' For i = 1 To m
' shuchu = ""
' For j = 1 To n + 1
' shuchu = shuchu & G(i, j) & Space(5)
' Next j
' Print #1, shuchu
' Next i
' Close 1
'求出法方程的系数矩阵和右端矩阵,

并求解系数矩阵
Dim a() As Single, d() As Single
ReDim a(n + 1, n + 1)
ReDim d(n + 1)
'a() = mattime(matT(G()), G(), False)
'd() = mattime(matT(G()), y(), True)
'aa() = answerX(a(), d())
aa() = answerX(G(), y())
answerA = aa
End Function

'定义n次多项式求值函数,a()为系数矩阵,xx为带入时刻(整数)
Private Function Pxn(a() As Single, xx As Single) As Single
Dim n As Single
Dim i As Integer
n = UBound(a)
Pxn = 0
For i = 1 To n
Pxn = Pxn + a(i) * xx ^ (i - 1)
Next i
End Function


'定义题目要求指数函数求值函数,按2次多项式求值,ae()为系数矩阵,x为带入时刻(整数)
Private Function Pxe(a() As Single, xx As Single) As Single
Pxe = a(1) * Exp(-a(2) * (xx - a(3)) ^ 2)
End Function

'定义n次多项式拟合求误差函数,n为次数
Private Function EN(n As Integer) As Single
Dim m As Integer
Dim aa() As Single
aa = answerA(x(), y(), n)
m = UBound(x())
EN = 0
Dim z As Integer
For z = 1 To m
EN = EN + (Pxn(aa, x(z)) - y(z)) ^ 2
Next z
EN = Sqr(EN)
End Function

'定义多项式函数形式生成函数,n为次数
Private Function Pnstr(n As Integer) As String
Dim m As Integer
Dim aa() As Single
aa = answerA(x(), y(), n)
m = UBound(x())
Pnstr = ""
Dim z As Integer
For z = 1 To n + 1
Pnstr = Pnstr & IIf(aa(z) > 0, "+" & aa(z), aa(z)) & IIf((z - 1) = 0, "", "t^" & CStr(z - 1))
Next z
End Function

'定义题目要求指数函数拟合误差函数
Private Function EE() As Single
Dim m As Integer
Dim aa() As Single
aa = answerA(x(), ye(), 2)
ae(1) = Exp(aa(1) - aa(2) ^ 2 / (4 * aa(3)))
ae(2) = -aa(3)
ae(3) = -aa(2) / (2 * aa(3))
m = UBound(x)
EE = 0
Dim z As Integer
For z = 1 To m
EE = EE + (Pxe(ae(), x(z)) - ye(z)) ^ 2
Next z
EE = Sqr(EE)
End Function

Private Sub cmd1_Click()
txte2 = EN(2)
txte3 = EN(3)
txte4 = EN(4)
txtee = EE
txtpe = ae(1) & "Exp(" & -ae(2) & "(t-" & ae(3) & ")^2)"
txtp2 = Pnstr(2)
txtp3 = Pnstr(3)
txtp4 = Pnstr(4)
Command2.Enabled = True
End Sub

Private Sub Command1_Click()
End
End Sub

Private Sub Command2_Click()
Dim i As Single
Dim a2() As Single, a3() As Single, a4() As Single
a2() = answerA(x, y, 2)
a3() = answerA(x, y, 3)
a4() = answerA(x, y, 4)
pic1.DrawWidth = 2
For i = -2 To 26 Step 0.01
pic1.Line (i - 0.01, Pxn(a2, (i - 0.01)))-(i, Pxn(a2, i)), vbMagenta
Next i
For i = -2 To 26 Step 0.01
pic1.Line (i - 0.01, Pxn(a3, (i - 0.01)))-(i, Pxn(a3, i)), vbBlue
Next i
For i = -2 To 26 Step 0.01
pic1.Line (i - 0.01, Pxn(a4, (i - 0.01)))-(i, Pxn(a4, i)), vbCyan
Next i
For i = -2 To 26 Step 0.01
pic1.Line (i - 0.01, Pxe(ae(), (i - 0.01)))-(i,

Pxe(ae(), i)), vbGreen
Next i

End Sub

Private Sub Form_Load()
'窗口加载定义数组x,y,ye
Dim i As Integer
ReDim x(25)
ReDim y(25)
ReDim ye(25)
For i = 1 To 25
x(i) = i - 1
Next i
y(1) = 15
y(2) = 14
y(3) = 14
y(4) = 14
y(5) = 14
y(6) = 15
y(7) = 16
y(8) = 18
y(9) = 20
y(10) = 20
y(11) = 23
y(12) = 25
y(13) = 28
y(14) = 31
y(15) = 32
y(16) = 41
y(17) = 29
y(18) = 27
y(19) = 25
y(20) = 24
y(21) = 22
y(22) = 20
y(23) = 18
y(24) = 17
y(25) = 16
For i = 1 To 25
ye(i) = Log(y(i))
Next i
fram2.Top = fram1.Top
fram2.Left = fram1.Left
fram2.Visible = False
Command2.Enabled = False
' If Right(App.Path, 1) = "\" Then fullpath = App.Path Else fullpath = App.Path & "\"

End Sub

Private Sub Form_activate()
'绘制坐标线和各数据点
Dim i As Integer
pic1.Scale (-2, 45)-(26, 10)
pic1.DrawWidth = 1
For i = -2 To 26
pic1.Line (i, 10)-(i, 45), vbBlack
Next i
For i = 10 To 45 Step 5
pic1.Line (-2, i)-(26, i), vbBlack
Next i
pic1.DrawWidth = 8
For i = 1 To 25
pic1.PSet (x(i), y(i)), vbRed
Next i
End Sub
Private Sub TabStrip1_Click()
Select Case TabStrip1.SelectedItem.Index
Case 1
fram1.Visible = True
fram2.Visible = False
Case Else
fram1.Visible = False
fram2.Visible = True
End Select
End Sub

相关主题
相关文档 最新文档