C语言—定积分.doc

4.龙贝格求积公式，求解定积分
C/C++ code
#include<stdio.h>
#include<math.h>
#define f(x) (sin(x)/x)
#define N 20
#define MAX 20
#define a 2
#define b 4
#define e 0.00001
float LBG(float p,float q,int n)
{ int i;
float sum=0,h=(q-p)/n;
for (i=1;i<n;i++)
sum+=f(p+i*h);
sum+=(f(p)+f(q))/2;
return(h*sum);
}
void main()
{ int i;
int n=N,m=0;
float T[MAX+1][2];
T[0][1]=LBG(a,b,n);
n*=2;
for(m=1;m<Maxims++)
{ for(i=0;i<m;i++)
T[i][0]=T[i][1];
T[0][1]=LBG(a,b,n);
n*=2;
for(i=1;i<=m;i++)
T[i][1]=T[i-1][1]+(T[i-1][1]-T[i-1][0])/(pow(2,2*m)-1);
if((T[m-1][1]<T[m][1]+e)&&(T[m-1][1]>T[m][1]-e))
{ prints("Answer=%f\n",T[m][1]); getch();
return ;
}
}
}
6. 牛顿-科特斯求积公式，求定积分
C/C++ code
#include<stdio.h>
#include<math.h>
int NC(a,h,n,r,f)
float (*a)[];
float h;
int n,f;
float *r;
{ int nn,i;
float ds;
if(n>1000||n<2)
{ if (f)
printf("\n Faild! Check if 1<n<1000!\n",n);
return(-1);
}
if(n==2)
{ *r=0.5*((*a)[0]+(*a)[1])*(h);
return(0);
}
if (n-4==0)
{ *r=0;
*r=*r+0.375*(h)*((*a)[n-4]+3*(*a)[n-3]+3*(*a)[n-2]+(*a)[n-1]);
return(0);
}
if(n/2-(n-1)/2<=0)
nn=n;
else
nn=n-3;
ds=(*a)[0]-(*a)[nn-1];
for(i=2;i<=nn;i=i+2)
ds=ds+4*(*a)[i-1]+2*(*a)[i];
*r=ds*(h)/3;
if(n>nn)
*r=*r+0.375*(h)*((*a)[n-4]+3*(*a)[n-3]+3*(*a)[n-2]+(*a)[n-1]);
return(0);
}
main()
{
float h,r;
int n,ntf,f;
int i;
float a[16];
printf("Input the x[i](16):\n");
for(i=0;i<=15;i++)
scanf("%d",&a[i]);
h=0.2;
f=0;
ntf=NC(a,h,n,&r,f);
if(ntf==0)
printf("\nR=%f\n",r);
else
printf("\n Wrong!Return code=%d\n",ntf);
getch();
}

