wyl.cpp

来自「算法包括:1.二分法求解 2.牛顿法求解 3.高斯消去法求解 4.雅可比迭代法求」· C++ 代码 · 共 333 行

#include <iostream>
#include <string>
#include <math.h>
using namespace std;
void twoFF()
{cout<<"所求方程为x*x*x-x*x-2*x+1=0\n";
double a,b;
cout<<"请输入求根区间：\n";
cout<<"a=";
cin>>a;
cout<<"b=";
cin>>b;
double c,d;
cout<<"请输入要求达到的精度：\n";
cout<<"精度为";
cin>>d;
cout<<"a=0\n";
cout<<"b=1\n";
int n=0;
while(c>d || c<-d)
{if((a*a*a-a*a-2*a+1)*((a+b)*(a+b)*(a+b)/8-(a+b)*(a+b)/4-a-b+1)>0){a=(a+b)/2;cout<<"a="<<a<<endl;}
else {b=(a+b)/2;cout<<"b="<<b<<endl;}
n++;
c=a-b;}
cout<<"迭代次数为"<<n<<endl;
}

void Newton()
{cout<<"求方程x*(x+1)*(x+1)-1=0在x=0.4附近的根\n";
double x=0.4,m;
double c,d;
cout<<"请输入要求达到的精度：\n";
cout<<"精度为";
cin>>d;
cout<<"x=0.4\n";
int n=0;
while(c>d || c<-d){m=x;
x = x -(x*(x+1)*(x+1)-1) /(3*x*x+4*x+1);
cout<<"x="<<x<<endl;
c=m-x;
n++;}
cout<<"迭代次数为"<<n<<endl;
}

int n;
double a[100][101];
double m[100];
void changerow(int g)
{int i,j;
if(a[g-1][g-1]==0)
{int k;
for(k=g;k<n;k++)
{if(a[k][g-1]!=0)break;
}
for(j=0;j<n+1;j++)
{m[j]=a[g-1][j];
a[g-1][j]=a[k][j];
a[k][j]=m[j];}
}
for(i=g;i<n;i++)
m[i]=a[i][g-1]/a[g-1][g-1];
for(i=g;i<n;i++)
for(j=g-1;j<n+1;j++)
a[i][j]=a[i][j]-m[i]*a[g-1][j];

for(i=0;i<n;i++)
{cout<<"\n";
	for(j=0;j<n+1;j++)
    printf("    %2f",a[i][j]);
    }
cout<<"\n";
}
void GSXU()
{
printf("输入方程个数:");
cin>>n;
printf("按行输入各个方程的系数:\n");
int i,j;
for(i=0;i<n;i++)
for(j=0;j<n+1;j++)
cin>>a[i][j];
for(i=0;i<n;i++)
{cout<<"\n";
	for(j=0;j<n+1;j++)
    printf("    %2f",a[i][j]);}
cout<<"\n";
for(i=1;i<n;i++)
{changerow(i);}

double x[100];
int r=n-1;
x[n-1]=a[n-1][n]/a[n-1][n-1];
double p=0,q;
for(i=n-2;i>=0;i--)
{
for(j=n-1;j>=r;j--)
{q=a[i][j]*x[j];
p=p+q;
}
x[i]=(a[i][n]-p)/a[i][i];
r--;
p=0;}
printf("方程组解为");
cout<<"X="<<"(";
for(i=0;i<n-1;i++)
cout<<x[i]<<", ";
cout<<x[n-1];
cout<<")"<<endl;
}

