| 网站首页 | 业界新闻 | 小组 | 威客 | 人才 | 下载频道 | 博客 | 代码贴 | 在线编程 | 编程论坛
欢迎加入我们,一同切磋技术
用户名:   
 
密 码:  
共有 423 人关注过本帖
标题:求助:关与C语言的计算程序,解卡尔丹(x**3+p*X+q=0)一元三次的实根,结果 ...
只看楼主 加入收藏
momokill
Rank: 1
等 级:新手上路
帖 子:5
专家分:0
注 册:2007-10-28
收藏
 问题点数:0 回复次数:0 
求助:关与C语言的计算程序,解卡尔丹(x**3+p*X+q=0)一元三次的实根,结果存入xa[

求助:关与C语言的计算程序,解卡尔丹(x**3+p*X+q=0)一元三次的实根,结果存入xa[],出不来正确结果。
/**************晚上**************/
#include<stdio.h>
#include<math.h>
#include<float.h>
#define Real double
#define PI 3.14159265354

#define Real double
#define Real_MIN DBL_MIN
/*double PI2(double xx){return 2*PI;}*/

/*解卡尔丹(x**3+p*X+q=0)一元三次的实根,结果存入xa[],实
根个数为*n*/
void getslovp3l(real p,real q,real xa,int n)
{real d,r,s;
p*=(1./.3);
q*=0.5;
d=q*q+p*p*p;
if(d>0.)
{*n=1;d=sqrt(d);
xa[0]=(((d-q)>=0.?pow(d-q,1./3.):-pow(q-d,1./3.))+((-q-d)>=0.?pow(-q-d,1./3.):-pow(q+d,1./3./)));

else if(d==0)
{*n=3;
if(p==q==0.)xa[0]=xa[1]=xa[2]=0.;
else {r=(q.0?-pow(q,1./3.):pow(-q,1./3.));
xa[0]=r+r;
xa[1]=xa[2]=-r; }

else
{*n=3;
r=pow(q*q-d,1./6.);
s=(1./3.)*atan2(sqrt(-d),-q);
xa[0]=r*cos(s);xa[0]+=xa[0];
xa[1]=r*cos(s+PI2/3.);xa[1]+=xa[1];
xa[2]=r*cos(s-PI2/3.);xa[2]+=xa[2];


/*解一元二次方程的实根,结果存入x[],实根个数为*n*/
void getslov2l(Real a,Real b,Real c,Real *x,int *n)
{Real d;
d=b*b-4.*a*c;
if(d>0.)
{*n=2;
if(b>=0)x[0]=(-b-sqrt(d))/(2.*a);
else x[0]=(-b+sqrt(d))/(2.*a);
if(fabs(x[0])>Real_MINx[1]=c/(a*x[0]);
else x[1]=-b/a;

else if(d==0){*n=2;x[0]=x[1]=-b/(2.*a);}
else {*n=0;return;}


/*解任意一元三次方程ax**3+b*X**2+c*x+d=0(a!=0)方程的实
根,其中*n取值为0-3,为实根个数*/
void getslovp3l(Real a,Real b,Real c,Real d,real *x,int *n)
{Real h;
b/=a;c/=a;d/=a;
h=b*(1./3.);
getslovp3l(c-b*b*(1./3.),2.*b*b*b*(1./27)-b*c*(1./3)+d,x,n);
x[0]-=h;x[1]-=h;x[2]-=h;

*解ax**3+b*X**2+c*x+d=0一元三次方程的实根,(其中包括a=0退化为低次方程的情况)
其中*n取值为0-3时为实根个数,取值为4时表示所有实数都是该方程的解*/
void getslovp3lall(Real a,Real b,Real c,Real d,real *x,int *n)
{Real h;
if(fabs(a)>Real_MIN)

b/=a;c/=a;d/=a;
h=b*(1./3.);
getslop3l(c-b*b*(1./3.),b*b*b*(2./27.)-b*c*(1./3.)+d,x,n);
x[0]-=h;x[1]-=h;x[2]-=h;

else if(fabs(b)>Real_MIN)


getslov2l(b,c,d,x,n);


else if(fabs(c)>Real_MIN)

*n=1;x[0]=d/c;

else if(d==0.)*n=4;
else *n=0;

mian[]
{double a,b,c,d,x[3];
char str1[10],str2[10];
float aa;
int i,n,tn,tmin,lmin.lmax,tmax;
FILE *fp;
/*scanf("%f%f%f%f,&a,&b,&c,&d);*/
/*scanf("%f,&a);*/
/*scanf("%f,&b);*/
float m0,s0,e0,a0,fpp,fmax,detp,d0,l,lm,af;
float c0,f1,f6,f,g,t,v,det;
float t1,t4,t6,v1,v4,v6;
Real g1,g4,g6;
if((fp=fopen("jcsj.txt","r"))==NULL{puts("cannot openfile jg.txt\n");}
fscanf(fp,"%f%f%f%f%f%f%f%f%f%d%d',&m0,&s0,&d0,&fmax,&e0,&a0,&af,&k,&lm,&tim,&tmax);
fclose(fp);
/*mo 导线计算单位质量(kg/km)*/
/*s0 导线计算截面积(mm~2)*/
/*d0 导线计算直径(mm)*/
/*fmax 导线拉断力(N)*/
/*e0 导线弹性系数(kg/mm~2)*/
/*a0 导线线膨胀系数(1/℃)*/
/*af 风速不均匀系数,绝缘导线af=1*/
/*lmin 导线需计算的最小档距(m)*/
/*lmax 导线需计算的最大档距(m)*/
/*tmin* 导线需计算的最低温度(℃)*/
/*tmax 导线需计算的最高温度(℃)*/
detp=0.;
t1=-5.;t4=20.;t6=20.;v1=0.;v4=25.;v6=0;/*气象条件*/
tn=0;
detp=(detp==0)?fmax/s0:detp;
fmax=(fmax==0)?detp*s0:fmax;
g1=9.80665*m0/(1000.*s0);
c0=(d0>=17)?1.1:1.2;
/*不考虑风压高度变化系数情况下*/
g4=0.6128*v4*v4*c0*d0*af/(s0*1000.);
g6=sqrt(g1*g1+g4*g4);
a=1.;
c=0.;
if((fp=fopen("jp.txt","w"))==NULL)(puts("cannot openfile jg.txt\n");}
fprintf(fp,"\n导线单位质量(kg/km)m=%g",m0);fprintf(fp,"\n 导线截面(mm~2) s=%g",s0);
fprintf(fp,"\n导线直径(mm)d=%g",d0);fprintf(fp,"\n 导线弹性系数(N/mm~2) E=%g",e0);
fprintf(fp,"\n导线线膨胀系数(1/℃)α=%g",a0);fpritf(fp,"\n导线破坏应力&p=%g",dept);
fprintf(fp,"\n最大使用应力detp=%g",detp/k);fprint(fp,"\n安全系数k=%g",k);
fprintf(fp,"\n自重比载个g1=%g",g1);fprint(fp,"\n风比载g4=%g,g4);
fprint(fp,"\n无冰综合比载g6=%g,g6);
for(tn=tmin;tn<=tman;tn=tn+5)
{printf("\n tn=%d",tn);
fprintf(fp,"\n温度tn=%d",tn);

f1=e0*g1*g1*lm*lm*(24*(detp/k)*(detp/k))-detp/k+0*e0*t1);
f4=e0*g6*g6*lm*lm*(24*(detp/k)*(detp/k))-detp/k+0*e0*t4);
f6=e0*g1*g1*lm*lm*(24*(detp/k)*(detp/k))-detp/4+0*e0*t6);
f=f1;t=t1;v=v1;g=g1;det=detp/k;
if(f1<f4){t=t4;v=v4;f=f4;g=g6;}
if(f<f6){t=t6;v=v6;f=f6;g=g1;det=detp/4;}
b=e0*g*g*lm*lm/(24.*det*det)-det+a0*e0)*(tn-t);
d=-e0*g1*lm*lm/24.;
/* forintf(fp,"\n a= %d",b);
fprintf(fp," b= %d",d)*/
getslov3lall(a,b,c,d,x,&n);
if(n!=4)for(i=0;i<n;++i)
{fpp=g1*lm*lm/(8*x[i]);
if(x[i]>0){
printf("\n1=%5.1f",lm);
printf(" det=%7.3f",x[i]);
printf(" f=%6.3f",fpp);
fprintf(fp,"\n=档距l=%5.1f",lm);
fprintf(fp," 应力&=%7.3f",lm);
fprintf(fp," 弧垂f=%6.3f",fpp);


else if(n==4)print("\n all Real numbers 解为所有实数");


fclose(fp);
getch0;

搜索更多相关主题的帖子: C语言 卡尔丹 结果 
2007-10-30 21:56
快速回复:求助:关与C语言的计算程序,解卡尔丹(x**3+p*X+q=0)一元三次的实根, ...
数据加载中...
 
   



关于我们 | 广告合作 | 编程中国 | 清除Cookies | TOP | 手机版

编程中国 版权所有,并保留所有权利。
Powered by Discuz, Processed in 0.036288 second(s), 7 queries.
Copyright©2004-2024, BCCN.NET, All Rights Reserved