原理
据说这个平方根计算法出自《九章算术》，成书在公元一世纪前后，总结了先秦到两汉的古代数学成就。古代数学家也真是鬼才，怎么想出来的？
操作过程：小数点前后2位2位的隔开，不足2位用0补。然后首位试开方，第二位开始反复用公式 (Qn+1 * 20 + Qn) * Qn来试商 ; 下一结果Qn+2 = 10*Qn+1 + Qn。具体流程请自行搜索了解，大致流程如下：
代码
本篇采用的计算方法既非二分法也非牛顿迭代法，而是把中国古代的手工计算平方根的算术方法转换成代码来完成。代码写得有点繁杂，算是抛砖引玉吧，期待高手们写出更好的代码来。
sqrt(x, w=20, Float=False)
x 为非负实数，允许科学计数法；或者是可以转成非负实数的字符串。
w 是当sqrt(x)为无理数或结果的小数部分很长时的最大位数。(默认w=20)
Float 当Float=True输出浮点数，w设置如超出float精度则失效。(默认输出字符串类型)
def sqrt(x,w=20,Float=False): if isinstance(x,str): try: x = float(x) except: x = -1 if x<0: assert x>=0 or float(x)>=0 if x==0 or x==1:return str(x*1.0) s,power = str(x),0 t = s.find('e') if t!=-1: s=s.split('e') s[1] = int(s[1]) if int(s[1])%2: tmp = str(s[0]) s = [tmp[0]+tmp[2]+tmp[1]+tmp[3:],s[1]-1] s,power = *s, t = s.find('.') if t!=-1: if t%2: s='0'+s if not len(s)%2: s+='0' else: if len(s)%2: s='0'+s s = s.split('.') t = t,len(s[0])//2 try: s = s[0]+s[1] except: s = s[0] s = [s[i:i+2] for i in range(len(s)) if not i%2][::-1] n = int(s.pop()) i = 1 while n>0 and i*i<=n: i+=1 i -= 1 m,n = str(n-i*i),i div = lambda x,y:[(i,y-(20*x+i)*i) for i in range(10) if y>=(20*x+i)*i] dot,ret = t[1],str(n) while w>0: dot -= 1; if dot==0: ret += '.' if dot<=0: w-=1 m += s.pop() if s else '00' q,d = div(n,int(m))[-1] m,n = str(d),10*n+q ret += str(q) if len(s)==0 and d==0: if ret.find('.')==-1: ret += '.0' break if power!=0: if power>0: ret += 'e+0'+str(power//2) else: ret += 'e-0'+str(power//2)[1:] if Float: ret = float(ret) return ret
测试
sqrt(n)与n**0.5
>>> for i in range(17): print(f'sqrt({i})={sqrt(i)}') print(f'{i}**0.5={i**0.5}') sqrt(0)=0.0 0**0.5=0.0 sqrt(1)=1.0 1**0.5=1.0 sqrt(2)=1.41421356237309504880 2**0.5=1.4142135623730951 sqrt(3)=1.73205080756887729352 3**0.5=1.7320508075688772 sqrt(4)=2.0 4**0.5=2.0 sqrt(5)=2.23606797749978969640 5**0.5=2.23606797749979 sqrt(6)=2.44948974278317809819 6**0.5=2.449489742783178 sqrt(7)=2.64575131106459059050 7**0.5=2.6457513110645907 sqrt(8)=2.82842712474619009760 8**0.5=2.8284271247461903 sqrt(9)=3.0 9**0.5=3.0 sqrt(10)=3.16227766016837933199 10**0.5=3.1622776601683795 sqrt(11)=3.31662479035539984911 11**0.5=3.3166247903554 sqrt(12)=3.46410161513775458705 12**0.5=3.4641016151377544 sqrt(13)=3.60555127546398929311 13**0.5=3.605551275463989 sqrt(14)=3.74165738677394138558 14**0.5=3.7416573867739413 sqrt(15)=3.87298334620741688517 15**0.5=3.872983346207417 sqrt(16)=4.0 16**0.5=4.0 >>> >>> >>> sqrt(2,Float=True) == 2**0.5 True >>> sqrt(3,Float=True) == 3**0.5 True >>> sqrt(123,Float=True) == 123**0.5 True >>>
科学计算法测试
>>> sqrt(1.23e+50) '1.10905365064094171620e+025' >>> 1.23e+50**0.5 1.1090536506409416e+25 >>> sqrt(1.23e+51) '3.50713558335003638336e+025' >>> 1.23e+51**0.5 3.5071355833500364e+25 >>> sqrt(1.23e-50) '1.10905365064094171620e-025' >>> 1.23e-50**0.5 1.1090536506409418e-25 >>> sqrt(1.23e-51) '3.50713558335003638336e-026' >>> 1.23e-51**0.5 3.5071355833500365e-26 >>>
精度测试
>>> sqrt(2,100) '1.41421356237309504880168872420969807856967187537694 80731766797379907324784621070388503875343276415727' >>> sqrt(10,500) '3.16227766016837933199889354443271853371955513932521 68268575048527925944386392382213442481083793002951 87347284152840055148548856030453880014690519596700 15390334492165717925994065915015347411333948412408 53169295770904715764610443692578790620378086099418 28371711548406328552999118596824564203326961604691 31433612894979189026652954361267617878135006138818 62785804636831349524780311437693346719738195131856 78403231241795402218308045872844614600253577579702 82864402902440797789603454398916334922265261206779' >>> sqrt(123,2000) '11.09053650640941716205160010260993291846337674245402 00228773128390850016331012896052334560795952104923 97609680678955280792187905933115292625456231839306 77251943912251383176574941199469582196976000438135 40867472202696805822192936674286399297485980945076 78295660716567970352602444301464684924479636459099 14867193001802834979182445692668356613065880869548 25999929108256938950212808340223106891096223696155 58407975630369894453553840958193501976809032496435 98351605065147790119568667920441106521130541775416 88403618227856731042118992185281429619080341919784 91236339758444278806715795552342614568993355758259 00257454016962686033789212494918951906742724387717 86816775058970553110606825772728456503631816127185 97236514403953501043370950197906664648656587402375 47456007486620199478701365589929806411967684259454 52983826062182085142259698311212942126801021825240 49746216571370198706996668818361845968235199845816 69902568378075870406180767774203205467306930309438 16832697339952984981648138573158982379175644441470 50866454616067044997958059525209028336061119869487 42330123683852090419874148218422184948658174234109 57266920859313915279433828848348895422861007776818 60008218306382698626583716405660864687800252861028 31434598682004467042622580768146827595778429365752 19552960734589394860581552514632178318401717998004 12904477515916582731378709757426533231624676879051 21277303475419893966421915682962327965202183939540 07041575389626966661084388718863067485822881867488 00698329255619212758441301451712107894934052042668 13620261632356745839351541875140890682676375675084 74266141686980135645708807796240849684849035273656 15026469136389861472142022294507047846522584450886 83823021978858519029304962564972861464819539784518 43062402845984165814550404820570869649363216162428 32571582506529019502195720558586845830492641836612 65156560044594234298289454649568716622838760718998 14917847213523793754337786561375901264942957378010 57481468928695787589958271831183026476218342365982 72020291818311722410861264031411990033164594115086' >>>
会被python强制转科学计数法的数建议用字符串形式：
>>> 12345679801234567890.123 1.2345679801234567e+19 >>> sqrt('12345679801234567890.123',30) '3513641956.892387732139011965897613363716' >>> sqrt(12345679801234567890.123,30) '3.513641956892387605472292516306e+09' >>> sqrt(1.2345679801234567e+19,30) '3.513641956892387605472292516306e+09'
1000位到1万位的耗时测试
>>> from time import time >>> for i in range(1,11): t = time() s = sqrt(2,i*1000) print(i*1000,':',time()-t) 1000 : 0.09099745750427246 2000 : 0.5609970092773438 3000 : 1.7489948272705078 4000 : 3.7850000858306885 5000 : 7.239997863769531 6000 : 12.108997106552124 7000 : 19.044997692108154 8000 : 28.065003395080566 9000 : 39.57198643684387 10000 : 54.04999876022339 >>>
其他方法
牛顿迭代法
牛顿迭代法的思想就是在一条曲线上从某一点的切线开始，首先求其与横轴的交点，之后再确定曲线上和该交点横坐标相同的点，并重复求该点的切线与横坐标的交点的方式，不断逼近真实解的过程。网上相关的讲解很多，我这里就简单总结一下牛顿迭代法的步骤：
1.确定需要求解的函数y=f(x)，在求平方根中该函数为f(x)=x*x-c;
2.假设给定初始点的横坐标为x0，则其对应的切线方程为Q(x)=f '(x0)*(x-x0) + f(x0),在求平方根的算法中该切线方程为Q(x) = 2*x0*(x - x0) + x0*x0-c;
3.根据切线方程与横坐标的交点得到下一个迭代点的横坐标，若前一迭代点的横坐标记为Xn，则下一迭代点的横坐标记为Xn+1，令第二条中的x0=Xn,Q(x) = 0可以得到Xn+1的表达式：
Xn+1 = (c/Xn + Xn) / 2
def newton_sqrt(n): x = n y = (x + n / x) / 2 while(x - y > 0.00001): x = y y = (x + n / x) / 2 return x print(newton_sqrt(9)) print(newton_sqrt(10)) print(10**0.5)
输出：
3.000000001396984
3.162277665175675
3.1622776601683795
二分法
二分法的理论基础是零点存在定理。先确定零点所在的区间范围[a,b]（区间范围端点值异号），通过不断缩小这个区间范围，达到逐渐逼近最终解的目的。而缩小区间范围的方法就是取[a,b]的中点m，判断f(m)的正负号与f(a)、f(b)的关系，改变a的值（f(b)*f(m)<0）或者b的值（f(a）*f(m)<0)为m,每次缩减一半的解空间。循环上述操作，直到最终的范围小于所给的精度范围。
def binary_sqrt(n): low = 0 high = n while low <= high: mid = (low + high) / 2 if mid * mid == n: return mid elif mid * mid > n: high = mid - 0.000001 else: low = mid + 0.000001 return mid print(binary_sqrt(9)) print(binary_sqrt(10)) print(10**0.5)
输出：
3.0000000242946148
3.162277569391012
3.1622776601683795