void Jacobi()
{double a[4][4], x1[100],x2[100],x3[100];
x1[0]=0;x2[0]=0;x3[0]=0;
printf("按行输入三元方程的各个系数:\n");
int i,j;
for(i=0;i<3;i++)
for(j=0;j<4;j++)
cin>>a[i][j];
x1[1]=a[0][3]/a[0][0];
x2[1]=a[1][3]/a[1][1];
x3[1]=a[2][3]/a[2][2];
double c=1;
for(i=1;c>0.000001;i++)
{
x1[i]=x1[1]-a[0][1]/a[0][0]*x2[i-1]-a[0][2]/a[0][0]*x3[i-1];
x2[i]=x2[1]-a[1][0]/a[1][1]*x1[i-1]-a[1][2]/a[1][1]*x3[i-1];
x3[i]=x3[1]-a[2][0]/a[2][2]*x1[i-1]-a[2][1]/a[2][2]*x2[i-1];
c=fabs(x1[i]-x1[i-1])+fabs(x2[i]-x2[i-1])+fabs(x3[i]-x3[i-1]);
cout<<i<<":  "<<x1[i]<<"  ,  "<<x2[i]<<"  ,  "<<x3[i]<<endl;
}
}

double x[100],y[100];
void input(int n)
{
int i;
for(i=0;i<n;i++)
{cout<<"x["<<i<<"]=";
scanf("%lf",&x[i]);
cout<<"y["<<i<<"]=";
scanf("%lf",&y[i]);}
}
void Lagrange()
{
int m;
printf("输入点的个数：\n");
scanf("%d",&m);
double v;
printf("输入所求X的值：\n");
scanf("%lf",&v);
input(m);
int i,j;
double a,b=1,c,d=1,e,f=1,g,h=1,p,q=0;
for(i=0;i<m;i++)
{  for(j=0;j<i;j++)
   {a=x[i]-x[j];
    b=b*a;
	c=v-x[j];
	d=d*c;}
   for(j=i+1;j<m;j++)
   {e=x[i]-x[j];
    f=f*e;
	g=v-x[j];
	h=h*g;}
p=d*h*y[i]/b/f;
cout<<"L"<<i<<"="<<p<<endl;
q=q+p;
b=1;
d=1;
f=1;
h=1;
}
cout<<"f("<<v<<")="<<q<<endl;
}
void NewtonCZ()
{int m;
printf("输入点的个数：\n");
scanf("%d",&m);
double v;
printf("输入所求X的值：\n");
scanf("%lf",&v);
input(m);
int i,j;
double z[100][100];
for(i=0;i<m;i++)
z[i][0]=y[i];
for(i=1;i<m;i++)
{cout<<i<<"阶差商：\n";
for(j=i;j<m;j++)
{z[j][i]=(z[j][i-1]-z[j-1][i-1])/(x[j]-x[j-i]);
cout<<z[j][i]<<endl;
}  }
double a,b=1,c;
c=y[0];
for(i=0;i<m-1;i++)
{
a=v-x[i];
b=b*a;
c=c+z[i+1][i+1]*b;
}
cout<<"f("<<v<<")="<<c<<endl;
}
double LJ1(double x[],int y,int z)
{double s=0;
for(int i=0;i<y;i++)
s=s+pow(x[i],z);
return s;
}
double LJ2(double x[],int y,int z,double q[])
{double s=0;
for(int i=0;i<y;i++)
s=s+pow(x[i],z)*q[i];
return s;
}
void zxecf()
{int m,n;
cout<<"输入变量个数和拟和阶数"<<endl;
cin>>m>>n;
double a[100],b[100];
for(int i=0;i<m;i++)
{cout<<"输入的x"<<i<<"值和y"<<i<<"的值"<<endl;
cin>>a[i]>>b[i];}
double c[100][100],d[100];
for( i=0;i<n+1;i++)
for(int j=0;j<n+1;j++)
{c[i][j]=LJ1(a,m,i+j);}
c[0][0]=9;
cout<<"最小二乘法法方程组系数为"<<endl;
for( i=0;i<n+1;i++)
{d[i]=LJ2(a,m,i,b);
for(int j=0;j<n+1;j++)
printf("    %.2f",c[i][j]);
cout<<"    "<<d[i]<<endl;}
}

