代码如下,加了大整数加法程序,因为最后变换得到的数组会超过15位的,有时候可能是有3*8=24位,所以8位一组不好,最好4位,最多弄5位的可能是就可以了,有可能可以省去大数加法程序。明天试试吧,今天的程序如下:
Dim l As Long, le As Long, le1 As Long, n As Long, r As Long, p As Long, q As Long, m As Byte
Dim wr As Double, w1 As Double, wlr As Double, wl1 As Double, tr As Double, t1 As Double
Dim pi As Double, t As Double, tr1 As Double
Private Sub Command1_Click()
Dim xr() As String, a As String
a = Trim(Text1)
b = Trim(Text3)
ts = Timer
x = Len(a) \ 8: Y = Len(b) \ 8
If Val(8 * x) = Len(a) Then
a = a
ElseIf Val(8 * Y) = Len(b) Then
b = b
Else
a = String(Val(x * 8 + 8 - Len(a)), "0") & a
b = String(Val(Y * 8 + 8 - Len(b)), "0") & b
x = x + 1: Y = Y + 1
End If
sb1 = x + Y
sb2 = Log(sb1) / Log(2)
If InStr(sb2, ".") = 0 Then
sb2 = sb2
Else
sb2 = Int(sb2) + 1
End If
sb = 2 ^ sb2
Print sb
a = String(Val(sb) * 8 - Len(a), "0") & a
b = String(Val(sb) * 8 - Len(b), "0") & b
a = dxcx0(Trim(a), Val(sb)): b = dxcx0(Trim(b), Val(sb))
Print a
ReDim xr(0 To (Len(a) - 8) \ 8): ReDim yr(0 To (Len(b) - 8) \ 8): ReDim zr(0 To (Len(b) - 8) \ 8)
If Len(a) = 8 Then
xr(0) = a: yr(0) = b
Else
For i1 = 0 To (Len(a) - 8) \ 8
xr(i1) = Mid(a, (i1 + 1) * 8 - 7, 8)
yr(i1) = Mid(b, (i1 + 1) * 8 - 7, 8)
Next
End If
Dim xi(): Dim yi(): Dim zi()
n = Len(a) \ 8 '求数组大小,其值必须是2的幂
m = 0
l = 2
pi = 3.14159265358979
Do
l = l + l
m = m + 1
Loop Until l > n
n = l / 2
ReDim xi(n - 1): ReDim yi(n - 1): ReDim zi(n - 1)
l = 1
Do
le = 2 ^ l
le1 = le / 2
wr = 1
wi = 0
If l = 1 Then
t = 0
Else
t = pi / le1
End If
w1r = Cos(t)
w1i = -Sin(t)
r = 0
Do
p = r
Do
q = p + le1
tr = xr(q) * wr - xi(q) * wi
ti = xr(q) * wi + xi(q) * wr
tr1 = yr(q) * wr - yi(q) * wi
ti1 = yr(q) * wi + yi(q) * wr
xr(q) = xr(p) - tr
xi(q) = xi(p) - ti
xr(p) = xr(p) + tr
xi(p) = xi(p) + ti
yr(q) = yr(p) - tr1
yi(q) = yi(p) - ti1
yr(p) = yr(p) + tr1
yi(p) = yi(p) + ti1
xr(p) = Format(Val(xr(p)), "0.000000"): xi(p) = Format(Val(xi(p)), "0.000000")
yr(p) = Format(Val(yr(p)), "0.000000"): yi(p) = Format(Val(yi(p)), "0.000000")
Print xr(p), xr(q)
p = p + le
Loop Until p > n - 1
wr2 = wr * w1r - wi * w1i
wi2 = wr * w1i + wi * w1r
wr = wr2
wi = wi2
r = r + 1
Loop Until r > le1 - 1
l = l + 1
Loop Until l > m
For I = 0 To n - 1 '仅输出模
zr(I) = xr(I) * yr(I) - xi(I) * yi(I): zi(I) = xr(I) * yi(I) + xi(I) * yr(I)
zr(I) = Format(Val(zr(I)), "0.000000"): zi(I) = Format(Val(zi(I)), "0.000000")
s = s & "/" & zr(I)
s1 = s1 & "/" & zi(I)
Next
s2 = nifft(dxcx1(Trim(s)), dxcx1(Trim(s1)), Trim(sb1))
Text2 = s2 & "有" & Len(s2) & "位,用时" & Timer - ts & "秒"
End Sub
Private Sub Command2_Click()
Text1 = ""
Text2 = ""
Text3 = ""
Form1.Cls
End Sub
Private Function qdqd0(sa As String) As String
a = sa
Do While Left(a, 1) = "0"
a = Mid(a, 2)
Loop
If a = "" Then
a = 0
Else
a = a
End If
qdqd0 = a
End Function
Private Function nifft(sa As String, sb As String, sb1 As String) As String
Dim xi(): Dim yi(): Dim zi()
Dim xr(), yr()
Dim zr() As String
s2 = Split(sa, "/")
s3 = Split(sb, "/")
J = UBound(s2)
n = J
For k = 1 To J
n1 = n1 + 1
ReDim Preserve xr(0 To n1 - 1)
ReDim Preserve yr(0 To n1 - 1)
xr(n1 - 1) = s2(n1): yr(n1 - 1) = s3(n1)
xr(n1 - 1) = Format(Val(xr(n1 - 1)), "0.000000"): yr(n1 - 1) = Format(Val(yr(n1 - 1)), "0.000000")
Next
ReDim zr(0 To J - 1)
m = 0
l = 2
pi = 3.14159265358979
Do
l = l + l
m = m + 1
Loop Until l > n
n = l / 2
ReDim xi(n - 1): ReDim yi(n - 1): ReDim zi(n - 1)
l = 1
Do
le = 2 ^ l
le1 = le / 2
wr = 1
wi = 0
If l = 1 Then
t = 0
Else
t = -1 * pi / le1
End If
w1r = Cos(t)
w1i = -Sin(t)
r = 0
Do
p = r
Do
q = p + le1
tr = xr(q) * wr - xi(q) * wi
ti = xr(q) * wi + xi(q) * wr
tr1 = yr(q) * wr - yi(q) * wi
ti1 = yr(q) * wi + yi(q) * wr
xr(q) = xr(p) - tr
xi(q) = xi(p) - ti
xr(p) = xr(p) + tr
xi(p) = xi(p) + ti
yr(q) = yr(p) - tr1
yi(q) = yi(p) - ti1
yr(p) = yr(p) + tr1
yi(p) = yi(p) + ti1
xr(p) = Format(Val(xr(p)), "0.000000"): xi(p) = Format(Val(xi(p)), "0.000000")
yr(p) = Format(Val(yr(p)), "0.000000"): yi(p) = Format(Val(yi(p)), "0.000000")
p = p + le
Loop Until p > n - 1
wr2 = wr * w1r - wi * w1i
wi2 = wr * w1i + wi * w1r
wr = wr2
wi = wi2
r = r + 1
Loop Until r > le1 - 1
l = l + 1
Loop Until l > m
For I = 0 To n - 1 '仅输出模
zr(I) = (xr(I) - yi(I)) / n
zr(I) = Format(Val(zr(I) + 0.5), "0.000000")
If InStr(zr(I), ".") = 0 Then
s1 = zr(I)
Else
s1 = Left(zr(I), InStr(zr(I), ".") - 1)
End If
s = "/" & s1 & s
zr(I) = s1
Next
For i1 = 1 To Val(J - sb1 + 1)
zr(sb1 + i1 - 2) = 0
Next
For i1 = 0 To n - 1
If zr(i1) < 0 Then
zr(i1) = "00000000"
Else
zr(i1) = zr(i1)
End If
s5 = s5 & "/" & zr(i1)
If i1 = 0 Then
If Len(zr(i1)) < 8 Then
zr(i1) = String(8 - Len(zr(i1)), "0") & zr(i1)
Else
zr(i1) = zr(i1)
End If
s6 = Val(Left(zr(i1), Len(zr(i1)) - 8))
If Len(s6) < 8 Then
s6 = String(8 - Len(s6), "0") & s6
Else
s6 = s6
End If
s8 = Right(zr(i1), 8)
ElseIf Val(zr(i1)) >= 0 Then
s7 = MPC1(Trim(zr(i1)), Trim(s6))
s10 = Right(s7, 8)
s11 = s10 & s11
If Len(s7) < 8 Then
s7 = String(8 - Len(s7), "0") & s7
ElseIf Len(s7) = 8 Then
s6 = "00000000"
Else
s7 = s7
s6 = Val(Left(s7, Len(s7) - 8))
End If
Else
s6 = s6
End If
Next
s9 = s6 & s11 & s8
nifft = qdqd0(Trim(s9))
End Function
Private Function dxcx0(sa As String, sb As String) As String
Dim x_() As String, a As String
a = Trim(sa)
ReDim x_(1 To sb)
For i1 = 1 To sb
x_(i1) = Mid(a, (sb - i1 + 1) * 8 - 7, 8)
If Len(x_(i1)) < 8 Then
x_(i1) = String(8 - Len(x_(i1)), "0") & x_(i1)
Else
x_(i1) = x_(i1)
End If
Next
Dim n As Integer, I As Long, J As Long, mn As Long, lh As Long, t As Double, k As Long
'位序倒置
n = sb '求数组大小,其值必须是2的幂
lh = n / 2
J = n / 2
For I = 1 To n - 2
Debug.Print I, J
k = lh '下面是向右进位算法
Do
If k > J Then Exit Do '高位是1吗
J = J - k '是的,高位置0
k = k / 2 '准备次高位的权
Loop Until k = 0 '次高位的权若非0,则检查新的次高位
J = J + k '非则若最高位是0,则置1
s = s & x_(J + 1)
Next
dxcx0 = x_(1) & x_(1 + sb / 2) & s
End Function
Private Function dxcx1(sa As String) As String
Dim x_() As Double, a As String
a = Trim(sa)
s2 = Split(sa, "/")
s3 = Split(sb, "/")
J = UBound(s2)
sb = J
ReDim x_(1 To sb)
For k = 1 To J
n1 = n1 + 1
ReDim Preserve x_(1 To n1)
x_(n1) = Format(Val(s2(n1)), "0.000000")
Next
Dim n As Integer, I As Long, mn As Long, lh As Long, t As Double
'位序倒置
n = sb '求数组大小,其值必须是2的幂
lh = n / 2
J = n / 2
For I = 1 To n - 2
Debug.Print I, J
k = lh '下面是向右进位算法
Do
If k > J Then Exit Do '高位是1吗
J = J - k '是的,高位置0
k = k / 2 '准备次高位的权
Loop Until k = 0 '次高位的权若非0,则检查新的次高位
J = J + k '非则若最高位是0,则置1
x_(J + 1) = Format(Val(x_(J + 1)), "0.000000")
s = s & "/" & x_(J + 1)
Next
x_(1) = Format(Val(x_(1)), "0.000000"): x_(1 + sb / 2) = Format(Val(x_(1 + sb / 2)), "0.000000")
dxcx1 = "/" & x_(1) & "/" & x_(1 + sb / 2) & s
End Function
Public Function MPC1(D1 As String, D2 As String) As String 'jiafa
Dim x, Y '两数长度
If Len(D1) >= Len(D2) Then
D4 = String(Len(D1) - Len(D2), "0") & D2
d3 = D1
Else
D4 = D2
d3 = String(Len(D2) - Len(D1), "0") & D1
End If
x = Len(d3): Y = Len(D4)
Dim a() As Integer, B1() As Integer, C1() As Integer, E1() As Integer
ReDim a(1 To x)
ReDim B1(1 To Y)
ReDim C1(1 To x)
ReDim E1(1 To x)
Dim I, J, C2, CJ, JW
For J = Y To 1 Step -1 'D2
JW = 0 '进位清0
B1(J) = Mid$(D4, J, 1) '每位数
For I = x To 1 Step -1 'D1
a(I) = Mid$(d3, I, 1) '每位数
C1(I) = a(I) + B1(I) + JW '计算jia
JW = C1(I) \ 10
E1(I) = C1(I) Mod 10
Next
Next
For r = 1 To x
If JW = 0 Then
MPC1 = MPC1 & E1(r)
Else
jc = jc & E1(r)
MPC1 = JW & jc
End If
Next
End Function
Dim l As Long, le As Long, le1 As Long, n As Long, r As Long, p As Long, q As Long, m As Byte
Dim wr As Double, w1 As Double, wlr As Double, wl1 As Double, tr As Double, t1 As Double
Dim pi As Double, t As Double, tr1 As Double
Private Sub Command1_Click()
Dim xr() As String, a As String
a = Trim(Text1)
b = Trim(Text3)
ts = Timer
x = Len(a) \ 8: Y = Len(b) \ 8
If Val(8 * x) = Len(a) Then
a = a
ElseIf Val(8 * Y) = Len(b) Then
b = b
Else
a = String(Val(x * 8 + 8 - Len(a)), "0") & a
b = String(Val(Y * 8 + 8 - Len(b)), "0") & b
x = x + 1: Y = Y + 1
End If
sb1 = x + Y
sb2 = Log(sb1) / Log(2)
If InStr(sb2, ".") = 0 Then
sb2 = sb2
Else
sb2 = Int(sb2) + 1
End If
sb = 2 ^ sb2
Print sb
a = String(Val(sb) * 8 - Len(a), "0") & a
b = String(Val(sb) * 8 - Len(b), "0") & b
a = dxcx0(Trim(a), Val(sb)): b = dxcx0(Trim(b), Val(sb))
Print a
ReDim xr(0 To (Len(a) - 8) \ 8): ReDim yr(0 To (Len(b) - 8) \ 8): ReDim zr(0 To (Len(b) - 8) \ 8)
If Len(a) = 8 Then
xr(0) = a: yr(0) = b
Else
For i1 = 0 To (Len(a) - 8) \ 8
xr(i1) = Mid(a, (i1 + 1) * 8 - 7, 8)
yr(i1) = Mid(b, (i1 + 1) * 8 - 7, 8)
Next
End If
Dim xi(): Dim yi(): Dim zi()
n = Len(a) \ 8 '求数组大小,其值必须是2的幂
m = 0
l = 2
pi = 3.14159265358979
Do
l = l + l
m = m + 1
Loop Until l > n
n = l / 2
ReDim xi(n - 1): ReDim yi(n - 1): ReDim zi(n - 1)
l = 1
Do
le = 2 ^ l
le1 = le / 2
wr = 1
wi = 0
If l = 1 Then
t = 0
Else
t = pi / le1
End If
w1r = Cos(t)
w1i = -Sin(t)
r = 0
Do
p = r
Do
q = p + le1
tr = xr(q) * wr - xi(q) * wi
ti = xr(q) * wi + xi(q) * wr
tr1 = yr(q) * wr - yi(q) * wi
ti1 = yr(q) * wi + yi(q) * wr
xr(q) = xr(p) - tr
xi(q) = xi(p) - ti
xr(p) = xr(p) + tr
xi(p) = xi(p) + ti
yr(q) = yr(p) - tr1
yi(q) = yi(p) - ti1
yr(p) = yr(p) + tr1
yi(p) = yi(p) + ti1
xr(p) = Format(Val(xr(p)), "0.000000"): xi(p) = Format(Val(xi(p)), "0.000000")
yr(p) = Format(Val(yr(p)), "0.000000"): yi(p) = Format(Val(yi(p)), "0.000000")
Print xr(p), xr(q)
p = p + le
Loop Until p > n - 1
wr2 = wr * w1r - wi * w1i
wi2 = wr * w1i + wi * w1r
wr = wr2
wi = wi2
r = r + 1
Loop Until r > le1 - 1
l = l + 1
Loop Until l > m
For I = 0 To n - 1 '仅输出模
zr(I) = xr(I) * yr(I) - xi(I) * yi(I): zi(I) = xr(I) * yi(I) + xi(I) * yr(I)
zr(I) = Format(Val(zr(I)), "0.000000"): zi(I) = Format(Val(zi(I)), "0.000000")
s = s & "/" & zr(I)
s1 = s1 & "/" & zi(I)
Next
s2 = nifft(dxcx1(Trim(s)), dxcx1(Trim(s1)), Trim(sb1))
Text2 = s2 & "有" & Len(s2) & "位,用时" & Timer - ts & "秒"
End Sub
Private Sub Command2_Click()
Text1 = ""
Text2 = ""
Text3 = ""
Form1.Cls
End Sub
Private Function qdqd0(sa As String) As String
a = sa
Do While Left(a, 1) = "0"
a = Mid(a, 2)
Loop
If a = "" Then
a = 0
Else
a = a
End If
qdqd0 = a
End Function
Private Function nifft(sa As String, sb As String, sb1 As String) As String
Dim xi(): Dim yi(): Dim zi()
Dim xr(), yr()
Dim zr() As String
s2 = Split(sa, "/")
s3 = Split(sb, "/")
J = UBound(s2)
n = J
For k = 1 To J
n1 = n1 + 1
ReDim Preserve xr(0 To n1 - 1)
ReDim Preserve yr(0 To n1 - 1)
xr(n1 - 1) = s2(n1): yr(n1 - 1) = s3(n1)
xr(n1 - 1) = Format(Val(xr(n1 - 1)), "0.000000"): yr(n1 - 1) = Format(Val(yr(n1 - 1)), "0.000000")
Next
ReDim zr(0 To J - 1)
m = 0
l = 2
pi = 3.14159265358979
Do
l = l + l
m = m + 1
Loop Until l > n
n = l / 2
ReDim xi(n - 1): ReDim yi(n - 1): ReDim zi(n - 1)
l = 1
Do
le = 2 ^ l
le1 = le / 2
wr = 1
wi = 0
If l = 1 Then
t = 0
Else
t = -1 * pi / le1
End If
w1r = Cos(t)
w1i = -Sin(t)
r = 0
Do
p = r
Do
q = p + le1
tr = xr(q) * wr - xi(q) * wi
ti = xr(q) * wi + xi(q) * wr
tr1 = yr(q) * wr - yi(q) * wi
ti1 = yr(q) * wi + yi(q) * wr
xr(q) = xr(p) - tr
xi(q) = xi(p) - ti
xr(p) = xr(p) + tr
xi(p) = xi(p) + ti
yr(q) = yr(p) - tr1
yi(q) = yi(p) - ti1
yr(p) = yr(p) + tr1
yi(p) = yi(p) + ti1
xr(p) = Format(Val(xr(p)), "0.000000"): xi(p) = Format(Val(xi(p)), "0.000000")
yr(p) = Format(Val(yr(p)), "0.000000"): yi(p) = Format(Val(yi(p)), "0.000000")
p = p + le
Loop Until p > n - 1
wr2 = wr * w1r - wi * w1i
wi2 = wr * w1i + wi * w1r
wr = wr2
wi = wi2
r = r + 1
Loop Until r > le1 - 1
l = l + 1
Loop Until l > m
For I = 0 To n - 1 '仅输出模
zr(I) = (xr(I) - yi(I)) / n
zr(I) = Format(Val(zr(I) + 0.5), "0.000000")
If InStr(zr(I), ".") = 0 Then
s1 = zr(I)
Else
s1 = Left(zr(I), InStr(zr(I), ".") - 1)
End If
s = "/" & s1 & s
zr(I) = s1
Next
For i1 = 1 To Val(J - sb1 + 1)
zr(sb1 + i1 - 2) = 0
Next
For i1 = 0 To n - 1
If zr(i1) < 0 Then
zr(i1) = "00000000"
Else
zr(i1) = zr(i1)
End If
s5 = s5 & "/" & zr(i1)
If i1 = 0 Then
If Len(zr(i1)) < 8 Then
zr(i1) = String(8 - Len(zr(i1)), "0") & zr(i1)
Else
zr(i1) = zr(i1)
End If
s6 = Val(Left(zr(i1), Len(zr(i1)) - 8))
If Len(s6) < 8 Then
s6 = String(8 - Len(s6), "0") & s6
Else
s6 = s6
End If
s8 = Right(zr(i1), 8)
ElseIf Val(zr(i1)) >= 0 Then
s7 = MPC1(Trim(zr(i1)), Trim(s6))
s10 = Right(s7, 8)
s11 = s10 & s11
If Len(s7) < 8 Then
s7 = String(8 - Len(s7), "0") & s7
ElseIf Len(s7) = 8 Then
s6 = "00000000"
Else
s7 = s7
s6 = Val(Left(s7, Len(s7) - 8))
End If
Else
s6 = s6
End If
Next
s9 = s6 & s11 & s8
nifft = qdqd0(Trim(s9))
End Function
Private Function dxcx0(sa As String, sb As String) As String
Dim x_() As String, a As String
a = Trim(sa)
ReDim x_(1 To sb)
For i1 = 1 To sb
x_(i1) = Mid(a, (sb - i1 + 1) * 8 - 7, 8)
If Len(x_(i1)) < 8 Then
x_(i1) = String(8 - Len(x_(i1)), "0") & x_(i1)
Else
x_(i1) = x_(i1)
End If
Next
Dim n As Integer, I As Long, J As Long, mn As Long, lh As Long, t As Double, k As Long
'位序倒置
n = sb '求数组大小,其值必须是2的幂
lh = n / 2
J = n / 2
For I = 1 To n - 2
Debug.Print I, J
k = lh '下面是向右进位算法
Do
If k > J Then Exit Do '高位是1吗
J = J - k '是的,高位置0
k = k / 2 '准备次高位的权
Loop Until k = 0 '次高位的权若非0,则检查新的次高位
J = J + k '非则若最高位是0,则置1
s = s & x_(J + 1)
Next
dxcx0 = x_(1) & x_(1 + sb / 2) & s
End Function
Private Function dxcx1(sa As String) As String
Dim x_() As Double, a As String
a = Trim(sa)
s2 = Split(sa, "/")
s3 = Split(sb, "/")
J = UBound(s2)
sb = J
ReDim x_(1 To sb)
For k = 1 To J
n1 = n1 + 1
ReDim Preserve x_(1 To n1)
x_(n1) = Format(Val(s2(n1)), "0.000000")
Next
Dim n As Integer, I As Long, mn As Long, lh As Long, t As Double
'位序倒置
n = sb '求数组大小,其值必须是2的幂
lh = n / 2
J = n / 2
For I = 1 To n - 2
Debug.Print I, J
k = lh '下面是向右进位算法
Do
If k > J Then Exit Do '高位是1吗
J = J - k '是的,高位置0
k = k / 2 '准备次高位的权
Loop Until k = 0 '次高位的权若非0,则检查新的次高位
J = J + k '非则若最高位是0,则置1
x_(J + 1) = Format(Val(x_(J + 1)), "0.000000")
s = s & "/" & x_(J + 1)
Next
x_(1) = Format(Val(x_(1)), "0.000000"): x_(1 + sb / 2) = Format(Val(x_(1 + sb / 2)), "0.000000")
dxcx1 = "/" & x_(1) & "/" & x_(1 + sb / 2) & s
End Function
Public Function MPC1(D1 As String, D2 As String) As String 'jiafa
Dim x, Y '两数长度
If Len(D1) >= Len(D2) Then
D4 = String(Len(D1) - Len(D2), "0") & D2
d3 = D1
Else
D4 = D2
d3 = String(Len(D2) - Len(D1), "0") & D1
End If
x = Len(d3): Y = Len(D4)
Dim a() As Integer, B1() As Integer, C1() As Integer, E1() As Integer
ReDim a(1 To x)
ReDim B1(1 To Y)
ReDim C1(1 To x)
ReDim E1(1 To x)
Dim I, J, C2, CJ, JW
For J = Y To 1 Step -1 'D2
JW = 0 '进位清0
B1(J) = Mid$(D4, J, 1) '每位数
For I = x To 1 Step -1 'D1
a(I) = Mid$(d3, I, 1) '每位数
C1(I) = a(I) + B1(I) + JW '计算jia
JW = C1(I) \ 10
E1(I) = C1(I) Mod 10
Next
Next
For r = 1 To x
If JW = 0 Then
MPC1 = MPC1 & E1(r)
Else
jc = jc & E1(r)
MPC1 = JW & jc
End If
Next
End Function