// [このプログラムの目的]卵の体積と表面積の 短径/長径 依存性のグラフを描くこと,2011年7月21日(木)
// 曲線の方程式:(x*x+y*y)**2=a*(x**3)+(a-b)*x*(y**2), 以下、a=1 として, b の値を 0=< b<=a に変えてグラフを計算する。
// file name: volume_&_surface.c
#include< stdio.h>
#include< math.h>
void main(void)
{
double v,s,ds,a,b,q;
double l1,l2;
double pi;
int i,imax;
double bmax,db;
double x,y,ypast;
double xmax,dx,dydx;
double l12[5001],vv[5001],ss[5001];
FILE *fp;
// 定数設定
pi=3.1415927;
a=1.0;
l1=a;// 長径
// 他のパラメータ設定
bmax=a;// b の最大値
db=bmax/2000;// b のプロット間隔
xmax=a;// x の最大値
dx=xmax/2000;// 卵型の横軸座標 x のプロット間隔
// 計算実行
i=0;
for(b=0;b<=bmax;b=b+db)
{
i++;
// 体積計算
if(b==0)
{
v=pi*a*a*a/6;
}
else
{
v=(pi/2)*((a+b)*(a+b)*(a+b)*a/(6*b)-a*a*a/6-a*a*b/2-((a+b)*(a+b)*(a+b)*(a+b)*(a+b)-(a-b)*(a-b)*(a-b)*(a-b)*(a-b))/(60*b*b));
}
// 表面積計算
if(b==0)
{
s=pi*a*a;
}
else
{
for(x=dx;x<=xmax;x=x+dx)
{
y=sqrt(((a-b-2*x)+sqrt(4*b*x+(a-b)*(a-b)))*x/2);// b=O のときの円の方程式 y=sqrt(a*a/2/2-(x-a/2)*(x-a/2)) を含んでいる.
dydx=(1/(2.0*sqrt(2.0)))*sqrt((a-b-2.0*x+sqrt(4.0*b*x+(a-b)*(a-b)))/x)+(sqrt(x/2.0))*((b/sqrt(4.0*b*x+(a-b)*(a-b))-1)/sqrt(a-b-2.0*x+sqrt(4.0*b*x+(a-b)*(a-b))));
ds=2*pi*y*sqrt(1+dydx*dydx)*dx;
s=s+ds;
}
}
// 短径計算
ypast=-1;
for(x=0;x<=xmax;x=x+dx)
{
y=sqrt(((a-b-2*x)+sqrt(4*b*x+(a-b)*(a-b)))*x/2);// b=O のときの円の方程式 y=sqrt(a*a/2/2-(x-a/2)*(x-a/2)) を含んでいる.
if(y< ypast)
{
q=ypast;
break;
}
ypast=y;
}
l2=2*q;// 短径
// データのまとめ
vv[i]=v;
ss[i]=s;
l12[i]=l2/l1;
printf("i=%d,l12=%f,v=%f,s=%f\n",i,l12[i],v,s);
s=0;
}
imax=i;
// 計算データのテキストファイルへの書き込み
fp=fopen("volume_&_surface.txt","w");
if(fp==NULL)
{
printf("FILE OPEN ERROR\n");
}
else
{
for(i=1;i<=imax;i++)
{
fprintf(fp,"%f,%f,%f\n",l12[i],vv[i],ss[i]);
}
fflush(fp);
fclose(fp);
}
printf("end\n");
}// the end of the programm
戻る