|
0
|
1 #include<stdio.h>
|
|
|
2 #include<stdlib.h>
|
|
|
3 #include<unistd.h>
|
|
|
4 #include<time.h>
|
|
|
5 #include<float.h>
|
|
3
|
6 #include<math.h>
|
|
0
|
7 #include<SDL.h>
|
|
|
8 #include<OpenGL/gl.h>
|
|
|
9 #include<OpenGL/glu.h>
|
|
|
10
|
|
3
|
11 #ifdef DEBUG
|
|
|
12 #define DEBUGlog(f, args...) \
|
|
|
13 fprintf(stderr, "in %s: "f, __FUNCTION__, ## args)
|
|
|
14 #else
|
|
|
15 #define DEBUGlog(f, args...) ;
|
|
|
16 #endif
|
|
0
|
17
|
|
|
18 #define W_HEIGHT 480
|
|
|
19 #define W_WIDTH 640
|
|
|
20
|
|
|
21 /*
|
|
|
22 N body problem
|
|
|
23 */
|
|
|
24
|
|
|
25 /* parameters */
|
|
|
26 static int NUM_BODY = 3;
|
|
1
|
27 static float Gravitation = 100.0f; // ?
|
|
|
28 static float delta = 0.05f; // 0.01 ~ 100 ?
|
|
|
29 static const float eps = 15.0f;
|
|
3
|
30 /*
|
|
|
31 * だいたい、平均距離をn倍にするんだったら、Gravitationoはn^2倍、
|
|
|
32 * epsはn倍
|
|
|
33 */
|
|
0
|
34
|
|
|
35 /* for OpenGL Utility. */
|
|
|
36 GLUquadricObj **sphere;
|
|
2
|
37 float DefDistance = 300.0f;
|
|
0
|
38
|
|
|
39 typedef struct
|
|
|
40 {
|
|
|
41 /* star's parameter. */
|
|
|
42 float weight;
|
|
|
43 float a[3]; /* acceleration */
|
|
|
44 float v[3]; /* velocity */
|
|
|
45 float r[3]; /* position */
|
|
|
46 /* for SDL. */
|
|
|
47 //SDL_Rect rect;
|
|
|
48 } body;
|
|
|
49 body *stars_old, *stars_new;
|
|
|
50
|
|
|
51
|
|
|
52 void start(void);
|
|
|
53 __code initialize(int num);
|
|
|
54 __code randomInit(SDL_Surface *screen, int num);
|
|
|
55 __code starsInit(SDL_Surface *screen, int num);
|
|
|
56 __code loop(int count, SDL_Surface *screen, int num);
|
|
|
57 __code compute(int count, SDL_Surface *screen, int num);
|
|
|
58 __code nextTurn(int count, SDL_Surface *screen, int num);
|
|
3
|
59 __code drawStars(SDL_Surface *screen, int num);
|
|
2
|
60 __code AdjustLooking(float distance, int count, SDL_Surface *screen, int num);
|
|
0
|
61
|
|
|
62 int main(int argc, char **argv)
|
|
|
63 {
|
|
|
64 int ch;
|
|
|
65 while ((ch = getopt(argc, argv, "s:g:")) != -1) {
|
|
|
66 switch (ch) {
|
|
|
67 case 's':
|
|
|
68 delta = atof(optarg);
|
|
|
69 break;
|
|
|
70 case 'g':
|
|
|
71 Gravitation = atof(optarg);
|
|
|
72 break;
|
|
|
73 case '?':
|
|
|
74 default:
|
|
|
75 break;
|
|
|
76 }
|
|
|
77 }
|
|
|
78 start();
|
|
|
79 return 0;
|
|
|
80 }
|
|
|
81 void start()
|
|
|
82 {
|
|
|
83 goto initialize(NUM_BODY);
|
|
|
84 }
|
|
|
85
|
|
|
86 __code finish(int ret, int num)
|
|
|
87 {
|
|
|
88 int i;
|
|
|
89 DEBUGlog("Gravity = %e\n", Gravitation);
|
|
|
90 DEBUGlog("FLT_MAX = %e\n", FLT_MAX);
|
|
|
91 DEBUGlog("FLT_MAX_EXP = %d\n", FLT_MAX_EXP);
|
|
|
92 DEBUGlog("FLT_MIN = %e\n", FLT_MIN);
|
|
|
93 DEBUGlog("FLT_MIN_EXP = %d\n", FLT_MIN_EXP);
|
|
|
94 DEBUGlog("FLT_EPSILON = %e\n", FLT_EPSILON);
|
|
|
95 free(stars_old);
|
|
|
96 free(stars_new);
|
|
|
97 for( i=0; i<NUM_BODY; i++){
|
|
|
98 gluDeleteQuadric(sphere[i]);
|
|
|
99 }
|
|
|
100 free(sphere);
|
|
|
101 exit(ret);
|
|
|
102 }
|
|
|
103
|
|
|
104
|
|
|
105 __code initialize(int num)
|
|
|
106 {
|
|
|
107 SDL_Surface *screen;
|
|
|
108 int i;
|
|
|
109
|
|
|
110 /* malloc. */
|
|
|
111 stars_old = malloc( sizeof(body)*num );
|
|
|
112 stars_new = malloc( sizeof(body)*num );
|
|
|
113 sphere = malloc( sizeof(GLUquadricObj*)*num );
|
|
|
114 if (stars_old==NULL||stars_new==NULL||sphere==NULL){
|
|
|
115 perror("malloc");
|
|
|
116 goto finish(1, num);
|
|
|
117 }
|
|
|
118
|
|
|
119 /* SDL initialization. */
|
|
|
120 if(SDL_Init(SDL_INIT_VIDEO) < 0){ //Could we start SDL_VIDEO?
|
|
|
121 fprintf(stderr,"Couldn't init SDL"); //Nope, output to stderr and quit
|
|
|
122 exit(1);
|
|
|
123 }
|
|
|
124 SDL_GL_SetAttribute( SDL_GL_RED_SIZE, 5);
|
|
|
125 SDL_GL_SetAttribute( SDL_GL_GREEN_SIZE, 5);
|
|
|
126 SDL_GL_SetAttribute( SDL_GL_BLUE_SIZE, 5);
|
|
|
127 SDL_GL_SetAttribute( SDL_GL_DEPTH_SIZE, 16);
|
|
|
128 SDL_GL_SetAttribute( SDL_GL_DOUBLEBUFFER, 1);
|
|
|
129
|
|
|
130 screen = SDL_SetVideoMode(W_WIDTH, W_HEIGHT, 32, SDL_HWSURFACE | SDL_RESIZABLE | SDL_OPENGL ); //Create a 640x480x32 resizable window
|
|
|
131 atexit(SDL_Quit); //Now that we're enabled, make sure we cleanup
|
|
|
132
|
|
|
133 /* initialize OpenGL. */
|
|
|
134 glViewport(0,0,screen->w,screen->h);
|
|
|
135 glMatrixMode(GL_PROJECTION);
|
|
|
136 glLoadIdentity();
|
|
|
137 gluPerspective( 60.0, (float)screen->w/(float)screen->h, 1.0, 1000.0);
|
|
1
|
138 gluLookAt( 300.0,300.0,300.0, 0.0,0.0,0.0, 1.0,0.0,0.0);
|
|
0
|
139 glClearColor(0.0, 0.0, 0.0, 0.0);
|
|
|
140 glMatrixMode(GL_MODELVIEW);
|
|
|
141
|
|
|
142 for( i=0; i<num; i++){
|
|
|
143 sphere[i] = gluNewQuadric();
|
|
|
144 gluQuadricDrawStyle(sphere[i], GLU_FILL);
|
|
|
145 }
|
|
|
146
|
|
|
147 goto starsInit(screen, num);
|
|
|
148 }
|
|
|
149
|
|
1
|
150 __code starsInitRandom(SDL_Surface *screen, int num)
|
|
0
|
151 {
|
|
|
152 int i;
|
|
|
153 srandom(time(NULL));
|
|
|
154 for (i=0; i<num; i++){ // this loop should be split into few code segment..
|
|
|
155 stars_old[i].weight = random()/(RAND_MAX+1.0)*5+5;
|
|
|
156 stars_old[i].v[0] = random()/(RAND_MAX+1.0)*5+5;
|
|
|
157 stars_old[i].v[1] = random()/(RAND_MAX+1.0)*5+5;
|
|
|
158 stars_old[i].v[2] = random()/(RAND_MAX+1.0)*5+5;
|
|
|
159 stars_old[i].r[0] = random()/(RAND_MAX+1.0)*5+5;
|
|
|
160 stars_old[i].r[1] = random()/(RAND_MAX+1.0)*5+5;
|
|
|
161 stars_old[i].r[2] = random()/(RAND_MAX+1.0)*5+5;
|
|
|
162 stars_new[i].weight = stars_old[i].weight;
|
|
|
163 }
|
|
|
164
|
|
|
165 goto loop(0, screen, num);
|
|
|
166 }
|
|
|
167 __code starsInit(SDL_Surface *screen, int num)
|
|
|
168 {
|
|
|
169 int i;
|
|
|
170 /* */
|
|
1
|
171 stars_old[0].weight = 1000;
|
|
0
|
172 stars_old[0].v[0] = 0.0;
|
|
|
173 stars_old[0].v[1] = 0.0;
|
|
|
174 stars_old[0].v[2] = 0.0;
|
|
|
175 stars_old[0].r[0] = 0.0;
|
|
|
176 stars_old[0].r[1] = 0.0;
|
|
|
177 stars_old[0].r[2] = 0.0;
|
|
1
|
178 /* */
|
|
|
179 stars_old[1].weight = 5;
|
|
|
180 stars_old[1].v[0] = 0.1;
|
|
|
181 stars_old[1].v[1] = 5.0;
|
|
0
|
182 stars_old[1].v[2] = 0.0;
|
|
1
|
183 stars_old[1].r[0] = 100.0;
|
|
0
|
184 stars_old[1].r[1] = 0.0;
|
|
|
185 stars_old[1].r[2] = 0.0;
|
|
1
|
186 /* */
|
|
|
187 stars_old[2].weight = 5;
|
|
0
|
188 stars_old[2].v[0] = 0.0;
|
|
1
|
189 stars_old[2].v[1] = -5.0;
|
|
|
190 stars_old[2].v[2] = 0.1;
|
|
|
191 stars_old[2].r[0] = -100.0;
|
|
0
|
192 stars_old[2].r[1] = 0.0;
|
|
|
193 stars_old[2].r[2] = 0.0;
|
|
|
194
|
|
|
195 for( i=0; i<num; i++){
|
|
|
196 stars_new[i].weight = stars_old[i].weight;
|
|
|
197 }
|
|
|
198
|
|
|
199 goto loop(0, screen, num);
|
|
|
200 }
|
|
|
201
|
|
3
|
202
|
|
|
203 /* MAIN LOOP ROUTIN
|
|
|
204 *
|
|
|
205 +------+
|
|
|
206 | loop |<-----------------------+
|
|
|
207 +---+--+ |
|
|
|
208 +-----------+ |
|
|
|
209 |count=num | count<num |
|
|
|
210 +---v------+ +--v-------------+ |
|
|
|
211 | nextTurn | | COMPUTE ROUTIN | |
|
|
|
212 +---+------+ +-------+--------+ |
|
|
|
213 | | |
|
|
|
214 +---v---------+ v |
|
|
|
215 | DRAW ROUTIN +----->+----------+
|
|
|
216 +-------------+ back to loop
|
|
|
217 */
|
|
0
|
218 __code loop(int count, SDL_Surface *screen, int num)
|
|
|
219 {
|
|
|
220 SDL_Event event;
|
|
|
221
|
|
3
|
222 usleep(1000);
|
|
0
|
223 /* check SDL event. */
|
|
|
224 while(SDL_PollEvent(&event)){ //Poll events
|
|
|
225 switch(event.type){ //Check event type
|
|
|
226 case SDL_QUIT: //User hit the X (or equivelent)
|
|
|
227 goto finish(1, num);
|
|
|
228 case SDL_KEYDOWN:
|
|
2
|
229 if (event.key.state==SDL_PRESSED){
|
|
|
230 switch(event.key.keysym.sym){
|
|
|
231 case SDLK_UP:
|
|
|
232 goto AdjustLooking(DefDistance/=1.2, count, screen, num);
|
|
|
233 case SDLK_DOWN:
|
|
|
234 goto AdjustLooking(DefDistance*=1.2, count, screen, num);
|
|
|
235 case SDLK_l:
|
|
|
236 goto AdjustLooking(DefDistance, count, screen, num);
|
|
|
237 case SDLK_q:
|
|
|
238 goto finish(0, num);
|
|
|
239 default:
|
|
|
240 break;
|
|
|
241 }
|
|
|
242 }
|
|
0
|
243 break;
|
|
|
244 default:
|
|
|
245 break;
|
|
|
246 } //Finished with current event
|
|
|
247 } //Done with all events for now
|
|
|
248
|
|
|
249 if ( count<num ){
|
|
|
250 DEBUGlog("count %d, goto commpute().\n", count);
|
|
|
251 goto compute(count, screen, num);
|
|
|
252 }else{
|
|
|
253 DEBUGlog("count %d, goto nextTurn()\n", count);
|
|
|
254 goto nextTurn(count, screen, num);
|
|
|
255 }
|
|
|
256 }
|
|
|
257
|
|
2
|
258 __code AdjustLooking(float distance, int count, SDL_Surface *screen, int num)
|
|
0
|
259 {
|
|
|
260 glMatrixMode(GL_PROJECTION);
|
|
|
261 glLoadIdentity();
|
|
2
|
262 gluPerspective( 60.0, (float)screen->w/(float)screen->h, 1.0, distance*2.5);
|
|
|
263 gluLookAt( distance,distance,distance, 0.0,0.0,0.0, 1.0,0.0,0.0);
|
|
0
|
264 glClearColor(0.0, 0.0, 0.0, 0.0);
|
|
|
265 glMatrixMode(GL_MODELVIEW);
|
|
|
266
|
|
|
267 goto loop(0, screen, num);
|
|
|
268 }
|
|
|
269
|
|
|
270 /*
|
|
3
|
271 * COMPUTING ROUTINs
|
|
|
272 *
|
|
0
|
273 * next x = x+dx
|
|
|
274 * dx = v*dt => dx/dt = v
|
|
|
275 * dv = a*dt => dv/dt = a
|
|
|
276 * a = F/m;
|
|
|
277 * F = F1 + F2 + ... + Fn
|
|
|
278 * Fi = G m1m2/r^2
|
|
3
|
279 |
|
|
|
280 +-----v-----+
|
|
|
281 | compute |
|
|
|
282 +-----+-----+
|
|
|
283 |<----------+
|
|
|
284 +-----v---------+ | i++
|
|
|
285 | computeForce +-+
|
|
|
286 +-----+---------+ i<num
|
|
|
287 | i==num
|
|
|
288 +-----+---------+
|
|
|
289 | com|ute_out |
|
|
|
290 +-----+---------+
|
|
|
291 v to loop
|
|
0
|
292 */
|
|
3
|
293 __code compute_out(int count, SDL_Surface *screen, int num);
|
|
|
294 __code computeForce(int i, int count, SDL_Surface *screen, int num);
|
|
0
|
295 __code compute(int count, SDL_Surface *screen, int num)
|
|
|
296 {
|
|
|
297 /* a is accel this planet receive now. */
|
|
|
298 stars_old[count].a[0]=stars_old[count].a[1]=stars_old[count].a[2]=0.0f;
|
|
|
299 DEBUGlog("count=%d\n", count);
|
|
3
|
300
|
|
|
301 goto computeForce( 0, count, screen, num);
|
|
|
302 }
|
|
|
303 __code computeForce(int i, int count, SDL_Surface *screen, int num)
|
|
|
304 {
|
|
|
305 //float F;
|
|
|
306 float d0, d1, d2, d;
|
|
|
307 float a;
|
|
|
308 //body *o, *m;
|
|
|
309 //o = &stars_old[i]; m = &stars_old[count];
|
|
|
310
|
|
|
311 if (i>=num) goto compute_out(count, screen, num);
|
|
|
312 /* skip compute with itself. */
|
|
|
313 else if (i==count) goto computeForce(i+1, count, screen, num);
|
|
0
|
314
|
|
3
|
315 /* compute distance between two i-th planet and itself. */
|
|
|
316 d0 = stars_old[i].r[0] - stars_old[count].r[0];
|
|
|
317 d1 = stars_old[i].r[1] - stars_old[count].r[1];
|
|
|
318 d2 = stars_old[i].r[2] - stars_old[count].r[2];
|
|
|
319 d = ( d0*d0+d1*d1+d2*d2+eps*eps );
|
|
|
320 /* compute force it receive from i-th planet. */
|
|
|
321 //F = Gravitation * stars_old[i].weight * stars_old[count].weight / d;
|
|
|
322 /* and accel. */
|
|
|
323 //a = F/stars_old[count].weight;
|
|
|
324 a = Gravitation/d * stars_old[i].weight ;
|
|
|
325 stars_old[count].a[0] += a*d0/sqrtf(d);
|
|
|
326 stars_old[count].a[1] += a*d1/sqrtf(d);
|
|
|
327 stars_old[count].a[2] += a*d2/sqrtf(d);
|
|
|
328 DEBUGlog("a=%e, d=%e, d0=%e\n", a, d, d0);
|
|
|
329 DEBUGlog("g*w=%e\n", Gravitation*stars_old[i].weight);
|
|
|
330
|
|
|
331 goto computeForce( i+1, count, screen, num);
|
|
|
332 }
|
|
|
333 __code compute_out(int count, SDL_Surface *screen, int num)
|
|
|
334 {
|
|
0
|
335
|
|
|
336 stars_new[count].v[0] = stars_old[count].v[0]+stars_old[count].a[0]*delta;
|
|
|
337 stars_new[count].v[1] = stars_old[count].v[1]+stars_old[count].a[1]*delta;
|
|
|
338 stars_new[count].v[2] = stars_old[count].v[2]+stars_old[count].a[2]*delta;
|
|
|
339 stars_new[count].r[0] = stars_old[count].r[0]+stars_new[count].v[0]*delta/*+stars_old[count].a[0]*delta*delta/2.0*/;
|
|
|
340 stars_new[count].r[1] = stars_old[count].r[1]+stars_new[count].v[1]*delta/*+stars_old[count].a[1]*delta*delta/2.0*/;
|
|
|
341 stars_new[count].r[2] = stars_old[count].r[2]+stars_new[count].v[2]*delta/*+stars_old[count].a[2]*delta*delta/2.0*/;
|
|
|
342
|
|
|
343 DEBUGlog("x = %e\n", stars_new[count].r[0]);
|
|
|
344 DEBUGlog("y = %e\n", stars_new[count].r[1]);
|
|
|
345 DEBUGlog("z = %e\n", stars_new[count].r[2]);
|
|
|
346 DEBUGlog("vx = %e\n", stars_new[count].v[0]);
|
|
|
347 DEBUGlog("vy = %e\n", stars_new[count].v[1]);
|
|
|
348 DEBUGlog("vz = %e\n", stars_new[count].v[2]);
|
|
|
349 DEBUGlog("a0=%e,a1=%e,a2=%e\n", stars_new[count].a[0], stars_new[count].a[1], stars_new[count].a[2]);
|
|
|
350
|
|
3
|
351 goto loop(1+count, screen, num);
|
|
0
|
352 }
|
|
|
353
|
|
|
354 /* it is executed after all bodies parameter is computed. */
|
|
|
355 __code nextTurn(int count, SDL_Surface *screen, int num)
|
|
|
356 {
|
|
|
357 body *tmp;
|
|
|
358 tmp = stars_old;
|
|
|
359 stars_old = stars_new;
|
|
|
360 stars_new = tmp;
|
|
3
|
361 //SDL_UpdateRect(screen, 0, 0, 0, 0);
|
|
0
|
362
|
|
3
|
363 goto drawStars(screen, num);
|
|
0
|
364 }
|
|
|
365
|
|
3
|
366
|
|
|
367
|
|
|
368 /*
|
|
|
369 * DRAWING ROUTINs
|
|
|
370 *
|
|
|
371 |
|
|
|
372 +----+------+
|
|
|
373 | drawStars |
|
|
|
374 +----+------+
|
|
|
375 |<--------+
|
|
|
376 +----v-----+ | count++
|
|
|
377 | drawStar +---+
|
|
|
378 +----+-----+ count<num
|
|
|
379 | count==num
|
|
|
380 +----v----------+
|
|
|
381 | drawStars_out |
|
|
|
382 +----+----------+
|
|
|
383 | to loop
|
|
|
384 */
|
|
|
385 __code drawStar(int i, SDL_Surface *screen, int num);
|
|
|
386 __code drawStars_out(SDL_Surface *screen, int num);
|
|
|
387
|
|
|
388 __code drawStars(SDL_Surface *screen, int num)
|
|
0
|
389 {
|
|
|
390 DEBUGlog("draw\n");
|
|
|
391 glClear( GL_COLOR_BUFFER_BIT );
|
|
|
392 glColor3d( 1.0, 1.0, 1.0);
|
|
|
393 glBegin(GL_LINES);
|
|
|
394 glVertex3d( 0.0, 0.0, 0.0);
|
|
2
|
395 glVertex3d( 200.0, 0.0, 0.0);
|
|
0
|
396 glEnd();
|
|
|
397 glBegin(GL_LINES);
|
|
|
398 glVertex3d( 0.0, 0.0, 0.0);
|
|
2
|
399 glVertex3d( 0.0, 200.0, 0.0);
|
|
0
|
400 glEnd();
|
|
|
401 glBegin(GL_LINES);
|
|
|
402 glVertex3d( 0.0, 0.0, 0.0);
|
|
2
|
403 glVertex3d( 0.0, 0.0, 200.0);
|
|
0
|
404 glEnd();
|
|
|
405
|
|
3
|
406 goto drawStar(0, screen, num);
|
|
|
407 }
|
|
|
408
|
|
|
409 __code drawStars_out(SDL_Surface *screen, int num)
|
|
|
410 {
|
|
|
411 DEBUGlog("draw axis: screen=%p, num=%d\n", screen, num);
|
|
0
|
412 SDL_GL_SwapBuffers();
|
|
3
|
413 DEBUGlog("go loop\n");
|
|
0
|
414 goto loop(0, screen, num);
|
|
|
415 }
|
|
|
416
|
|
3
|
417 __code drawStar(int i, SDL_Surface *screen, int num)
|
|
|
418 {
|
|
|
419 if (i==num)
|
|
|
420 goto drawStars_out(screen, num);
|
|
|
421
|
|
|
422 DEBUGlog("i=%d, num=%d\n", i, num);
|
|
|
423 glPushMatrix();
|
|
|
424 glTranslatef( stars_new[i].r[0], stars_new[i].r[1], stars_new[i].r[2]);
|
|
|
425 gluSphere( sphere[i], 2.0, 32.0, 32.0 );
|
|
|
426 glPopMatrix();
|
|
|
427
|
|
|
428 DEBUGlog("to myself\n");
|
|
|
429 goto drawStar(i+1, screen, num);
|
|
|
430 }
|
|
|
431
|
