・3次元空間に引力、斥力を他者に及ぼす点を配置し、その中を質点（持っている力ー１）が動き回りどのような動きをするか確認するプログラム
困っていること
・動き回れる質点を現在の1個から複数にしたいが計算結果がすべて０になりきちんとした結果が得られないこと
・今は無限の3次元空間となっているが有限の3次元空間にしたいができないこと （したい大きさ16*16*16）
・現在は手入力で質点等のデータを入れているがこれを疑似乱数で自動入力方式にしたいこと
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <float.h>
#include <time.h>
#define DT 0.001
#define TIMELIMIT 20.0
#define R 0.1
#define BUFSIZE 256
#define N 256
#define M 3
#define LIMIT 2
#define ULONG_MAX 50
struct coordinate {
double x;/*x座標*/
double y;/*y座標*/
double z;/*z座標*/
};
struct charge {
struct coordinate qxy;/*電荷の位置*/
double q;/*電荷*/
struct coordinate v,x;
};
int getdouble(double *x);/*実数の読み込み*/
int inputq(struct charge qi[]);/*電荷入力*/
int inputn(struct charge qp[]);
int main()
{
struct coordinate v,x;/*質点の速度と位置*/
struct charge qi[N];
struct charge qp[M];
int nofq;
int nofn;
double t=0;
double h=DT;
double rx,ry,rz,r,rmin;
int i;
clock_t start,end;/*処理時間計測関数*/
//入力終了
nofn=inputn(qp);
nofq=inputq(qi);
printf("%f\t%f\t%f\t%f\t%f\t%f\t%f\n",t,x.x,x.y,x.z,v.x,v.y,v.z);
start = clock();
printf("開始時間：%f\n",start);
while(t<=TIMELIMIT){
t+=h;
rmin=DBL_MAX;
for(i=0;i<nofq;++i){
rx=qi[i].qxy.x-x.x;//位置x座標計算
ry=qi[i].qxy.y-x.y;//位置y座標計算
rz=qi[i].qxy.z-x.z;//位置z座標計算
r=sqrt(rx*rx+ry*ry+rz*rz);
if(r<rmin) rmin=r;
qp[i].v.x+=(rx/r/r/r*qi[i].q)*h;//速度x座標計算
qp[i].v.y+=(ry/r/r/r*qi[i].q)*h;//速度y座標計算
qp[i].v.z+=(rz/r/r/r*qi[i].q)*h;//速度z座標計算
}
end = clock();
printf("処理時間=%f\n",(double)(end-start)/CLOCKS_PER_SEC);/*処理時間出力*/
x.x+=v.x*h;
x.y+=v.y*h;
x.z+=v.z*h;
printf("刻み時間=%f\t\n質点の位置(%f\t,%f\t,%f)\n質点の速度(%f\t,%f\t,%f)\n",t,x.x,x.y,x.z,v.x,v.y,v.z);
if(rmin<R) break;
}
return 0;
}
int inputq(struct charge qi[])
{
int i;
for(i=0;i<N;++i){
fprintf(stderr,"電荷%d\n",i);
fprintf(stderr,"電荷の配置箇所qxy.x\n");
if(getdouble(&qi[i].qxy.x)==EOF)break;
fprintf(stderr,"電荷の配置箇所qxy.y\n");
if(getdouble(&qi[i].qxy.y)==EOF)break;
fprintf(stderr,"電荷の配置箇所qxy.z\n");
if(getdouble(&qi[i].qxy.z)==EOF)break;
fprintf(stderr,"電荷の値q\n");
if(getdouble(&qi[i].q)==EOF)break;
}
return i;
}
int inputn(struct charge qp[])
{
int i;
for(i=0;i<M;i++){
//入力開始
fprintf(stderr,"初速度%dv0x\n",i);
if(getdouble(&qp[i].v.x)==EOF)exit(1);
fprintf(stderr,"初速度%dv0y\n",i);
if(getdouble(&qp[i].v.y)==EOF)exit(1);
fprintf(stderr,"初速度%dv0z\n",i);
if(getdouble(&qp[i].v.z)==EOF)exit(1);
fprintf(stderr,"初期位置%dx.x\n",i);
if(getdouble(&qp[i].x.x)==EOF)exit(1);
fprintf(stderr,"初期位置%dx.y\n",i);
if(getdouble(&qp[i].x.y)==EOF)exit(1);
fprintf(stderr,"初期位置%dx.z\n",i);
if(getdouble(&qp[i].x.z)==EOF)exit(1);
}
return i;
}
int getdouble(double *x)
{
char linebuf[BUFSIZE];
int result=0;
if(fgets(linebuf,BUFSIZE,stdin)!=NULL){
if(sscanf(linebuf,"%lf",x)<=0)
result=EOF;
}else{
result=EOF;
}
return result ;
}
メニュー
カレンダー
Clock
最新メンバー
投稿数
統計データ
Link to us
カウンター