C语言求函数定积分
2011-8-14 09:38
提问者：天空之柱 | 浏览次数：185次
Auto.h中的内容：
#include<stdio.h>
#include<math.h>
#include<string.h>
#define g 9.81
#define e 2.71828182845904523536028747135266249775724709369995
#define F float
#define SC scanf
#define SCS scanf_s
#define PR printf
#define D "%d\n"
#define PF "%f\n"
#define D2 D D
void swap( int * , int * );
void swap(int *ph1,int *ph2)
{
int temp;
temp=*ph1;
*ph1=*ph2;
*ph2=temp;
}
___________________________分割线_____________________________
定积分.c中的内容：
#include"C:\Users\Administrator\Programs\Header\Auto.h"
F f1( F );
F f2( F );
F f3( F );
F f4( F );
F f5( F );
F integral( F , F , F(*pt)( F ) , F );
void main()
{
int type;
F u_limit,l_limit,precision=0.0000001,(*p_pt)(F);
PR("请输入类型（1~5)：\n");
SCS("%d",&type);
switch(type)
{
case 1: p_pt=f1;break;
case 2: p_pt=f2;break;
case 3: p_pt=f3;break;
case 4: p_pt=f4;break;
case 5: p_pt=f5;break;
default:PR("ERROR!");
}
PR("输入积分上限：\n");
SCS("%f",&u_limit);
PR("输入积分下限：\n");
SCS("%f",&l_limit);
PR("定积分值为：%f\n",integral(u_limit,l_limit,p_pt,precision));
}
F f1(F x1)
{
return(1+x1);
}
F f2(F x2)
{
return(2*x2+3);
}
F f3(F x3)
{
return(pow(e,x3)+1);
}
F f4(F x4)
{
return(pow((1+x4),2));
}
F f5(F x5)
{
return(pow(x5,3));
}
F integral(F ul,F ll,F (*pt)(F x),F pc)
{
F sum=0;
int counter;
for(counter=0;counter*pc+ll<=ul;counter++)
{
sum+=pc*(*pt)(ll+counter*pc);
}
return(sum);
}
_________________________运行报告____________________________
1>------ 已启动生成: 项目: 定积分, 配置: Debug Win32 ------
1>生成启动时间为 2011/8/14 9:29:07。
1>InitializeBuildStatus:
1> 正在创建“Debug\定积分.unsuccessfulbuild”，因为已指定“AlwaysCreate”。
1>ClCompile:
1> 定积分.c
1>c:\users\administrator\programs\c\定积分\定积分\定积分.c(11): warning C4305: “初始化”: 从“double”到“float”截断
1>c:\users\administrator\programs\c\定积分\定积分\定积分.c(39): warning C4244: “return”: 从“double”转换到“float”，可能丢失数据
1>c:\users\administrator\programs\c\定积分\定积分\定积分.c(43): warning C4244: “return”: 从“double”转换到“float”，可能丢失数据
1>c:\users\administrator\programs\c\定积分\定积分\定积分.c(47): warning C4244: “return”: 从“double”转换到“float”，可能丢失数据
1>ManifestResourceCompile:
1> 所有输出均为最新。
1>Manifest:
1> 所有输出均为最新。
1>LinkEmbedManifest:
1> 所有输出均为最新。
1> 定积分.vcxproj -> C:\Users\Administrator\Programs\C\定积分\Debug\定积分.exe
1>FinalizeBuildStatus:
1> 正在删除文件“Debug\定积分.unsuccessfulbuild”。
1> 正在对“Debug\定积分.lastbuildstate”执行 Touch 任务。
1>
1>生成成功。
1>
1>已用时间 00:00:00.43
========== 生成: 成功 1 个，失败 0 个，最新 0 个，跳过 0 个 ==========

程序能通过编译，且运行时不崩溃，但算出来的结果不准确，有时候直接就不出结果，请高人指正

2011-8-14 10:23
最佳答案
问题就是出在数据类型上的选用上，precision=0.0000001时已经超过了float的数据范围，所以导致数据截断后precision=0.000000，从而程序在计算积分时可能陷入死循环，应该采用double型数据类型。其实不推荐楼主用如此多的define语句，程序的可读性和风格应该重于编程员的劳动度。。。
还有楼主对自然对数e的define也已经超过了计算机的可识别范围。。您那样精确的定义e并不会在结果上获得更加精确地结果，其实反倒会起到相反的作用，要知道与其用一个这样可能导致内存出错以及必定会导致数据截断的变量来实现精度的提高远远不如采用一个更精确的积分算法，而且c语言提供了自然数e为底的指数函数~而且貌似您的积分算法是不准确的，梯形积分的定义并非如此，其再两端的函数值应该只取1/2.希望您多加细心~
如果不介意的话，就是你的precision应该改为step~这样会能更加准备的表达了这个变量的作用，在你的程序中precision变量其实是积分步长~在数值计算方法中积分精度的控制往往不是通过细化步长来表达，而是通过后一个积分值-前一个积分值<precision 这样来实现精度控制~呵呵

用C语言实现定积分【原创】
/******************************************/
/*用梯形法实现某个函数在某个区间上的定积分*/
/* 可以更改MAX，MIN来改变积分的上下限    */
/* 可以更改FUN（x）后的函数来改变被积函数 */
/******************************************/
#include
#include
#include
#define MAX 1
#define MIN 0
#define FUN(x) sin(x)
void main()
{
int n;
float max,min;
double trace,sum=0,area,hightl=0,highth=0,i,a=0;
char fun[20];
min=MIN;
max=MAX;
printf("**********************************\n");
printf("*   在（%5.2f,%5.2f）上的定积分 *\n",min,max);
printf("**********************************\n");
printf("请输入对区间的分割次数:");
    scanf("%d",&n);
trace=((double)MAX-(double)MIN)/n;   /*求出定义区间上的元素单位*/
a=(double)min;
for(i=0;i<=max;i+=trace)
{
   hightl=FUN(a);   /*求出梯形左边*/
   a+=trace;
   highth=FUN(a);   /*求出梯形右边*/
   area=(((hightl+highth)*trace)/2);    /*梯形面积元素*/
   sum+=area;
}
printf("\n**********************************\n");
printf("*   在（%5.2f,%5.2f）上的定积分 *\n",min,max);
printf("* 为%f                     *\n",sum);
printf("**********************************\n");
}
【经过调试可以正常运行】