编程 | 边际效应 - 杨文博的个人博客 - 页码 5

Something About Scanf() In C (1)

一. scanf 函数输入格式中的字符串

scanf 函数输入格式中也可以含有普通字符串，但他的含义是这些字符必须在输入中出现，例如:

int num;
scanf("hello %d", &num);

他的含义是首先要求输入一个hello字符串，然后再输入一个十进制数。注意在等待输入时忽略hello与要输入的数之间的空格，制表符，回车。

因此这两种输入都是正确的: hello 1234 hello1234

二. scanf函数的返回值

程序:

{
int num, result=0;
printf("please input the student’s score: ");
while(result==0) {
    /* 清空输入缓冲区。 */
    fflush(stdin);
    if(result!=1) printf("Please input a digital score: ");
    result=scanf("%d",&num);
}
}

一切OK!

三. scanf函数中一个参数的应用

在 scanf 函数中，我们可以使用 %c 来读取一个字符，使用 %s 读取一个字符串。但是读取字符串时不忽略空格，读字符串时忽略开始的空格，并且读到空格为止，因此我们只能读取一个单词，而不是整行字符串。因此一般使用 fgets 来读取一个字符串。

其实 scanf 函数也可完成这样的功能，而且还更强大。这里主要介绍一个参数：%[] ，这个参数的意义是读入一个字符集合。 [] 是个集合的标志，因此 %[] 特指读入此集合所限定的那些字符，比如 %[A-Z] 是输入大写字母，一旦遇到不在此集合的字符便停止。如果集合的第一个字符是"^"，这说明读取不在 "^" 后面集合的字符，既遇到 "^" 后面集合的字符便停止。注意此时读入的字符串是可以含有空格的。 Eg: 输入一个字符串，这个字符串只含有小写字符。遇到第一个不是小写字符时停止, scanf("%[a-z], str); Eg: 想输入一个字符串，遇到 "." 停止，可设计如下: scanf("%[^.]", str); 使用这个参数，你可以完成许多强大的功能呦!

通常来讲，scanf 函数和他的一些参数并不是很常用，主要是因为:

1. 许多系统的 scanf 函数都有漏洞。 (典型的就是TC再输入浮点型时有时会出错)。

2. 用法复杂，容易出错。

3. 编译器作语法分析时会很困难，从而影响目标代码的质量和执行效率。

About "double" in C

昨天发现一个很有趣的现象，在 Turbo C 里 double 类型的变量无法用通常模式进行输入操作，即无法用 scanf() 进行赋值，程序举例：

void main()
{
double a,b,c;
scanf("%f%f%f",&a,&b,&c);
printf("%f%f%f",a,b,c);
}

输出结果均是 0.000000，猜想 C 语言应该没有默认初始值的功能，之所以是 0.000000 可能是保留六位小数的原因，由于数值较小且跟地址有关，输入的数值没有传入到地址中。如果将 double 型改为 float 型，则能正常操作。是何道理？难道 C 里输入函数对 double 不兼容？还是另有其他输入方法？回去要查一下 MSDN，不知道 C 有没有专门的帮助文档。

还有要注意的是 TC2 和 TC3 之间对程序要求的变化，在 TC3 中如果不包含标准输入输出的头文件 stdio.h，程序中使用 scanf 和 printf 会报错，而且主函数也必须声明类型，不然会有警告。难道是 ANSI C 的标准变了？还是 TC 为了完善自己只适用于 TC？看来还是要注意编译软件的兼容性问题。 C 语言在好多地方还是不如 C++ 呀！

PS: 输出结果是 0.000000 的原因是因为 scanf() 如果取不到值就会把变量赋零。

2009年6月11日15:40
只因当时年纪小呀，这篇博客现在让我看都快要笑死了。为了不误导读者，特更正如下：
1. 之所以读取错误，是因为本文中程序是错的，%f 匹配浮点类型，%lf 才能匹配 double 类型，这不是 scanf 的错。
2. TC3 对程序的检查是正确的，C 标准确实有变化，我当时是受一本错误的 C 语言教科书误导。
3. C 虽然很多地方不如 C++，但是 C 有简洁性的长处，C++ 是不能比的。

Newton均匀插值多项式

#include <iostream.h>
#include <fstream.h>
#include <stdlib.h>

#define Q(a,b) Q[(a)*(a+1)/2+b]

ifstream infile("in.txt");
ofstream outfile("out.txt");

int n,i,j;
double *x,*root,*Q,*b;

void input()
{
infile>>n;
x=new double[n];
root=new double[n-1];
b=new double[n];
Q=new double[n*(n+1)/2];

for(i=0;i<n;i++)
{
infile>>x[i]>>Q(i,0);
}
}

void setroot()
{
double t;
for(i=0;i<n;i++)
for(j=i;j<n;j++)
if(x[j]<x[i])
{
t=x[j];
x[j]=x[i];
x[i]=t;
t=Q(j,0);
Q(j,0)=Q(i,0);
Q(i,0)=t;
}
for(i=0;i<n-1;i++)
{
root[i]=(x[i]+x[i+1])/2;
}

}

void Newton()
{
for(i=1;i<n;i++)
for(j=1;j<=i;j++)
Q(i,j)=(Q(i,j-1)-Q(i-1,j-1))/(x[i]-x[i-j]);
b[n-1]=Q(n-1,n-1);
}

void f(double x0)
{
for(i=n-1;i>0;i--)
b[i-1]=Q(i-1,i-1)+b[i]*(x0-x[i-1]);
outfile<<"f("<<x0<<")的近似值为:"<<b[0]<<endl;
}

void main(void)
{
input();
setroot();
Newton();
outfile<<"差分表为："<<endl;
for(i=0;i<n;i++)
{ outfile<<i<<’ ’<<x[i]<<’ ’;
for(j=0;j<=i;j++)
outfile<<Q(i,j)<<’ ’;
outfile<<endl;
}

for(j=0;j<n-1;j++)
f(root[j]);

delete x,b,Q,root;
infile.close();
outfile.close();
}

input===============

7
-3 -4.5
6 2.1
-1 -3.2
5 2.8
-2 -1.8
3 1.7
0 0.5

output==============

差分表为：
0 -3 -4.5
1 -2 -1.8 2.7
2 -1 -3.2 -1.4 -2.05
3 0 0.5 3.7 2.55 1.53333
4 3 1.7 0.4 -0.825 -0.675 -0.368056
5 5 2.8 0.55 0.03 0.1425 0.116786 0.0606052
6 6 2.1 -0.7 -0.416667 -0.0744444 -0.0309921 -0.0184722 -0.00878638
f(-2.5)的近似值为:-1.06517
f(-1.5)的近似值为:-3.0675
f(-0.5)的近似值为:-1.86002
f(1.5)的近似值为:5.90552
f(4)的近似值为:-0.577778
f(5.5)的近似值为:4.59237

逐次线性

#include<iostream.h>
#include<fstream.H>
#include<math.h>
#include<stdlib.h>

#define q(a,b) q[(a)*(a+1)/2+b]
#define p(a,b) q[a+1+(b)*(b+1)/2]

ifstream infile("in.txt");
ofstream outfile("out.txt");

double *x,*f,*q,root;
int n,i,j,k;

void indata()
{
infile>>n;
x=new double[n];
f=new double[n];
q=new double[(n+2)*(n+1)/2];
for(i=0;i<=n;i++) infile>>x[i]>>f[i];
infile>>root;
}

void neville()
{
for(i=0;i<=n;i++)
{
q(i,0)=f[i];
for(k=1;k<=i;k++)
q(i,k)=((root-x[i-k])*q(i,k-1)-(root-x[i])*q(i-1,k-1))/(x[i]-x[i-k]);
if(fabs(q(i,i)-q(i-1,i-1))<1e-10)
break;
}
if(i>n) i--;
outfile<<"Neville:"<<q(i,i)<<endl;
}
void aitken()
{
for(i=1;i<=n;i++)
{
p(0,i)=f[0]+(root-x[0])*(f[i]-f[0])/(x[i]-x[0]);
for(j=1;j<i;j++)
p(j,i)=p(j-1,j)+(root-x[j])*(p(j-1,i)-p(j-1,j))/(x[i]-x[j]);
if(fabs(p(i,i)-p(i-1,i-1))<1e-10)
break;
}
if(i>n) i--;
outfile<<endl<<"Aitken:"<<p(i-1,i)<<endl;
}

void mine()
{

q(0,0)=f[0];
for(i=1;i<=n;i++)
{ q(i,0)=f[i];
q(i,1)=1/(x[0]-x[i])*((root-x[i])*f[0]-(root-x[0])*f[i]);

for(k=2;k<=i;k++)
q(i,k)=-1/(x[i-k+1]-x[i])*((root-x[i-k+1])*q(i,k-1)-(root-x[i])*q(i-1,k-1));
if(fabs(q(i,i)-q(i-1,i-1))<1e-10)
break;
}
if(i>n) i--;
outfile<<endl<<"mine:"<<q(i,i)<<endl;
}

int main()
{
indata();
neville();
aitken();
mine();
return 0;
}
input===============

4
-2 0.04
-1 0.2
0 1
1 5
2 25
0.5

output==============

Neville:2.04

Aitken:2.04

mine:2.04

Gauss-Jordan列主元消去法

#include <iostream.h>
#include <fstream.h>
#include <stdlib.h>
#include <math.h>
#define a(u,v) a[(u)*n+(v)-w]

void main()
{//-----------------------文件流----------------------------
char filename1[256],filename2[256];
cout<<"Enter the path to input file:";
cin>>filename1;
cout<<"Enter the path to output file:";
cin>>filename2;
ifstream infile(filename1,ios::in|ios::nocreate);
ofstream outfile(filename2);
if(!infile)
{
cout<<"Can not open input file:"<<filename1<<’
’;
exit(1);
}
//----------------------变量定义----------------------------
int n,w,k,i,j;
int *p;

double *a;
double t;
//------------------------输入------------------------------
infile>>n;
w=n+1;

p=new int[w];
a=new double[n*n];

for(i=0;i<n*n;i++)
infile>>a[i];

//------------------------算法------------------------------
for(k=1;k<w;k++) //step 1
{ p[k]=k;
for(i=k;i<w;i++) //step 2
{
if(fabs(a(i,k))>fabs(a(p[k],k))) p[k]=i;
}

if(a(p[k],k)<1e-10) //step 3
{
cout<<"A is singular!"<<endl;
exit(1);
}
if(p[k]!=k) //step 4
for(i=1;i<w;i++)
{
t=a(k,i);
a(k,i)=a(p[k],i);
a(p[k],i)=t;
}
a(k,k)=1/a(k,k); //step 5
for(i=1;i<w;i++)
if(i!=k)
a(i,k)*=(-a(k,k)); //step 6
for(i=1;i<w;i++) //step 7
if(i!=k)
for(j=1;j<w;j++)
if(j!=k)
a(i,j)+=a(i,k)*a(k,j);
for(j=1;j<w;j++) //step 8
if(j!=k)
a(k,j)*=a(k,k);
}
for(k=n;k>0;k--) //step 9
if(k!=p[k])
for(i=1;i<w;i++)
{
t=a(i,k);
a(i,k)=a(i,p[k]);
a(i,p[k])=t;
}
outfile<<"各步主行行号依次为：:";
for(i=1;i<w;i++)
outfile<<p[i]<<’ ’;

outfile<<"
A的逆矩阵为：
"; //step 10
for(i=1;i<w;i++)
{
for(j=1;j<w;j++)
outfile<<a(i,j)<<’ ’;
outfile<<’
’;
}

delete a,p;

infile.close();
outfile.close();
}
input=================

10
1 2 3 4 5 6 7 8 9 10
2 3 4 5 6 7 8 9 10 1
3 4 5 6 7 8 9 10 1 2
4 5 6 7 8 9 10 1 2 3
5 6 7 8 9 10 1 2 3 4
6 7 8 9 10 1 2 3 4 5
7 8 9 10 1 2 3 4 5 6
8 9 10 1 2 3 4 5 6 7
9 10 1 2 3 4 5 6 7 8
10 1 2 3 4 5 6 7 8 9

output=================

各步主行行号依次为：:10 9 8 7 6 6 7 8 9 10
A的逆矩阵为：
-0.0981818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.101818
0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.101818 -0.0981818
0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.101818 -0.0981818 0.00181818
0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.101818 -0.0981818 0.00181818 0.00181818
0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.101818 -0.0981818 0.00181818 0.00181818 0.00181818
0.00181818 0.00181818 0.00181818 0.00181818 0.101818 -0.0981818 0.00181818 0.00181818 0.00181818 0.00181818
0.00181818 0.00181818 0.00181818 0.101818 -0.0981818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818
0.00181818 0.00181818 0.101818 -0.0981818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818
0.00181818 0.101818 -0.0981818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818
0.101818 -0.0981818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818 0.00181818

三对角方程组的追赶法求解方程组的解

#include <iostream.h>
#include <fstream.h>
#include <stdlib.h>

int main()
{
char /*filename1[256],*/filename2[256];
// cout<<"Enter the path to input file:";
// cin>>filename1;
cout<<"Enter the path to output file:";
cin>>filename2;
// ifstream infile(filename1,ios::in|ios::nocreate);
ofstream outfile(filename2);
// if(!infile)
// {
// cout<<"Can not open input file:"<<filename1<<’
’;
// exit(1);
// }
int n,w,k;
// infile>>n;
n=50;
w=n+1;
double *a,*b,*d,*c,*p,*q,*x,*y,*Iy,*Ix;

a=new double[w];
d=new double[w];
c=new double[w];
b=new double[w];
p=new double[w];
q=new double[w];
x=new double[w];
y=new double[w];
Iy=new double[w];
Ix=new double[w];

//输入======================================================================================

int i;
infile>>d[1]>>c[1];
// infile.seekg((n-3)*sizeof(double));
// infile.seekg((n-3)*sizeof(int),ios::cur);
for(k=0;k<n-3;k++)
infile>>a[0];
infile>>a[1];

for(k=2;k<n;k++)
{
// infile.seekg((k-2)*sizeof(int),ios::cur);
for(i=0;i<k-2;i++)
infile>>a[0];
infile>>a[k]>>d[k]>>c[k];
// infile.seekg((n-k-1)*sizeof(int),ios::cur);
for(i=0;i<n-k-1;i++)
infile>>a[0];
}
infile>>c[n];
// infile.seekg((n-3)*sizeof(int),ios::cur);
for(i=0;i<n-3;i++)
infile>>a[0];
infile>>a[n]>>d[n];
for(k=1;k<w;k++)
infile>>b[k];
*/
for(k=1;k<w;k++)
{
a[k]=1;
c[k]=1;
d[k]=4;
b[k]=k+1;
}
d[1]=2;d[n]=4;

//算法===========================================================================
if(d[1]<1e-8)
{
cout<<"Method failed!";
exit(1);
}
p[1]=d[1];
q[1]=c[1]/d[1];
for(k=2;k<n-1;k++)
{
p[k]=d[k]-a[k]*q[k-1];
if(p[k]<1e-8)
{
cout<<"Method failed!";
exit(2);
}
q[k]=c[k]/p[k];
}
p[n-1]=d[n-1]-a[n-1]*q[n-2];
if(p[n-1]<1e-8)
{
cout<<"Method failed!";
exit(3);
}
y[1]=b[1]/d[1]; Iy[1]=-a[1]/d[1];
for(k=2;k<n-1;k++)
{
y[k]=(b[k]-a[k]*y[k-1])/p[k];
Iy[k]=(-a[k]*Iy[k-1])/p[k];
}
y[n-1]=(b[n-1]-a[n-1]*y[n-2])/p[n-1];
Iy[k]=(-c[n-1]-a[k]*Iy[k-1])/p[k];
x[n-1]=y[n-1]; Ix[n-1]=Iy[n-1];
for(k=n-2;k>0;k--)
{
x[k]=y[k]-q[k]*x[k+1];
Ix[k]=Iy[k]-q[k]*Ix[k+1];
}
x[n]=(b[n]-c[n]*x[1]-a[n]*x[n-1])/(c[n]*Ix[1]+a[n]*Ix[n-1]+d[n]);
for(k=1;k<n;k++)
x[k]+=Ix[k]*x[n];
outfile<<"The roots are:(";
for(k=1;k<n;k++)
outfile<<x[k]<<’,’;
outfile<<x[n]<<");
";
delete a,d,c,p,q,x,y,Iy,Ix;
// infile.close();
outfile.close();
return 0;
}
output============

The roots are:(-6.4282,2.31175,0.18121,0.963411,0.965146,1.17601,1.33083,1.50067,1.66649,1.83338,1.99999,2.16667,2.33333,2.5,2.66667,2.83333,3,3.16667,3.33333,3.5,3.66667,3.83333,4,4.16667,4.33333,4.5,4.66667,4.83333,5,5.16667,5.33333,5.5,5.66667,5.83333,6,6.16667,6.33333,6.5,6.66666,6.83334,6.99997,7.16677,7.33293,7.5015,7.66108,7.85418,7.92219,8.45706,7.24957,12.5447);

修改的LDLt分解

#include <iostream.h>
#include <fstream.h>
#include <stdlib.h>

int main()
{
char filename1[256],filename2[256];
cout<<"Enter the path to input file:";
cin>>filename1;
cout<<"Enter the path to output file:";
cin>>filename2;
ifstream infile(filename1,ios::in|ios::nocreate);
ofstream outfile(filename2);
if(!infile)
{
cout<<"Can not open input file:"<<filename1<<’
’;
exit(1);
}

int n,w,k,i,j;
infile>>n;
w=n+1;

double *a,*x,*b,*g;
double t;

a=new double[n*n];
g=new double[w];
x=new double[w];
b=new double[w];
for(i=0;i<n*n;i++)
infile>>a[i];
for(i=1;i<w;i++)
infile>>b[i];

if(a[0]<1e-8) //Step1
{
cout<<"Method failed!";
exit(2);
}
/*其实这里可以定义：#define a(int u,int v) a[u*n-n+v-1]
这样引用起来比较直观!但是运算是一样的.*/
g[1]=a[n]; //Step2
a[n]=g[1]/a[0];
a[n+1]-=g[1]*a[n];
if(a[n+1]<1e-8) //Step3
{
cout<<"Method failed!";
exit(3);
}
for(i=3;i<w;i++) //Step4
{
g[1]=a[i*n-n]; //Step5
for(j=2;j<i;j++) //Step6
{
for(t=0,k=1;k<j;k++)
t+=g[k]*a[j*n-n+k-1];
g[j]=a[i*n-n+j-1]-t;
}
a[i*n-n]=g[1]/a[0]; //Step7
for(j=2;j<i;j++) //Step8
a[i*n-n+j-1]=g[j]/a[j*n-n+j-1];
for(t=0,k=1;k<i;k++) //Step9
t+=g[k]*a[i*n-n+k-1];
a[i*n-n+i-1]-=t;
if(a[i*n-n+i-1]<1e-8) //Step10
{
cout<<"Method failed!";
exit(4);
}
}
for(i=2;i<w;i++) //Step11
{
for(t=0,k=1;k<i;k++)
t+=a[i*n-n+k-1]*b[k];
b[i]-=t;
}
x[n]=b[n]/a[n*n-1]; //Step12
for(i=n-1;i>0;i--) //Step13
{
for(t=0,k=i+1;k<w;k++)
t+=a[k*n-n+i-1]*x[k];
x[i]=b[i]/a[i*n-n+i-1]-t;
}

for(i=1;i<w;i++)
outfile<<"g("<<i<<")="<<g[i]<<’
’;

outfile<<"y=(";
for(i=1;i<n;i++)
outfile<<b[i]<<’,’;
outfile<<b[n]<<");
";

outfile<<"diag(D)=(";
for(i=1;i<n;i++)
outfile<<a[i*n-n+i-1]<<’,’;
outfile<<a[n*n-1]<<");
";

outfile<<"L=
";
for(i=1;i<w;i++)
{ for(j=1;j<w;j++)
{ if(j<=i)
outfile<<’ ’<<a[i*n-n+j-1];
else outfile<<’ ’<<’0’;
}
outfile<<’
’;
}

outfile<<"The roots are:("; //Step14
for(i=1;i<n;i++)
outfile<<x[i]<<’,’;
outfile<<x[n]<<");";
delete a;
delete g;
delete b;
delete x;
infile.close();
outfile.close();
return 0;
}
input=============

10
1 1 1 1 1 1 1 1 1 1
1 2 2 2 2 2 2 2 2 2
1 2 3 3 3 3 3 3 3 3
1 2 3 4 4 4 4 4 4 4
1 2 3 4 5 5 5 5 5 5
1 2 3 4 5 6 6 6 6 6
1 2 3 4 5 6 7 7 7 7
1 2 3 4 5 6 7 8 8 8
1 2 3 4 5 6 7 8 9 9
1 2 3 4 5 6 7 8 9 10
8 2 -12 4 -45 6 -5 9 10 30

output=============

g(1)=1
g(2)=1
g(3)=1
g(4)=1
g(5)=1
g(6)=1
g(7)=1
g(8)=1
g(9)=1
g(10)=-6.27744e+066
y=(8,-6,-14,16,-49,51,-11,14,1,20);
diag(D)=(1,1,1,1,1,1,1,1,1,1);
L=
1 0 0 0 0 0 0 0 0 0
1 1 0 0 0 0 0 0 0 0
1 1 1 0 0 0 0 0 0 0
1 1 1 1 0 0 0 0 0 0
1 1 1 1 1 0 0 0 0 0
1 1 1 1 1 1 0 0 0 0
1 1 1 1 1 1 1 0 0 0
1 1 1 1 1 1 1 1 0 0
1 1 1 1 1 1 1 1 1 0
1 1 1 1 1 1 1 1 1 1
The roots are:(14,8,-30,65,-100,62,-25,13,-19,20);

按比例列选主元Gauss消去法求线性方程组的解

实习二

用按比例列选主元Gauss消去法求下列线性方程组Ax=b的解。

系数矩阵A：
15 7 6 5 1
7 10 8 7 2
6 8 10 9 3
5 7 9 10 4
1 2 3 4 5

右端向量b：

-4.0000
21.0000
13.2000
5.0000
-9.0000

要求：输出每次选取的主元的行号，解向量x.
输出时要有一定的文字说明，以使输出结果一目了然。

#include <iostream.h>
#include <fstream.h>
#include <stdlib.h>
#include <math.h>
void main()
{
char filename1[256],filename2[256];
cout<<"Enter the path to input file:";
cin>>filename1;
cout<<"Enter the path to output file:";
cin>>filename2;
ifstream infile(filename1,ios::in|ios::nocreate);
ofstream outfile(filename2);
if(!infile){
cout<<"Can’t open input file:"<<filename1<<’
’;
exit(1);
}

int k,m,i,j,n,w;
infile>>n;
w=n+1;

double *p,*x;
double t;

p=new double[w*w];
x=new double[w];
for(i=0;i<w;i++) p[i]=0;
for(i=w;i<w*w;i++)
infile>>p[i];

for(k=1;k<n;k++) //step 1

{
t=p[k*w+k-1];
m=k;
for(i=k;i<w;i++) //step 2
{ if(fabs(p[i*w+k-1])>fabs(t))
{ t=p[i*w+k-1];
m=i;
}
}
if(t<1e-8) //step 3
{ cout<<"A is singular"<<endl;
exit(2);
}
outfile<<"第"<<k<<"列的主元行号为："<<m<<’
’;
//----------------------------------------
if(m!=k) //step 4
for(j=k;j<w+1;j++)
{ t=p[k*w+j-1];
p[k*w+j-1]=p[m*w+j-1];
p[m*w+j-1]=t;
}
for(i=k+1;i<w;i++) //step 5
{ p[i*w+k-1]=p[i*w+k-1]/p[k*w+k-1]; //step 6
for(j=k+1;j<w+1;j++) //step 7
p[i*w+j-1]=p[i*w+j-1]-p[i*w+k-1]*p[k*w+j-1];
}

}
x[n]=p[n*w+n]/p[n*w+n-1]; //step 8

for(k=n-1;k>0;k--) //step 9
{ t=0;
for(j=k+1;j<n+1;j++) t+=p[k*w+j-1]*x[j];
x[k]=(p[k*w+n]-t)/p[k*w+k-1];
}
outfile<<"方程组的根为：("; //step 10
for(i=1;i<n;i++)
outfile<<x[i]<<’,’;
outfile<<x[n]<<’)’;
infile.close();
outfile.close();

}
input===========

