用C语言计算微积分？

我给一楼加的注释以及修改：

#include<stdio.h>

#include<math.h>

#define ARRAYBOUND 10001

void main()

{

int i = 0 //辅助变量，最常见那种

int n = 0 //将所求定积分函数曲线在x轴方向，平均分成n等分n越大,结果越精确不过限于此算法限制n<ARRAYBOUND,否则溢出.

float x[ARRAYBOUND]//ARRAYBOUND维浮点数组，存放离散的x坐标值

float y[ARRAYBOUND]//ARRAYBOUND维浮点数组，存放每个x坐标对应的函数值x[i],y[i]满足y[i]=f(x[i]),f是你要求定积分的函数

float x0 = 0.0 //定积分下限

float xn = 0.0 //定积分上限

float h = 0.0 //面积微元宽度

float J = 0.0 //辅助变量

/*f=x^3*/ //这里说明要求定积分的是函数f(x)=x*x*x(y等于x的立方,x^3是vb的写法)

// printf("input x0,xn,n:")

printf("请分别输入下限(x0),上限(xn),精度(n):")

scanf("%f",&x0)

scanf("%f",&xn)

scanf("%d",&n)

h=(xn-x0)/n//将函数图形在x方向平分成n份,h是每个面积微元的宽度

x[0]=x0 //将积分下限赋值给x[0]

for(i=0i<=n &&n<ARRAYBOUNDi++)

{

x[i]=x[0]+i*h //计算n个离散的横坐标值，存入x[]数组

y[i]=(float)pow(x[i],3)//计算n个横坐标对应的函数值，存入y[]数组。在此可以改变要求积分的函数

}

// J=0.0

for(i=0i<ni++)

{

//J=J+y[i]+y[i+1]

J+=y[i]//将所有纵坐标值代数相加，存入J

}

//J=J*h/2.0

J=J*h//所有微元面积一次求解，因为∑h*y[i]=h*∑y[i]

printf("\nn=%d \n所求定积分值是: %f\n",n,J)

}

我将//J=J+y[i]+y[i+1]改为J+=y[i]将//J=J*h/2.0改为J=J*h只是帮助lz理解

其实，这两种表达在理论上是等价的，不过我发现修改后，在n同样大小的情况下，结果的精度有一点点下降，还真不知为什么？？？

这样的话lz应该能理解了吧，其实一楼的算法还有不少值得改进的地方，希望lz能有所突破！！

#include <stdio.h>

main()

{

float x,y,ji_fen,wei_fen

for(x=0x<2x=x+0.01)

{

y=x*x

ji_fen=ji_fen+y*0.01

wei_fen=(y-(x-0.01)*(x-0.01))/0.01

printf("x=%f y=%f ji_fen=%f wei_fen=%f\n",x,y,ji_fen,wei_fen)

}

}

这是C语言的程序用的的积分分离，增量式算法你可以看看~在我用的片子是89C52

/*********************************************************/

/* 名称: float PIDprocess1*/

/* 功能: PID adjust */

/* 说明: */

/* 调用: */

/* 输入: float xdata *Yn, float xdata *Rn */

/* 返回值: deltaPn */

/*********************************************************/

float PIDprocess1()

{

int data E_0

float data deltaPn,deltaPi,deltaPp,deltaPd,PsumCopy

E_0=SetTemperature1-CurrentTemperature1

if(abs(E_0)>Emax)

{

deltaPp=(float)Kp*(E_0-E_11)

deltaPd=(float)Kd*(E_0-2*E_11+E_21)

// if(deltaPd>=dPdmax) deltaPd=0

deltaPn=deltaPp+deltaPd

}

else

{

if(abs(E_0)>E0)

{

deltaPi=(float)Ki*E_0*(Emax-abs(E_0))/(Emax-E0)

}

else

{//小误差时的处理

deltaPi=(float)Ki*E_0

// if(fabs(deltaPi)<dPimin) deltaPi=0

}

deltaPp=(float)Kp*(E_0-E_11)

deltaPd=(float)Kd*(E_0-2*E_11+E_21)

// if(deltaPd>=dPdmax) deltaPd=0

deltaPn=deltaPp+deltaPi+deltaPd

}

/* if(fabs(deltaPn)>dPmax)

{

if(deltaPn>0) deltaPn=dPmax

else deltaPn=-dPmax

} */

Psum1+=deltaPn

PsumCopy=Psum1

if(PsumCopy>Pmax) PsumCopy=Pmax

if(PsumCopy<Pmin) PsumCopy=Pmin

E_21=E_11

E_11=E_0

return(PsumCopy)

}

/*********************************************************/

/* 名称: float PIDprocess2 */

/* 功能: PID adjust */

/* 说明: */

/* 调用: */

/* 输入: float xdata Yn, float xdata Rn*/

/* 返回值: Psum2 */

/*********************************************************/

float PIDprocess2()

{

int data E_0

float data deltaPn,deltaPi,PsumCopy

E_0=SetTemperature2-CurrentTemperature2

if(abs(E_0)>Emax)

{

deltaPn=(float)Kp*(E_0-E_12)+(float)Kd*(E_0-2*E_12+E_22)

}

else

{

if(abs(E_0)>E0)

{

deltaPi=(float)Ki*E_0*(Emax-abs(E_0))/(Emax-E0)

}

else

{//小误差时的处理

deltaPi=(float)Ki*E_0

}

deltaPn=(float)Kp*(E_0-E_12)+deltaPi+(float)Kd*(E_0-2*E_12+E_22)

}

/* if(fabs(deltaPn)>dPmax)

{

if(deltaPn>0) deltaPn=dPmax

else deltaPn=-dPmax

}*/

Psum2+=deltaPn

PsumCopy=Psum2

if(PsumCopy>Pmax) PsumCopy=Pmax

if(PsumCopy<Pmin) PsumCopy=Pmin

E_22=E_12

E_12=E_0

return(PsumCopy)

}

---