double f(double x)
{double v;
 v=exp(-x*x);
 return v;
}
void Romberge()
{double a,b,t[9],s[5],c[3],r;
cout<<"输入区间两端点"<<endl;
cin>>a>>b;
t[1]=(b-a)/2*(f(a)+f(b));
t[2]=0.5*t[1]+(b-a)/2*f(a+(b-a)/2);
s[1]=4.0/3*t[2]-1.0/3*t[1];
t[4]=0.5*t[2]+(b-a)/4*(f(a+(b-a)/4)+f(a+3*(b-a)/4));
s[2]=4.0/3*t[4]-1.0/3*t[2];
c[1]=16.0/15*s[2]-1.0/15*s[1];
t[8]=0.5*t[4]+(b-a)/8*(f(a+(b-a)/8)+f(a+3*(b-a)/8)+f(a+5*(b-a)/8)+f(a+7*(b-a)/8));
s[4]=4.0/3*t[8]-1.0/3*t[4];
c[2]=16.0/15*s[4]-1.0/15*s[2];
r=64.0/63*c[2]-c[1]/63;
cout<<"t[1]="<<t[1]<<endl;
cout<<"t[2]="<<t[2]<<endl;
cout<<"t[4]="<<t[4]<<endl;
cout<<"t[8]="<<t[8]<<endl;
cout<<"s[1]="<<s[1]<<endl;
cout<<"s[2]="<<s[2]<<endl;
cout<<"s[4]="<<s[4]<<endl;
cout<<"c[1]="<<c[1]<<endl;
cout<<"c[2]="<<c[2]<<endl;
cout<<"R="<<r<<endl;
}

double ff(double x,double y)
{return (-2*x*y);
}
void Euler()
{int n;
double a,b,h,c,k;
cout<<"输入变量x的取值范围"<<endl;
cin>>a>>b;
cout<<"步长h=";
cin>>h;
cout<<"y(0)=";
cin>>c;
n=(b-a)/h;
double x[100],y[100];
x[0]=a;
y[0]=c;
cout<<"x[0]="<<x[0]<<"    "<<"y[0]="<<y[0]<<endl;
for(int i=0;i<=n;i++)
{   k=ff(x[i],y[i]);
	x[i+1]=x[i]+h;
    y[i+1]=y[i]+h*k;
	cout<<"x["<<i+1<<"]="<<x[i+1]<<"    "<<"y["<<i+1<<"]="<<y[i+1]<<endl;
}
}


int main()
{string strChoice;
cout<<"\t\t***数值分析程序,欢迎您使用!***\n";
do{cout<<"\t\t作者04信科1班王宇亮,学号I04630123\n";
cout<<"\t\t1.二分法求解；\n";
cout<<"\t\t2.牛顿法求解；\n";
cout<<"\t\t3.高斯消去法求解；\n";
cout<<"\t\t4.雅可比迭代法求解；\n";
cout<<"\t\t5.拉格朗日插值；\n";
cout<<"\t\t6.牛顿插值；\n";
cout<<"\t\t7.最小二乘法拟合；\n";
cout<<"\t\t8.龙贝格方法计算积分；\n";
cout<<"\t\t9.欧拉方法求解初值问题；\n";
cout<<"\t\t0.退出程序；\n";
cout<<"请您选择操作：\n";
	cin>>strChoice;
		if(strChoice=="1")
		twoFF();
		else if(strChoice=="2")
	    Newton();
		else if(strChoice=="3")
	    GSXU();
        else if(strChoice=="4")
        Jacobi();
        else if(strChoice=="5")
        Lagrange();
        else if(strChoice=="6")
        NewtonCZ();
        else if(strChoice=="7")
        zxecf();
        else if(strChoice=="8")
        Romberge();
        else if(strChoice=="9")
        Euler();
		else if(strChoice!="0")
			cout<<"输入错误,请重新选择: \n";
	}while(strChoice!="0");
cout<<"\n\n\t\t谢谢使用,下次再见!\n";
return 0;
}