5
15 7 6 5 1 -4.0000
7 10 8 7 2 21.0000
6 8 10 9 3 13.2000
5 7 9 10 4 5.0000
1 2 3 4 5 -9.0000

output==========

第1列的主元行号为：1
第2列的主元行号为：2
第3列的主元行号为：3
第4列的主元行号为：4
方程组的根为：(-2.01506,3.73891,2.00518,-1.87662,-2.59436)

三次方程求法

#include <iostream.h>
#include <math.h>
#include <iomanip>
#include <fstream.h>
#include <stdlib.h>

#define f(x) (1*x*x*x-3*x*x-6*x+8)

double bisection(double x,double y,double z)
{
double a=x,b=y,p;
while(((b-a)/2)>=z)
{ p=(a+b)/2;
if(f(p)*f(b)<=0)
a=p;
else b=p;
}
return p;

}

void main()
{
double a[5];
char filename[256];
cout.setf(ios::fixed);
cout<<"请输入误差容限:";
cin>>a[0];
cout<<"输入目的文件全路径：";
cin>>filename;
ofstream outfile(filename);
a[2]=(3-sqrt(21))/2;
a[3]=(3+sqrt(21))/2;
a[1]=a[2];
a[4]=a[3];
for(int i=1; ;i++)
{
if(f(a[1])>=0) a[1]--;

else if(f(a[4])<=0) a[4]++;
else break;
}
for(i=1;i<4;i++)

outfile<<"函数f(x)在区间("<<a[i]<<’,’<<a[i+1]<<")上的根为"<<bisection(a[i],a[i+1],a[0])<<endl;

cout<<"已经将结果写入："<<filename<<endl;
cout<<"误差在"<<a[0]<<"之内."<<endl;
outfile.close();
}

