常見的四種單片機開方根算法分享

C語言中要求平方根，可以在頭文件中加入#include 《math.h》。然后調(diào)用sqrt（n）;函數(shù)即可。但在單片機中調(diào)用此函數(shù)無疑會耗費大量資源和時間，是極不合適的。在此，總結(jié)下網(wǎng)上常見的四種單片機常用開方根算法：

對于擁有專門的乘除法指令的單片機，可采用以下兩種方法：

1、二分法

對于一個非負(fù)數(shù)n，它的平方根不會小于大于（n/2+1）（謝謝@linzhi-cs提醒）。在［0， n/2+1］這個范圍內(nèi)可以進行二分搜索，求出n的平方根。

-------------------------------------------------------------------------------

1 int sqrt（int x） { 2 long long i = 0; 3 long long j = x / 2 + 1; 4 while （i 《= j） 5 { 6 long long mid = （i + j） / 2; 7 long long sq = mid * mid; 8 if （sq == x） return mid; 9 else if （sq 《 x） i = mid + 1;10 else j = mid - 1;11 }12 return j;13 }

-------------------------------------------------------------------------------

2、更為常用的牛頓迭代法

-------------------------------------------------------------------------------

1 int sqrt（int x） { 2 if （x == 0） return 0; 3 double last = 0; 4 double res = 1; 5 while （res ！= last） 6 { 7 last = res; 8 res = （res + x / res） / 2; 9 }10 return int（res）;11 }

-------------------------------------------------------------------------------

牛頓迭代法也可以求解多次方程。

對于不帶乘除法指令的單片機，可采取以下兩種算法：

算法3：

本算法只采用移位、加減法、判斷和循環(huán)實現(xiàn)，因為它不需要浮點運算，也不需要乘除運算，因此可以很方便地運用到各種芯片上去。

我們先來看看10進制下是如何手工計算開方的：

先看下面兩個算式：

x = 10*p + q （1）

公式（1）左右平方之后得：

x^2 = 100*p^2 + 20pq + q^2 （2）

現(xiàn)在假設(shè)我們知道x^2和p，希望求出q來，求出了q也就求出了x^2的開方x了。

我們把公式（2）改寫為如下格式：

q = （x^2 - 100*p^2）/（20*p+q） （3）

這個算式左右都有q，因此無法直接計算出q來，因此手工的開方算法和手工除法算法一樣有一步需要猜值。

我們來一個手工計算的例子：計算1234567890的開方

首先我們把這個數(shù)兩位兩位一組分開，計算出最高位為3。也就是（3）中的p，最下面一行的334為余數(shù)，也就是公式（3）中的（x^2 - 100*p^2）近似值

3 --------------- | 12 34 56 78 90 9 --------------- | 3 34

下面我們要找到一個0-9的數(shù)q使它最接近滿足公式（3）。我們先把p乘以20寫在334左邊：

3 q --------------- | 12 34 56 78 90 9 --------------- 6q| 3 34

我們看到q為5時（60+q*q）的值最接近334，而且不超過334。于是我們得到：

3 5 --------------- | 12 34 56 78 90 9 --------------- 65| 3 34 | 3 25 --------------- 9 56

接下來就是重復(fù)上面的步驟了，這里就不再啰嗦了。

這個手工算法其實和10進制關(guān)系不大，因此我們可以很容易的把它改為二進制，改為二進制之后，公式（3）就變成了：

q = （x^2 - 4*p^2）/（4*p+q） （4）

我們來看一個例子，計算100（二進制1100100）的開方：

1 0 1 0 --------------- | 1 10 01 00 1 --------------- 100| 0 10 | 0 00 --------------- | 10 011001| 10 01 --------------- 0 00

這里每一步不再是把p乘以20了，而是把p乘以4，也就是把p右移兩位，而由于q的值只能為0或者1，所以我們只需要判斷余數(shù)（x^2 - 4*p^2）和（4*p+1）的大小關(guān)系，如果余數(shù)大于等于（4*p+q）那么該上一個1，否則該上一個0。

下面給出完成的C語言程序，其中root表示p，rem表示每步計算之后的余數(shù)，divisor表示（4*p+1），通過a》》30取a的最高 2位，通過a《《=2將計算后的最高2位剔除。其中root的兩次《《1相當(dāng)于4*p。程序完全是按照手工計算改寫的，應(yīng)該不難理解。

-------------------------------------------------------------------------------

unsigned short sqrt（unsigned long a）{ unsigned long rem = 0; unsigned long root = 0; unsigned long divisor = 0; for（int i=0; i《16; i++）{ root 《《= 1; rem = （（rem 《《 2） + （a 》》 30））; a 《《= 2; divisor = （root《《1） + 1; if（divisor 《= rem）{ rem -= divisor; root++; } } return （unsigned short）（root）; }

