2678145490787694710138406683675886762424440548647327525594133873038347266950777818028529804877891353347244255046176220880989519574806598935727647635506939520211209794328712328126947740602255036264392628020291100698432900022333877728824270241063030019536831205264726059774673235680796770108334187808344470141093691232926021244583649747591446066921535333854157441330339217364515328375889937313870826987520723914810358707312468857074455945610038118763873992681688857200178452843740327143990247797893489008808544535536081089994800692174759199578596435124153437689121057209523619989593900355037400039597165241729693472001589594642879166693443994939359639600558001937547209387175085119128575833396037478944824322315197771750730768383267287661855880237549141571520374996029107650341071074081501632143441365461865431362023790463469688071768208320226561675630348708531708110373889241713511809029065048200924851629507265095151843687852114948736130517024578482059771284515915597919026074110811882700307150650858080695566672584075215926450360584738705106755945308199980660478183563586945902705730231409948938478910176292179410952310592864745025373902249844607078729230098512959498522851990040833671698873183076996740726479481352429756295966452788332789259316653805673821716832204368483794501474641036177256005810416507649770210480747056764697918800240201155395628561684853684747063434500664449114760179367086926249813006056735202795763742271472136500817313178050955623517148960973128433455798444839070215800699075256371341531892841721456473301085595261371522078416420399702683504171019191363372883472415693921509069484448792009627325245658854979159002549542769741504596589556474585723101403486046952706847309188488390899096649210853675228464611999630360164703837131620319658277117300208326944026347808645954842782309633205645471959842390075986572806483755223497552549948718545978420793425036666212191832327909389117021539426917904222971515424022606157676379417252939026404363841301595348953340657745965591468370831371449860774296757556463049237248472613472436666799615694672489600751763956185893123579218797539971253828850163437053164723845196089466298000486042246621634263064752906475198879616484926366053597071510886820878425907736453843197228282179943365206039514850477377542684361730162273246709577442085159313947668552462343596819299335060863266863059421286080919061657175757802930378698062290163092494909661982639303946271889619023417245433500782685534406291391527022657608608130680134939380126048353657648559735671562230714451258764762350451448379438849319710473631625533527602085050858633346936315115309333355986556546113134787213610042970405680162084678021281774169131725296977361816657904233850711837826724976507876061043779177845550646612327663283839673950405472039432309339110648301347672597370034346896334064135099613033528785924976781865447506814375627647103020761819104412269362847848884053177378762432543285984193091289854949394023586451139393716370633759317698908834326333339634557238526607966375935960004339535161203546369910853171524516945554862019134684985868547068709986010719120259421723270478586992153989283736650948010894951434857207445141886625290777550599700738909777829412358941157674866550887608199469723859806486384262574182811908513849905040209776132986875560077482527819093584150609347560543627976011009349905362597935257895287412087233810106693141683912683308201407373637490835862452911394718399632093441745957478079754723679119343225931665785380948807114301164623717144211285944154289093317610519159281033583963015810556464180819928345256832572962092994993950655913067768151113345566231088759030148164710287289335782462899525700500822021472132532771242198844498579102733028256523810973405350508253183399668290339844052874956694050028526577775764032598023622192221226501334606994447870277619118208162133562541678250491647671042612088902952091163094215645902682086024337098910110924115381788866213335441814530141911754446270296630426675020592640816713810936850961517475097707179713115856513919631998747192606018612535436009162735908155778165736923393260822217415131632296810207091918499370220533980666393617256227781299273913099896641458499949922475383775878771011361417556799406151773618626138749998464514570482965286184387350245762317271666329142667987131139265167425606863654107927682419958167046956962340419295394058625063862272314570914511355913547470055072878835349142049813368019960713974852351833144227544851625677753419141348585495260230041111323576432417502998832485027071744397105338790364666312325483218209150368759120749671419554589928540715813531144241705501211397506662448276056294585639532800634727428900861538486753011006379264618338782036690340859441153926551954112364946638768223012718907190897165613382555036872673785444132692292217184974818914501886689736449912131900762592363447945036085692999163815947845753538999001633767986447907729462453017601用时0.9394531秒,有4888位(这是利用快速傅里叶变换的乘法结果,改进了一下,5位一组,速度稍有提高,但仍然不如前面的模仿手工计算的快速乘法)
Private Sub Command1_Click()
Dim m, n
m = Trim(Text1): n = Trim(Text2)
ts = Timer
c = MbC4(Trim(m), Trim(n))
Text3 = c & "用时" & Timer - ts & "秒,有" & Len(c) & "位"
End Sub
Private Sub Command2_Click()
Text1 = ""
Text2 = ""
Text3 = ""
End Sub
Public Function MbC4(D1 As String, D2 As String) As String '快速乘法
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
Dim xr() As Double, a As String
a = Trim(D1)
B = Trim(D2)
X = Len(a) \ 5: Y = Len(B) \ 5
a = String(Val(X * 5 + 5 - Len(a)), "0") & a
B = String(Val(Y * 5 + 5 - Len(B)), "0") & B
X = X + 1: Y = Y + 1
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
a = String(Val(sb) * 5 - Len(a), "0") & a
B = String(Val(sb) * 5 - Len(B), "0") & B
ReDim x_(1 To sb): ReDim y_(1 To sb)
For i1 = 1 To sb
x_(i1) = Mid(a, (sb - i1 + 1) * 5 - 4, 5): y_(i1) = Mid(B, (sb - i1 + 1) * 5 - 4, 5)
If Len(x_(i1)) < 5 Then
x_(i1) = String(5 - Len(x_(i1)), "0") & x_(i1)
ElseIf Len(y_(i1)) < 5 Then
y_(i1) = String(5 - Len(y_(i1)), "0") & y_(i1)
Else
x_(i1) = x_(i1): y_(i1) = y_(i1)
End If
Next
ReDim xr(0 To (Len(a) - 5) \ 5): ReDim yr(0 To (Len(B) - 5) \ 5): ReDim zr(0 To (Len(B) - 5) \ 5)
If Len(a) = 5 Then
xr(0) = a: yr(0) = B
Else
Dim I As Long, J As Long, mn As Long, lh As Long, 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
xr(I + 1) = x_(J + 1): yr(I + 1) = y_(J + 1)
Next
xr(0) = x_(1): xr(1) = x_(1 + sb / 2)
yr(0) = y_(1): yr(1) = y_(1 + sb / 2)
End If
Dim xi(): Dim yi(): Dim zi()
n = sb '求数组大小,其值必须是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")
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
J = sb
ReDim x_(1 To sb): ReDim y_(1 To sb)
For k = 1 To J
n1 = n1 + 1
ReDim Preserve x_(1 To n1)
x_(n1) = zr(n1 - 1): y_(n1) = zi(n1 - 1)
x_(n1) = Format(Val(x_(n1)), "0.000000"): y_(n1) = Format(Val(y_(n1)), "0.000000")
Next
'位序倒置
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
xr(I + 1) = x_(J + 1): yr(I + 1) = y_(J + 1)
'js = js & "/" & x_(J + 1)
'js1 = js1 & "/" & y_(J + 1)
Next
'sx1 = "/" & x_(1) & "/" & x_(1 + sb / 2) & js
'sy1 = "/" & y_(1) & "/" & y_(1 + sb / 2) & js1
xr(0) = x_(1): xr(1) = x_(1 + sb / 2)
yr(0) = y_(1): yr(1) = y_(1 + sb / 2)
ns = Len(a) \ 5: Jn = ns
ReDim zr(0 To ns - 1)
m = 0
l = 2
pi = 3.14159265358979
Do
l = l + l
m = m + 1
Loop Until l > ns
ns = l / 2
ReDim xi(ns - 1): ReDim yi(ns - 1): ReDim zi(ns - 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 > ns - 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 ns - 1 '仅输出模
zr(I) = (xr(I) - yi(I)) / n
zr(I) = Format(Val(zr(I) + 0.5), "0.000000")
If InStr(zr(I), ".") = 0 Then
s121 = zr(I)
Else
s121 = Left(zr(I), InStr(zr(I), ".") - 1)
End If
's0 = "/" & s121 & s0
zr(I) = s121
Next
For i1 = 1 To Val(Jn - sb1 + 1)
zr(sb1 + i1 - 2) = 0
Next
For i1 = 0 To n - 1
If zr(i1) < 0 Then
zr(i1) = "00000"
ElseIf Len(zr(i1)) < 5 Then
zr(i1) = String(5 - Len(zr(i1)), "0") & zr(i1)
Else
zr(i1) = zr(i1)
End If
's5 = s5 & "/" & zr(i1)
If i1 = 0 Then
s6 = Val(Left(zr(i1), Len(zr(i1)) - 5))
If Len(s6) < 5 Then
s6 = String(5 - Len(s6), "0") & s6
Else
s6 = s6
End If
s8 = Right(zr(i1), 5)
ElseIf Val(zr(i1)) >= 0 Then
s7 = Val(zr(i1)) + Val(s6)
If Len(s7) = 5 Or Len(s7) = 10 Or Len(s7) = 15 Then
s7 = s7
Else
s7 = String(20 - Len(s7), "0") & s7
End If
s10 = Right(s7, 5)
s11 = s10 & s11
If Len(s7) < 5 Then
s7 = String(5 - Len(s7), "0") & s7
ElseIf Len(s7) = 5 Then
s6 = "00000"
Else
s7 = s7
s6 = Val(Left(s7, Len(s7) - 5))
End If
Else
s6 = s6
End If
Next
s9 = s6 & s11 & s8
s9 = qdqd0(Trim(s9))
's2 = nifft(dxcx1(Trim(s)), dxcx1(Trim(s1)), Trim(sb1))
's3 = nifft(Trim(sx1), Trim(sy1), Trim(sb1))
MbC4 = s9
End Function
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
[此贴子已经被作者于2021-5-4 18:41编辑过]