实验5 常微分方程初值问题和矩阵特征值的计算
思考问题：欧拉折线法与泰勒级数法的关系是什么？如何把乘幂法变为反幂法来求解按模最小特征值？
答：欧拉折线法就是用一系列折线近似的代替曲线，欧拉折线法思想就是泰勒级数展开，其解得精度是差分步长的一阶精度。

泰勒公式展开是一种高阶显式一步法，理论上它具有任意阶精度；
把乘幂法变为反幂法来求解按模最小特征值的方法为该求方阵A 的逆矩阵
1-A ，乘幂法对应的模最大特征值λ取倒数为λ
1，变为模最小特征值，由此得到反幂法，可反复迭代：)1()(-=k k Ay y ，)1()(/-k k y y 收敛于A 的主特征值，)(k y 为对应的特征向量。

5.1 取步长2.0=h ，用显示欧拉法求⎩⎨⎧=--='1
)0(2
y xy y y 在6.0=x 处y 的近似值。

提示：循环3次，用3段折线来逼近函数y 。

运行结果为：
5.2 已知矩阵⎥⎥⎥⎦
⎤⎢⎢⎢⎣⎡=361641593642A ，用改进后的幂乘法求A 得主特征值及其对应的特征向量。

运行结果为：
程序代码
5.1
#include<stdio.h>
#include <conio.h>
#define MAXSIZE 50
double f(double x,double y);
void main(void)
{
double a,b,h,x[MAXSIZE],y[MAXSIZE];
long i,n;
printf("\nPlease input interval a,b:");
scanf("%lf,%lf",&a,&b);
printf("\nPlease input step length h:");
scanf("%lf",&h);
n=(long)((b-a)/h);
x[0]=a;
printf("\nPlease input the start point of x[0]=%lf ordinate y[0]:",x[0]); scanf("%lf",&y[0]);
for(i=0;i<n;i++)
{
x[i+1]=x[i]+h;
y[i+1]=y[i]+h * f(x[i],y[i]);
}
printf("\nThe calcuate result is:");
for(i=0;i<=n;i++)
printf("\nx[%ld]=%lf,y[%ld]=%lf",i,x[i],i,y[i]);
getch();
}
double f(double x,double y)
{
return(-y-x*y*y); /* Calcuate and return to function value f(x,y) */ }
5.2
#include<stdio.h>
#include <conio.h>
#include<math.h>
#define MAXSIZE 50
void input (double a[][MAXSIZE],double v[],long n);
void matrix_product (double a[][MAXSIZE],double u[],double v[],long n);
void normalization (double u[],double v[],long n,double * pm1);
void output (double v[],long n,double m1);
void main (void)
{
double a[MAXSIZE][MAXSIZE],u[MAXSIZE],v[MAXSIZE];
double epsilon,m0,m1;
long n,maxk,k;
printf("\nPlease input the square A's order number: ");
scanf("%ld",&n);
input(a,v,n);
printf("\nPlease input the max number of iteratings:");
scanf("%ld",&maxk);
printf("\nPlease input the main eigenvalue'sprecision:");
scanf("%lf",&epsilon);
matrix_product(a,u,v,n);
normalization(u,v,n,&m1);
for (k=1;k<=maxk;k++)
{
m0=m1;
matrix_product(a,u,v,n);
normalization(u,v,n,&m1);
if (fabs(m0-m1)<=epsilon)
break;
}
if (k<=maxk)
output (v,n,m1);
else
printf("\nThe number of iterations has the ultra limit."); }
/* A subroutine 1: Read in square A and the initial vector */
void input (double a[][MAXSIZE],double v[],long n)
{
long i,j;
printf("\nPlease input %ld order square A:\n",n);
for (i=0;i<=n-1;i++)
for (j=0;j<=n-1;j++)
scanf("%lf",&a[i][j]);
printf("\nPlease input the initial iterations vector:");
for (i=0;i<=n-1;i++)
scanf("%lf",&v[i]);
}
/* A subroutine 2: Calcuate U=AV */
void matrix_product (double a[][MAXSIZE],double u[],double v[],long n) {
long i,j;
for (i=0;i<=n-1;i++)
{
u[i]=0;
for (j=0;j<=n-1;j++)
u[i]+=a[i][j]*v[j];
}
}
/* A subroutine 3: To normalized the vector U’s length */
void normalization (double u[],double v[],long n,double * pm1)
{
long i;
*pm1=u[0];
for (i=1;i<=n-1;i++)
if(fabs(*pm1)<fabs(u[i]))
*pm1=u[i];
for (i=0;i<=n-1;i++)
v[i]=u[i]/(*pm1);
}
/* A subroutine 4: Output the result of calcuate */
void output (double v[],long n,double m1)
{
long i;
printf("\nThe square A's eigenvalue is about: %f",m1); printf("\nIts corresponding eigenvectors is about: \n"); for (i=0;i<=n-1;i++)
printf(" %lf",v[i]);
getch();
}。