迭代法(Iteration)也称辗转法,是一种不断用变量的旧值递推新值的过程,跟迭代法相对应的是直接法(或者称为一次解法),即一次性解决问题。迭代算法是用计算机解决问题的一种基本方法。它利用计算机运算速度快、适合做重复性操作的特点,让计算机对一组指令(或一定步骤)重复执行,在每次执行这组指令(或这些步骤)时,都从变量的原值推出它的一个新值。
// 自然数连加和int adds1(int n){int sum = 0;for(int i=1;i<=n; i++)sum+=i; // 迭代法return sum;}int adds2(int n){return (1+n)*n/2; // 直接法}
利用迭代算法解决问题,需要做好以下三个方面的工作:
① 确定迭代变量
在可以用迭代算法解决的问题中,至少存在一个可直接或间接地不断由旧值递推出新值的变量,这个变量就是迭代变量。
② 建立迭代关系式
所谓迭代关系式,指如何从变量的前一个值推出其下一个值的公式(或关系)。迭代关系式的建立是解决迭代问题的关键。
递推法也是一种利用旧值推新出新值的方法,但与迭代法不同的时,递推法往往会有一个初始条件,就如同分段函数一样,初始条件是递推表达式之一,而迭代法一般使用一个单一表达式。
③ 对迭代过程进行控制
在什么时候结束迭代过程?这是编写迭代程序必须考虑的问题。不能让迭代过程无休止地执行下去。迭代过程的控制通常可分为两种情况:一种是所需的迭代次数是个确定的值,可以计算出来;另一种是所需的迭代次数无法确定。对于前一种情况,可以构建一个固定次数的循环来实现对迭代过程的控制;对于后一种情况,需要进一步分析得出可用来结束迭代过程的条件。
迭代法又分为精确迭代和近似迭代。前者如辗转相除法,后者如“二分法”、"牛顿迭代法”来近似求函数值的问题。
1 二分法求函数值
对于函数y=f(x),设x的取值区间为[a,b]满足f(a)*f(b)<0k,也就是是f(x)在x的取值区间有0值,通过不断地把函数f(x)的零点所在的区间一分为二,使区间的两个端点逐步逼近零点,进而得到零点近似值的方法叫二分法。
用二分法求函数f(x)零点近似值的步骤如下:
1.1 确定区间[a,b],验证f(a)·f(b)<0,给定精确度ξ。
1.2 求区间(a,b)的中点c。
1.3 计算f(c)。
(1) 若f(c)=0,则c就是函数的零点;
(2) 若f(a)·f(c)<0,则令b=c;
(3) 若f(c)·f(b)<0,则令a=c。
(4) 判断是否达到精确度ξ:即若|a-b|
while (fabs(b-a) > e){ double c = (a+b)/2.0; if (f(a)* f(c) < 0) b = c; else a = c;}
测试代码:
#include #include #include #include double f(double x) // 测试函数f(x)=1+x-x*x*x{ return 1+x-x*x*x;}double bisect(double a,double b){double e = 1e-5; if (fabs(f(a)) <= e) printf("solution: %lg", a); else if (fabs(f(b)) <= e) printf("solution: %lg", b); else if (f(a)*f(b) > 0) printf("f(%lg)*f(%lg) > 0 ! need <= 0 !", a, b); else { while (fabs(b-a) > e) { double c = (a+b)/2.0; if (f(a)* f(c) < 0) b = c; else a = c; } return (a+b)/2.0; }}int main(){ printf("solution: %lg", bisect(1,2)); // solution: 1.32472 getchar(); return 0;}
二分法求平方根:
double sqrt1(double x) // 二分法 { double EPSINON = 1e-5; double low = 0.0; double high = x; double mid = (low + high) / 2; while ((high - low) > EPSINON){ if (mid*mid > x) high = mid; else low = mid; mid = (high + low) / 2; } return mid;}
2 牛顿迭代法求平方根
令f(x)=x^2-n,如图所示:
取x0,如果x0不是解,做一个经过(x0,f(x0))这个点的切线,与x轴的交点为x1。
同理,如果x1不是解,做一个经过(x1,f(x1))这个点的切线,与x轴的交点为x2。
依次类推。
以这样的方式得到的 xi 会无限趋近于 f(x)=0 的解。
(f(x)-f(xi))/(x-xi)=f’(x),f’(x)是斜率也是f(x)的导函数,即f’(x)=2x。
化简得:f(xi)=f(x)-f’(x)(x-xi),令f(xi)=0得:
(x²-n)-2x(x-xi)=0
持续化简得:
x² - n - 2x² + 2xxi=0
2xxi=x² + n
2xi=x+n/x
xi=(x+n/x)/2 // 迭代表达式
double sqrt2(double x) // 牛顿迭代法{ if (x == 0) { return 0; } double last = 0.0; double res = 1.0; while (res != last) // 判断前后两个解 xi 和 xi-1 是否无限接近 { last = res; res = (res + x / res) / 2; // 迭代表达式 } return res;}
3 欧几里德辗转相除法
最经典的迭代算法是欧几里德辗转相除法,用于计算两个整数a,b的最大公约数。枚举法虽然也可以求出最大公约数,但效率明显较差。
辗转相除法的计算原理:
gcd(a,b) = gcd(b,a % b)
a可以表示成a = kb + r,则r = a % b。假设d是a,b的一个公约数,则有 a%d==0,b%d==0,而r = a - kb,因此r%d==0(因为a%b==0,kb%d==0) ,因此d是(b,a % b)的公约数。
同理,假设d 是(b,a % b)的公约数,则 b%d==0,r%d==0 ,但是a = kb +r ,因此d也是(a,b)的公约数。
因此(a,b)和(b,a % b)的公约数是一样的,其最大公约数也必然相等。
欧几里德算法就是根据这个原理来做的,欧几里德算法又叫辗转相除法,它是一个反复迭代执行,直到余数等于0停止:
#include using namespace std;int gcdIter(int a,int b){ if (a<=0 || b<=0) return 0; int temp; while (b > 0) { temp = a % b; // 迭代关系式 a = b; b = temp; } return a;}int gcd(int a,int b){return b?gcd(b,a%b):a;} // 一次gcd()调用有a=b,b=a%bint main(){cout<
从上面的程序我们可以看到a,b是迭代变量,迭代关系是temp = a % b;根据迭代关系我们可以由旧值推出新值,然后循环执a = b; b = temp;直到迭代过程结束(余数为0)。
需要注意的是,一些问题可以用迭代法解决,也可以用递归法解决,如上面的辗转相除法。其迭代表达式隐含在函数递归调用的参数表达式中,通过隐式栈存储,其返回值也是如此。
附迭代法求平方根代码:
#include #include using namespace std;double sqrt1(double x) // 二分法 { double EPSINON = 1e-5; double low = 0.0; double high = x; double mid = (low + high) / 2; while ((high - low) > EPSINON){ if (mid*mid > x) high = mid; else low = mid; mid = (high + low) / 2; } return mid;}double sqrt2(double x) // 牛顿迭代法{ if (x == 0) { return 0; } double last = 0.0; double res = 1.0; while (res != last) // 判断前后两个解 xi 和 xi-1 是否无限接近 { last = res; res = (res + x / res) / 2; // 迭代表达式 } return res;}double sqrt3( double x ) { float y; float delta; float maxError; if ( x <= 0 ) return 0; y = x / 2; // initial guess maxError = x * 0.00001; // refine do { delta = ( y * y ) - x; y -= delta / ( 2 * y ); } while ( delta > maxError || delta < -maxError ); return y;}double sqrt4(double x) // 引自John Carmack 的Quake引擎中的源代码{ if(x == 0) return 0; float result = x; float xhalf = 0.5f*result; // float后加f转换成double类型 int i = *(int*)&result; i = 0x5f3759df - (i>>1); // 通过预先猜测的一个神奇数字 0x5f3759df, 来将迭代次数降到极限的一次 // (另外通过整形数的移位来进一步加速) result = *(float*)&i; result = result*(1.5f-xhalf*result*result); result = result*(1.5f-xhalf*result*result); return 1.0f/result; }int main(){ cout<
-End-