input=============

f:output.txt //可以是任何路径
1e-7

output============

函数f(x)在区间(-2.79129,-0.791288)上的根为-2
函数f(x)在区间(-0.791288,3.79129)上的根为1
函数f(x)在区间(3.79129,4.79129)上的根为4

变步长的Simpson求积法

#include<iostream.h>
#include<fstream.h>
#include<math.h>
#include<stdlib.h>

ifstream infile("in.dat");
ofstream outfile("out.out ");

double a[3],b[3];
int i,m[3],k;

void indata()
{
for(i=0;i<3;i++)
infile>>a[i]>>b[i]>>m[i];
m[i]/=2;
}

double cal0(double x)
{
return (6*exp(-x*x));
}

double cal1(double x)
{
return (3*log(x));
}

double cal2(double x)
{
return (sin(x)/x);
}

void simpson(double a,double b,int m,double (*g)(double))
{
double h[3],q[3],p[3],s[3],t;

h[0]=(b-a)/(2*m);
t=0;
for(i=1;i<m+1;i++)
t+=(g(a+(2*i-1)*h[0])+g(a+2*i*h[0]));
p[0]=g(a)-g(b)+2*t;
t=0;
for(i=1;i<m+1;i++)
t+=g(a+(2*i-1)*h[0]);
s[0]=h[0]*(p[0]+2*t)/3;

h[1]=h[0]/2;
q[0]=0;
for(i=1;i<2*m+1;i++)
q[0]+=g(a+(2*i-1)*h[1]);
s[1]=h[1]*(p[0]+4*q[0])/3;
h[2]=h[1]/2;
q[1]=0;
for(i=1;i<pow(2,2)*m+1;i++)
q[1]+=g(a+(2*i-1)*h[2]);
p[1]=p[0]+2*q[0];
s[2]=h[2]*(p[1]+4*q[1])/3;
for(k=2;;k++)
{
h[(k+1)%3]=(b-a)/(pow(2,k+2)*m);
q[k%3]=0;
for(i=1;i<pow(2,k+1)*m+1;i++)
q[k%3]+=g(a+(2*i-1)*h[(k+1)%3]);
p[k%3]=p[(k-1)%3]+2*q[(k-1)%3];
s[(k+1)%3]=h[(k+1)%3]*(p[k%3]+4*q[k%3])/3;
if(s[(k+1)%3]-s[k%3]<1e-7)
{ t=s[(k+1)%3];
break;
}
}
outfile<<t<<endl;
outfile<<"积分区间数目为:"<<pow(2,k+1)*m<<endl;
}

void main()
{
indata();

outfile<<"I("<<1<<")=";
simpson(a[0],b[0],m[0],cal0);

outfile<<"I("<<2<<")=";
simpson(a[1],b[1],m[1],cal1);

outfile<<"I("<<3<<")=";
simpson(a[2],b[2],m[2],cal2);

infile.close();
outfile.close();
}

input==================

0 1 2
1 2 2
1 4 2

output=================

I(1)=4.48094
积分区间数目为:16
I(2)=1.15888
积分区间数目为:32
I(3)=0.81212
积分区间数目为:32