-------------------------------------------------------------------------------

算法4

這種方法比牛頓迭代法更加快速的方法。

1.原理

下述用pow（X，Y）表示X的Y次冪，用B［0］，B［1］，。。。，B［m-1］表示一個序列，

其中［x］為下標(biāo)。

假設(shè)：

B［x］，b［x］都是二進制序列，取值0或1。

1、 M = B［m-1］*pow（2，m-1） + B［m-2］*pow（2，m-2） + 。。。 + B［1］*pow（2，1） + B［0］*pow

（2，0）

2、 N = b［n-1］*pow（2，n-1） + b［n-2］*pow（2，n-2） + 。。。 + b［1］*pow（2，1） + n［0］*pow

（2，0）

3、 pow（N，2） = M

（1） N的最高位b［n-1］可以根據(jù)M的最高位B［m-1］直接求得。

設(shè) m 已知，因為 pow（2， m-1） 《= M 《= pow（2， m），所以 pow（2， （m-1）/2） 《= N 《=

pow（2， m/2）

如果 m 是奇數(shù)，設(shè)m=2*k+1，

那么 pow（2，k） 《= N 《 pow（2， 1/2+k） 《 pow（2， k+1），

n-1=k， n=k+1=（m+1）/2

如果 m 是偶數(shù)，設(shè)m=2k，

那么 pow（2，k） 》 N 》= pow（2， k-1/2） 》 pow（2， k-1），

n-1=k-1，n=k=m/2

所以b［n-1］完全由B［m-1］決定。

余數(shù) M［1］ = M - b［n-1］*pow（2， 2*n-2）

（2） N的次高位b［n-2］可以采用試探法來確定。

因為b［n-1］=1，假設(shè)b［n-2］=1，則 pow（b［n-1］*pow（2，n-1） + b［n-1］*pow（2，n-2），

2） = b［n-1］*pow（2，2*n-2） + （b［n-1］*pow（2，2*n-2） + b［n-2］*pow（2，2*n-4）），

然后比較余數(shù)M［1］是否大于等于 （pow（2，2）*b［n-1］ + b［n-2］） * pow（2，2*n-4）。這種

比較只須根據(jù)B［m-1］、B［m-2］、。。。、B［2*n-4］便可做出判斷，其余低位不做比較。

若 M［1］ 》= （pow（2，2）*b［n-1］ + b［n-2］） * pow（2，2*n-4）， 則假設(shè)有效，b［n-2］ =

1；

余數(shù) M［2］ = M［1］ - pow（pow（2，n-1）*b［n-1］ + pow（2，n-2）*b［n-2］， 2） = M［1］ -

（pow（2，2）+1）*pow（2，2*n-4）；

若 M［1］ 《 （pow（2，2）*b［n-1］ + b［n-2］） * pow（2，2*n-4）， 則假設(shè)無效，b［n-2］ =

0；余數(shù) M［2］ = M［1］。

（3） 同理，可以從高位到低位逐位求出M的平方根N的各位。

使用這種算法計算32位數(shù)的平方根時最多只須比較16次，而且每次比較時不必把M的各位逐

一比較，尤其是開始時比較的位數(shù)很少，所以消耗的時間遠低于牛頓迭代法。

2. 實現(xiàn)代碼

這里給出實現(xiàn)32位無符號整數(shù)開方得到16位無符號整數(shù)的C語言代碼。

-------------------------------------------------------------------------------

/****************************************/ /*Function： 開根號處理 */ /*入口參數(shù)：被開方數(shù)，長整型 */ /*出口參數(shù)：開方結(jié)果，整型 */ /****************************************/ unsigned int sqrt_16（unsigned long M） { unsigned int N， i; unsigned long tmp， ttp; // 結(jié)果、循環(huán)計數(shù) if （M == 0） // 被開方數(shù)，開方結(jié)果也為0 return 0; N = 0; tmp = （M 》》 30）; // 獲取最高位：B［m-1］ M 《《= 2; if （tmp 》 1） // 最高位為1 { N ++; // 結(jié)果當(dāng)前位為1，否則為默認(rèn)的0 tmp -= N; } for （i=15; i》0; i--） // 求剩余的15位 { N 《《= 1; // 左移一位 tmp 《《= 2; tmp += （M 》》 30）; // 假設(shè) ttp = N; ttp = （ttp《《1）+1; M 《《= 2; if （tmp 》= ttp） // 假設(shè)成立 { tmp -= ttp; N ++; } } return N; }

-------------------------------------------------------------------------------

以上算法結(jié)尾網(wǎng)上收集所得，雖然原理可能比較難懂，但都可在單片機中實際運用。

責(zé)任編輯人：CC