1 /*****************************************************************************
2 * Copyright (c) 2019 FrontISTR Commons
3 * This software is released under the MIT License, see LICENSE.txt
4 *****************************************************************************/
5
6 /**
7 * @brief 並列計算された結果を読込み処理するためのユーティリティ
8 */
9
10 #include "fstr_rmerge_util.h"
11
12 extern FILE* log_fp;
13 extern int nrank;
14 extern int strid;
15 extern int endid;
16 extern int intid;
17
out_log(const char * fmt,...)18 static void out_log(const char* fmt, ...) {
19 va_list arg;
20 va_start(arg, fmt);
21 vfprintf(log_fp, fmt, arg);
22 va_end(arg);
23 }
24
25 /**
26 * @brief 全分散メッシュの読込み
27 */
28
get_dist_fname(char * name_ID,char * fheader,int * fg_single,int * refine,int irank)29 static int get_dist_fname(char* name_ID, char* fheader, int* fg_single,
30 int* refine, int irank) {
31 struct hecmw_ctrl_meshfiles* files;
32
33 files = HECMW_ctrl_get_meshfiles_header_sub(name_ID, nrank, irank);
34 if (!files) return -1;
35
36 if (files->n_mesh == 1) {
37 strcpy(fheader, files->meshfiles[0].filename);
38 if (files->meshfiles[0].type == HECMW_CTRL_FTYPE_HECMW_DIST) {
39 *fg_single = 0;
40 } else {
41 *fg_single = 1;
42 }
43 *refine = files->meshfiles[0].refine;
44 out_log("refine number is %d\n", *refine);
45 } else {
46 HECMW_ctrl_free_meshfiles(files);
47 return -1;
48 }
49
50 HECMW_ctrl_free_meshfiles(files);
51 return 0;
52 }
53
get_area_n(char * fheader)54 static int get_area_n(char* fheader) {
55 FILE* fp;
56 char buff[HECMW_FILENAME_LEN + 1];
57 int area = 0;
58
59 while (1) {
60 sprintf(buff, "%s.%d", fheader, area);
61 out_log("try open : %s ... ", buff);
62 fp = fopen(buff, "r");
63 if (!fp) {
64 out_log("fail\n");
65 out_log("area number is %d\n", area);
66 return area;
67 } else {
68 out_log("success\n");
69 fclose(fp);
70 }
71 area++;
72 }
73 }
74
75 /**
76 * @brief 全分散メッシュの読込み
77 */
78
fstr_get_all_local_mesh(char * name_ID,int * area_number,int * refine)79 struct hecmwST_local_mesh** fstr_get_all_local_mesh(char* name_ID,
80 int* area_number,
81 int* refine) {
82 int i;
83 int area_n;
84 int fg_single;
85 char fheader[HECMW_HEADER_LEN + 1];
86 char fname[HECMW_FILENAME_LEN + 1];
87 struct hecmwST_local_mesh** mesh;
88
89 if (get_dist_fname(name_ID, fheader, &fg_single, refine, 0)) return NULL;
90
91 if (fg_single) {
92 out_log("mesh file type is NOT HECMW_DIST.\n");
93 area_n = 1;
94 out_log("area number is %d\n", area_n);
95 mesh = HECMW_malloc(area_n * sizeof(struct hecmwST_local_mesh*));
96 mesh[0] = HECMW_get_mesh(name_ID);
97 if (!mesh[0]) return NULL;
98 } else {
99 out_log("mesh file type is HECMW_DIST.\n");
100 if (nrank == 0) {
101 area_n = get_area_n(fheader);
102 } else {
103 area_n = nrank;
104 }
105 if (area_n == 0) return NULL;
106 mesh = HECMW_malloc(area_n * sizeof(struct hecmwST_local_mesh*));
107 for (i = 0; i < area_n; i++) {
108 if (nrank == 0) {
109 sprintf(fname, "%s.%d", fheader, i);
110 } else {
111 get_dist_fname(name_ID, fheader, &fg_single, refine, i);
112 sprintf(fname, "%s.%d", fheader, i);
113 }
114 out_log("loading dist mesh from %s\n", fname);
115 mesh[i] = HECMW_get_dist_mesh(fname);
116 if (!mesh[i]) return NULL;
117 }
118 }
119 *area_number = area_n;
120
121 return mesh;
122 }
123
124 /**
125 * @brief メッシュの削除
126 */
127
fstr_free_mesh(struct hecmwST_local_mesh ** mesh,int area_n)128 void fstr_free_mesh(struct hecmwST_local_mesh** mesh, int area_n) {
129 int i;
130
131 if (!mesh) return;
132
133 for (i = 0; i < area_n; i++) {
134 HECMW_dist_free(mesh[i]);
135 }
136 HECMW_free(mesh);
137 return;
138 }
139
140 /**
141 * @brief ステップ数を調べる(ファイルの存在を調べる)
142 */
143
fstr_get_step_n(char * name_ID)144 int fstr_get_step_n(char* name_ID) {
145 FILE* fp;
146 int step, fg_text;
147 char* fheader;
148 char fname[HECMW_FILENAME_LEN + 1];
149
150 if (endid > -1) return endid;
151
152 if (nrank == 0) {
153 if ((fheader = HECMW_ctrl_get_result_fileheader(name_ID, 1, &fg_text)) ==
154 NULL)
155 return 0;
156 } else {
157 if ((fheader = HECMW_ctrl_get_result_fileheader_sub(name_ID, 1, nrank, 0,
158 &fg_text)) == NULL)
159 return 0;
160 }
161
162 step = 1;
163 while (1) {
164 sprintf(fname, "%s.0.%d", fheader, step);
165 out_log("try open : %s ... ", fname);
166 fp = fopen(fname, "r");
167 if (!fp) {
168 out_log("fail\n");
169 out_log("step number is %d\n", step - 1);
170 return step - 1;
171 } else {
172 out_log("success\n");
173 fclose(fp);
174 }
175 step++;
176 }
177 }
178
179 /**
180 * @brief ステップの全領域データの読み込み
181 */
182
fstr_get_all_result(char * name_ID,int step,int area_n,int refine)183 fstr_res_info** fstr_get_all_result(char* name_ID, int step, int area_n,
184 int refine) {
185 char* fheader;
186 char fname[HECMW_FILENAME_LEN + 1];
187 fstr_res_info** res;
188 struct hecmwST_result_data* data;
189 int i, j, k, count, num, nnode, nelem, flag, fg_text;
190 int *node_gid, *elem_gid;
191 int refine_nnode = 0;
192
193 if (nrank == 0) {
194 if ((fheader = HECMW_ctrl_get_result_fileheader(name_ID, step,
195 &fg_text)) == NULL)
196 return 0;
197 }
198
199 res = HECMW_malloc(area_n * sizeof(fstr_res_info*));
200 if (!res) return NULL;
201
202 for (i = 0; i < area_n; i++) {
203 res[i] = HECMW_malloc(sizeof(fstr_res_info));
204 if (!res[i]) return NULL;
205
206 if (nrank != 0) {
207 if ((fheader = HECMW_ctrl_get_result_fileheader_sub(
208 name_ID, step, nrank, i, &fg_text)) == NULL)
209 return 0;
210 }
211 sprintf(fname, "%s.%d.%d", fheader, i, step);
212 data = HECMW_result_read_by_fname(fname);
213 if (!data) return NULL;
214 nnode = HECMW_result_get_nnode();
215 node_gid = HECMW_malloc(nnode * sizeof(int));
216 if (!node_gid) return NULL;
217 HECMW_result_get_nodeID(node_gid);
218 nelem = HECMW_result_get_nelem();
219 elem_gid = HECMW_malloc(nelem * sizeof(int));
220 if (!elem_gid) return NULL;
221 if (data->ne_component) HECMW_result_get_elemID(elem_gid);
222
223 if (refine) {
224 res[i]->result = HECMW_malloc(sizeof(struct hecmwST_result_data));
225 if (!res[i]->result) return NULL;
226
227 res[i]->result->ng_component = data->ng_component;
228 res[i]->result->ng_dof = HECMW_malloc(data->ng_component * sizeof(int));
229 if (!res[i]->result->ng_dof) return NULL;
230 res[i]->result->global_label =
231 HECMW_malloc(data->ng_component * sizeof(char*));
232 if (!res[i]->result->global_label) return NULL;
233 for (j = 0; j < data->ng_component; j++) {
234 res[i]->result->ng_dof[j] = data->ng_dof[j];
235 res[i]->result->global_label[j] = HECMW_strdup(data->global_label[j]);
236 }
237
238 res[i]->result->nn_component = data->nn_component;
239 res[i]->result->nn_dof = HECMW_malloc(data->nn_component * sizeof(int));
240 if (!res[i]->result->nn_dof) return NULL;
241 res[i]->result->node_label =
242 HECMW_malloc(data->nn_component * sizeof(char*));
243 if (!res[i]->result->node_label) return NULL;
244 num = 0;
245 for (j = 0; j < data->nn_component; j++) {
246 res[i]->result->nn_dof[j] = data->nn_dof[j];
247 res[i]->result->node_label[j] = HECMW_strdup(data->node_label[j]);
248 num += data->nn_dof[j];
249 }
250
251 count = 1;
252 flag = 0;
253 for (j = 1; j < nnode; j++) {
254 if (flag == refine) break;
255 count++;
256 if (node_gid[j] > 0 && node_gid[j - 1] < 0) flag++;
257 if (node_gid[j] < 0 && node_gid[j] > node_gid[j - 1]) flag++;
258 }
259 count--;
260 out_log("\narea:%d -- refined_nn_internal:%d", i, count);
261 refine_nnode += count;
262
263 count = 0;
264 for (j = 0; j < nnode; j++)
265 if (node_gid[j] > 0) count++;
266 res[i]->nnode_gid = count;
267 res[i]->result->node_val_item =
268 HECMW_malloc(num * count * sizeof(double));
269 if (!res[i]->result->node_val_item) return NULL;
270 count = 0;
271 for (j = 0; j < nnode; j++) {
272 if (node_gid[j] > 0) {
273 for (k = 0; k < num; k++) {
274 res[i]->result->node_val_item[count++] =
275 data->node_val_item[num * j + k];
276 }
277 }
278 }
279
280 res[i]->result->ne_component = data->ne_component;
281 res[i]->result->ne_dof = HECMW_malloc(data->ne_component * sizeof(int));
282 if (!res[i]->result->ne_dof) return NULL;
283 res[i]->result->elem_label =
284 HECMW_malloc(data->ne_component * sizeof(char*));
285 if (!res[i]->result->elem_label) return NULL;
286 num = 0;
287 for (j = 0; j < data->ne_component; j++) {
288 res[i]->result->ne_dof[j] = data->ne_dof[j];
289 res[i]->result->elem_label[j] = HECMW_strdup(data->elem_label[j]);
290 num += data->ne_dof[j];
291 }
292
293 count = 0;
294 for (j = 0; j < nelem; j++)
295 if (elem_gid[j] > 0) count++;
296 out_log("\narea:%d -- ne_original from result:%d", i, count);
297 res[i]->nelem_gid = count;
298 res[i]->result->elem_val_item =
299 HECMW_malloc(num * count * sizeof(double));
300 if (!res[i]->result->elem_val_item) return NULL;
301 count = 0;
302 for (j = 0; j < nelem; j++) {
303 if (elem_gid[j] > 0) {
304 for (k = 0; k < num; k++) {
305 res[i]->result->elem_val_item[count++] =
306 data->elem_val_item[num * j + k];
307 }
308 }
309 }
310
311 HECMW_result_free(data);
312 HECMW_result_free_nodeID();
313 HECMW_result_free_elemID();
314
315 res[i]->node_gid = HECMW_malloc(res[i]->nnode_gid * sizeof(int));
316 if (!res[i]->node_gid) return NULL;
317 count = 0;
318 for (j = 0; j < nnode; j++) {
319 if (node_gid[j] > 0) res[i]->node_gid[count++] = node_gid[j];
320 }
321 free(node_gid);
322
323 res[i]->elem_gid = HECMW_malloc(res[i]->nelem_gid * sizeof(int));
324 if (!res[i]->elem_gid) return NULL;
325 count = 0;
326 for (j = 0; j < nelem; j++) {
327 if (elem_gid[j] > 0) res[i]->elem_gid[count++] = elem_gid[j];
328 }
329 free(elem_gid);
330 } else {
331 res[i]->result = data;
332 res[i]->nnode_gid = nnode;
333 res[i]->node_gid = node_gid;
334 res[i]->nelem_gid = nelem;
335 res[i]->elem_gid = elem_gid;
336 }
337 }
338
339 if (refine) out_log("\ntotal refined_nn_internal:%d\n", refine_nnode);
340
341 return res;
342 }
343
344 /**
345 * @brief ステップの全領域データの結合
346 */
347
fstr_all_result(fstr_glt * glt,fstr_res_info ** res,int refine)348 struct hecmwST_result_data* fstr_all_result(fstr_glt* glt, fstr_res_info** res,
349 int refine) {
350 int i, j, k, l_id, area;
351 int gitem, nitem, eitem, count, irec;
352 struct hecmwST_result_data* data;
353
354 gitem = 0;
355 for (i = 0; i < res[0]->result->ng_component; i++)
356 gitem += res[0]->result->ng_dof[i];
357 nitem = 0;
358 for (i = 0; i < res[0]->result->nn_component; i++)
359 nitem += res[0]->result->nn_dof[i];
360 eitem = 0;
361 for (i = 0; i < res[0]->result->ne_component; i++)
362 eitem += res[0]->result->ne_dof[i];
363
364 data = HECMW_malloc(sizeof(struct hecmwST_result_data));
365 data->ng_dof = HECMW_malloc(res[0]->result->ng_component * sizeof(int));
366 data->global_label = HECMW_malloc(res[0]->result->ng_component * sizeof(char*));
367 data->global_val_item = HECMW_malloc(gitem * sizeof(double));
368 data->nn_dof = HECMW_malloc(res[0]->result->nn_component * sizeof(int));
369 data->node_label = HECMW_malloc(res[0]->result->nn_component * sizeof(char*));
370 data->node_val_item = HECMW_malloc(nitem * glt->node_n * sizeof(double));
371 data->ne_dof = HECMW_malloc(res[0]->result->ne_component * sizeof(int));
372 data->elem_label = HECMW_malloc(res[0]->result->ne_component * sizeof(char*));
373 data->elem_val_item = HECMW_malloc(eitem * glt->elem_n * sizeof(double));
374
375 data->ng_component = res[0]->result->ng_component;
376 for (i = 0; i < res[0]->result->ng_component; i++) {
377 data->ng_dof[i] = res[0]->result->ng_dof[i];
378 data->global_label[i] = HECMW_strdup(res[0]->result->global_label[i]);
379 }
380 for (i = 0; i < gitem; i++) {
381 data->global_val_item[i] = res[0]->result->global_val_item[i];
382 }
383 count = 0;
384 for (i = 0; i < glt->node_n; i++) {
385 l_id = glt->nrec[i].local;
386 area = glt->nrec[i].area;
387 irec = nitem * l_id;
388 for (j = 0; j < res[0]->result->nn_component; j++) {
389 for (k = 0; k < res[0]->result->nn_dof[j]; k++) {
390 data->node_val_item[count++] = res[area]->result->node_val_item[irec++];
391 }
392 }
393 }
394 data->nn_component = res[0]->result->nn_component;
395 for (i = 0; i < res[0]->result->nn_component; i++) {
396 data->nn_dof[i] = res[0]->result->nn_dof[i];
397 data->node_label[i] = HECMW_strdup(res[0]->result->node_label[i]);
398 }
399
400 count = 0;
401 for (i = 0; i < glt->elem_n; i++) {
402 l_id = glt->erec[i].local;
403 area = glt->erec[i].area;
404 irec = eitem * l_id;
405 for (j = 0; j < res[0]->result->ne_component; j++) {
406 for (k = 0; k < res[0]->result->ne_dof[j]; k++) {
407 if (refine) {
408 data->elem_val_item[count++] = 0.0;
409 } else {
410 data->elem_val_item[count++] =
411 res[area]->result->elem_val_item[irec++];
412 }
413 }
414 }
415 }
416 data->ne_component = res[0]->result->ne_component;
417 for (i = 0; i < res[0]->result->ne_component; i++) {
418 data->ne_dof[i] = res[0]->result->ne_dof[i];
419 data->elem_label[i] = HECMW_strdup(res[0]->result->elem_label[i]);
420 }
421
422 return data;
423 }
424
425 /**
426 * @brief fstr_res_info の削除
427 */
428
fstr_free_result(fstr_res_info ** res,int area_n)429 void fstr_free_result(fstr_res_info** res, int area_n) {
430 int i;
431
432 if (!res) return;
433
434 for (i = 0; i < area_n; i++) {
435 HECMW_result_free(res[i]->result);
436 HECMW_free(res[i]->node_gid);
437 HECMW_free(res[i]->elem_gid);
438 HECMW_free(res[i]);
439 }
440 HECMW_free(res);
441 return;
442 }
443
444 /**
445 * @brief グローバルとローカル、所属領域のテーブル fstr_glt の作成
446 */
447
cmp_global_glt(const fstr_gl_rec * g1,const fstr_gl_rec * g2)448 static int cmp_global_glt(const fstr_gl_rec* g1, const fstr_gl_rec* g2) {
449 return (g1->global - g2->global);
450 }
451
452 typedef int (*cmp_func)(const void*, const void*);
453
fstr_create_glt(struct hecmwST_local_mesh ** mesh,int area_n)454 fstr_glt* fstr_create_glt(struct hecmwST_local_mesh** mesh, int area_n) {
455 int i, j, k, eid, count;
456 int all_n, all_e;
457 int area_e;
458 fstr_gl_rec* nrec;
459 fstr_gl_rec* erec;
460 fstr_glt* glt;
461 int* area_etype_list;
462
463 all_n = 0;
464 for (i = 0; i < area_n; i++) {
465 out_log("area:%d -- nn_internal:%d\n", i, mesh[i]->nn_internal);
466 all_n += mesh[i]->nn_internal;
467 }
468 out_log("total nn_internal:%d\n", all_n);
469
470 nrec = HECMW_malloc(sizeof(fstr_gl_rec) * all_n);
471 if (!nrec) return NULL;
472
473 count = 0;
474 for (i = 0; i < area_n; i++) {
475 for (j = 0; j < mesh[i]->nn_internal; j++) {
476 nrec[count].global = mesh[i]->global_node_ID[j];
477 nrec[count].local = j;
478 nrec[count].area = i;
479 count++;
480 }
481 }
482
483 qsort(nrec, all_n, sizeof(fstr_gl_rec), (cmp_func)cmp_global_glt);
484
485 all_e = 0;
486 for (i = 0; i < area_n; i++) {
487 area_e = 0;
488 area_etype_list = HECMW_malloc(sizeof(int) * mesh[i]->n_elem);
489 for (j = 0; j < mesh[i]->n_elem_type; j++) {
490 for (k = mesh[i]->elem_type_index[j]; k < mesh[i]->elem_type_index[j+1]; k++) {
491 area_etype_list[k] = mesh[i]->elem_type_item[j];
492 }
493 }
494
495 for (j = 0; j < mesh[i]->ne_internal; j++) {
496 eid = mesh[i]->elem_internal_list[j] - 1;
497 if ( HECMW_is_etype_patch(area_etype_list[eid]) ) continue;
498 area_e++;
499 }
500
501 HECMW_free(area_etype_list);
502 all_e += area_e;
503 out_log("area:%d -- ne_internal:%d\n", i, area_e);
504 }
505 out_log("total ne_internal:%d\n", all_e);
506
507 erec = HECMW_malloc(sizeof(fstr_gl_rec) * all_e);
508 if (!erec) return NULL;
509
510 count = 0;
511 for (i = 0; i < area_n; i++) {
512
513 area_etype_list = HECMW_malloc(sizeof(int) * mesh[i]->n_elem);
514 for (j = 0; j < mesh[i]->n_elem_type; j++) {
515 for (k = mesh[i]->elem_type_index[j]; k < mesh[i]->elem_type_index[j+1]; k++) {
516 area_etype_list[k] = mesh[i]->elem_type_item[j];
517 }
518 }
519
520 for (j = 0; j < mesh[i]->ne_internal; j++) {
521 eid = mesh[i]->elem_internal_list[j] - 1;
522 if ( HECMW_is_etype_patch(area_etype_list[eid]) ) continue;
523 if ( HECMW_is_etype_link(area_etype_list[eid]) ) continue;
524
525 erec[count].global = mesh[i]->global_elem_ID[eid];
526 erec[count].local = eid;
527 erec[count].area = i;
528 count++;
529 }
530
531 HECMW_free(area_etype_list);
532 }
533
534 qsort(erec, all_e, sizeof(fstr_gl_rec), (cmp_func)cmp_global_glt);
535
536 glt = HECMW_malloc(sizeof(fstr_glt));
537 if (!glt) return NULL;
538 glt->nrec = nrec;
539 glt->erec = erec;
540 glt->node_n = all_n;
541 glt->elem_n = all_e;
542
543 return glt;
544 }
545
546 /**
547 * @brief fstr_glt の削除
548 */
549
fstr_free_glt(fstr_glt * glt)550 void fstr_free_glt(fstr_glt* glt) {
551 if (!glt) return;
552
553 HECMW_free(glt->nrec);
554 HECMW_free(glt->erec);
555 HECMW_free(glt);
556 return;
557 }
558
559 /**
560 * @brief 単一領域メッシュの作成
561 */
562
fstr_create_glmesh(fstr_glt * glt)563 struct hecmwST_local_mesh* fstr_create_glmesh(fstr_glt* glt) {
564 struct hecmwST_local_mesh* mesh;
565 int i;
566
567 mesh = HECMW_calloc(1, sizeof(struct hecmwST_local_mesh));
568 mesh->global_node_ID = HECMW_malloc(glt->node_n * sizeof(int));
569 mesh->global_elem_ID = HECMW_malloc(glt->elem_n * sizeof(int));
570
571 for (i = 0; i < glt->node_n; i++) {
572 mesh->global_node_ID[i] = glt->nrec[i].global;
573 }
574 mesh->n_node = glt->node_n;
575
576 for (i = 0; i < glt->elem_n; i++) {
577 mesh->global_elem_ID[i] = glt->erec[i].global;
578 }
579 mesh->n_elem = glt->elem_n;
580
581 return mesh;
582 }
583
584 /**
585 * @brief 単一領域メッシュの削除
586 */
587
fstr_free_glmesh(struct hecmwST_local_mesh * mesh)588 void fstr_free_glmesh(struct hecmwST_local_mesh* mesh) {
589 if (!mesh) return;
590
591 HECMW_free(mesh->global_node_ID);
592 HECMW_free(mesh->global_elem_ID);
593 HECMW_free(mesh);
594 return;
595 }
596