FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_dist_refine.c
Go to the documentation of this file.
1/*****************************************************************************
2 * Copyright (c) 2019 FrontISTR Commons
3 * This software is released under the MIT License, see LICENSE.txt
4 *****************************************************************************/
5
6#include <stdio.h>
7#include <stdlib.h>
8#include <errno.h>
9#include <string.h>
10#include "hecmw_struct.h"
11#include "hecmw_util.h"
12#include "hecmw_common_define.h"
13#include "hecmw_etype.h"
14#include "hecmw_bit_array.h"
15#include "hecmw_set_int.h"
16#include "hecmw_varray_int.h"
17#include "hecmw_dist_alloc.h"
18#include "hecmw_dist_free.h"
19#include "hecmw_dist_refine.h"
20
21#ifdef HECMW_WITH_REFINER
22
23#include "rcapRefiner.h"
24
25/*============================================================================*/
26/* */
27/* convert elem_node_item */
28/* */
29/*============================================================================*/
30static const int *get_enode_h2r(int etype, int *ierror) {
31 /* tetra, prism, hexa: no need to convert */
32 static const int enode_pyr_h2r[13] = {4, 0, 1, 2, 3, 9, 10,
33 11, 12, 5, 6, 7, 8};
34 static const int enode_rod2_h2r[3] = {0, 2, 1};
35 const int *enode_h2r;
36
37 *ierror = 0;
38
39 if (HECMW_is_etype_link(etype)) return NULL;
40 if (HECMW_is_etype_truss(etype)) return NULL;
41
42 switch (etype) {
44 /* case HECMW_ETYPE_ROD31: */
69 enode_h2r = NULL;
70 break;
73 enode_h2r = enode_pyr_h2r;
74 break;
77 enode_h2r = enode_rod2_h2r;
78 break;
79 default:
81 "Element type %d not supported for rerinement.\n", etype);
82 *ierror = 1;
83 enode_h2r = NULL;
84 break;
85 }
86 return enode_h2r;
87}
88
89static const int *get_enode_r2h(int etype, int *ierror) {
90 /* tetra, prism, hexa: no need to convert */
91 static const int enode_pyr_r2h[13] = {1, 2, 3, 4, 0, 9, 10,
92 11, 12, 5, 6, 7, 8};
93 static const int enode_rod2_r2h[3] = {0, 2, 1};
94 const int *enode_r2h;
95
96 *ierror = 0;
97
98 if (HECMW_is_etype_link(etype)) return NULL;
99 if (HECMW_is_etype_truss(etype)) return NULL;
100
101 switch (etype) {
102 case HECMW_ETYPE_ROD1:
104 case HECMW_ETYPE_TRI1:
105 case HECMW_ETYPE_TRI2:
106 case HECMW_ETYPE_QUA1:
107 case HECMW_ETYPE_QUA2:
108 case HECMW_ETYPE_TET1:
110 case HECMW_ETYPE_TET2:
111 case HECMW_ETYPE_PRI1:
112 case HECMW_ETYPE_PRI2:
113 case HECMW_ETYPE_HEX1:
115 case HECMW_ETYPE_HEX2:
116 case HECMW_ETYPE_BEM1:
117 case HECMW_ETYPE_BEM3:
118 case HECMW_ETYPE_SHT1:
119 case HECMW_ETYPE_SHT2:
120 case HECMW_ETYPE_SHT6:
121 case HECMW_ETYPE_SHQ1:
122 case HECMW_ETYPE_SHQ2:
123 case HECMW_ETYPE_SHQ8:
124 case HECMW_ETYPE_PTT1:
125 case HECMW_ETYPE_PTT2:
126 case HECMW_ETYPE_PTQ1:
127 case HECMW_ETYPE_PTQ2:
128 enode_r2h = NULL;
129 break;
130 case HECMW_ETYPE_PYR1:
131 case HECMW_ETYPE_PYR2:
132 enode_r2h = enode_pyr_r2h;
133 break;
134 case HECMW_ETYPE_ROD2:
135 case HECMW_ETYPE_BEM2:
136 enode_r2h = enode_rod2_r2h;
137 break;
138 default:
140 "Element type %d not supported for refinement.\n", etype);
141 *ierror = 1;
142 enode_r2h = NULL;
143 }
144 return enode_r2h;
145}
146
147static int elem_node_item_hecmw2rcap(struct hecmwST_local_mesh *mesh) {
148 int i, j;
149 int etype, ierror;
150 const int *enode_h2r;
151 int *enode;
152 int nn, istart, iend;
153 int temp_enode[HECMW_MAX_NODE_MAX];
154
155 for (i = 0; i < mesh->n_elem_type; i++) {
156 etype = mesh->elem_type_item[i];
157 enode_h2r = get_enode_h2r(etype, &ierror);
158 if (ierror) return HECMW_ERROR;
159 if (enode_h2r == NULL) continue;
160 nn = HECMW_get_max_node(etype);
161 istart = mesh->elem_type_index[i];
162 iend = mesh->elem_type_index[i + 1];
163 for (j = istart; j < iend; j++) {
165 for (j = 0; j < nn; j++) temp_enode[j] = enode[enode_h2r[j]];
166 for (j = 0; j < nn; j++) enode[j] = temp_enode[j];
167 }
168 }
169 return HECMW_SUCCESS;
170}
171
172static int elem_node_item_rcap2hecmw(struct hecmwST_local_mesh *mesh) {
173 int i, j;
174 int etype, ierror;
175 const int *enode_r2h;
176 int *enode;
177 int nn, istart, iend;
178 int temp_enode[HECMW_MAX_NODE_MAX];
179
180 for (i = 0; i < mesh->n_elem_type; i++) {
181 etype = mesh->elem_type_item[i];
182 enode_r2h = get_enode_r2h(etype, &ierror);
183 if (ierror) return HECMW_ERROR;
184 if (enode_r2h == NULL) continue;
185 nn = HECMW_get_max_node(etype);
186 istart = mesh->elem_type_index[i];
187 iend = mesh->elem_type_index[i + 1];
188 for (j = istart; j < iend; j++) {
190 for (j = 0; j < nn; j++) temp_enode[j] = enode[enode_r2h[j]];
191 for (j = 0; j < nn; j++) enode[j] = temp_enode[j];
192 }
193 }
194 return HECMW_SUCCESS;
195}
196
197/*============================================================================*/
198/* */
199/* convert surface IDs */
200/* */
201/*============================================================================*/
202static const int *get_sid_h2r(int etype, int *ierror) {
203 /* HECMW to RCAP */
204 /*
205 static const int sid_tri_h2r[4] = {-1, 0, 1, 2};
206 static const int sid_qua_h2r[5] = {-1, 3, 1, 0, 2};
207 static const int sid_tet_h2r[5] = {-1, 0, 1, 2, 3};
208 static const int sid_pri_h2r[6] = {-1, 3, 4, 2, 0, 1};
209 static const int sid_hex_h2r[7] = {-1, 5, 3, 2, 4, 0, 1};
210 static const int sid_pyr_h2r[6] = {-1, 3, 1, 0, 2, 4};
211 */
212 /* FSTR to RCAP*/
213 static const int sid_tri_h2r[4] = {-1, 2, 0, 1};
214 static const int sid_qua_h2r[5] = {-1, 0, 1, 2, 3};
215 static const int sid_tet_h2r[5] = {-1, 3, 2, 0, 1};
216 static const int sid_pri_h2r[6] = {-1, 0, 1, 2, 3, 4};
217 static const int sid_hex_h2r[7] = {-1, 0, 1, 2, 3, 4, 5};
218 static const int sid_pyr_h2r[6] = {-1, 3, 1,
219 0, 2, 4}; /* TEMPORARY: same as HECMW */
220
221 const int *sid_h2r;
222 *ierror = 0;
223 switch (etype) {
224 case HECMW_ETYPE_TRI1:
225 case HECMW_ETYPE_TRI2:
226 sid_h2r = sid_tri_h2r;
227 break;
228 case HECMW_ETYPE_QUA1:
229 case HECMW_ETYPE_QUA2:
230 sid_h2r = sid_qua_h2r;
231 break;
232 case HECMW_ETYPE_TET1:
234 case HECMW_ETYPE_TET2:
235 sid_h2r = sid_tet_h2r;
236 break;
237 case HECMW_ETYPE_PRI1:
238 case HECMW_ETYPE_PRI2:
239 sid_h2r = sid_pri_h2r;
240 break;
241 case HECMW_ETYPE_HEX1:
243 case HECMW_ETYPE_HEX2:
244 sid_h2r = sid_hex_h2r;
245 break;
246 case HECMW_ETYPE_PYR1:
247 case HECMW_ETYPE_PYR2:
248 sid_h2r = sid_pyr_h2r;
249 break;
250 case HECMW_ETYPE_SHT1:
251 case HECMW_ETYPE_SHT2:
252 case HECMW_ETYPE_SHT6:
253 case HECMW_ETYPE_SHQ1:
254 case HECMW_ETYPE_SHQ2:
255 case HECMW_ETYPE_SHQ8:
256 case HECMW_ETYPE_PTT1:
257 case HECMW_ETYPE_PTT2:
258 case HECMW_ETYPE_PTQ1:
259 case HECMW_ETYPE_PTQ2:
260 sid_h2r = NULL;
261 break;
262 default:
264 "Element type %d not supported for refinement.\n", etype);
265 *ierror = 1;
266 sid_h2r = NULL;
267 }
268 return sid_h2r;
269}
270
271static const int *get_sid_r2h(int etype, int *ierror) {
272 /* RCAP to HECMW */
273 /*
274 static const int sid_tri_r2h[3] = {1, 2, 3};
275 static const int sid_qua_r2h[4] = {3, 2, 4, 1};
276 static const int sid_tet_r2h[4] = {1, 2, 3, 4};
277 static const int sid_pri_r2h[5] = {4, 5, 3, 1, 2};
278 static const int sid_hex_r2h[6] = {1, 2, 3, 4, 5, 6};
279 static const int sid_pyr_r2h[5] = {3, 2, 4, 1, 5};
280 */
281 /* RCAP to FSTR */
282 static const int sid_tri_r2h[3] = {2, 3, 1};
283 static const int sid_qua_r2h[4] = {1, 2, 3, 4};
284 static const int sid_tet_r2h[4] = {3, 4, 2, 1};
285 static const int sid_pri_r2h[5] = {1, 2, 3, 4, 5};
286 static const int sid_hex_r2h[6] = {1, 2, 3, 4, 5, 6};
287 static const int sid_pyr_r2h[5] = {3, 2, 4, 1,
288 5}; /* TEMPORARY: same as HECMW */
289
290 const int *sid_r2h;
291 *ierror = 0;
292 switch (etype) {
293 case HECMW_ETYPE_TRI1:
294 case HECMW_ETYPE_TRI2:
295 sid_r2h = sid_tri_r2h;
296 break;
297 case HECMW_ETYPE_QUA1:
298 case HECMW_ETYPE_QUA2:
299 sid_r2h = sid_qua_r2h;
300 break;
301 case HECMW_ETYPE_TET1:
303 case HECMW_ETYPE_TET2:
304 sid_r2h = sid_tet_r2h;
305 break;
306 case HECMW_ETYPE_PRI1:
307 case HECMW_ETYPE_PRI2:
308 sid_r2h = sid_pri_r2h;
309 break;
310 case HECMW_ETYPE_HEX1:
312 case HECMW_ETYPE_HEX2:
313 sid_r2h = sid_hex_r2h;
314 break;
315 case HECMW_ETYPE_PYR1:
316 case HECMW_ETYPE_PYR2:
317 sid_r2h = sid_pyr_r2h;
318 break;
319 case HECMW_ETYPE_SHT1:
320 case HECMW_ETYPE_SHT2:
321 case HECMW_ETYPE_SHT6:
322 case HECMW_ETYPE_SHQ1:
323 case HECMW_ETYPE_SHQ2:
324 case HECMW_ETYPE_SHQ8:
325 case HECMW_ETYPE_PTT1:
326 case HECMW_ETYPE_PTT2:
327 case HECMW_ETYPE_PTQ1:
328 case HECMW_ETYPE_PTQ2:
329 sid_r2h = NULL;
330 break;
331 default:
333 "Element type %d not supported for refinement.\n", etype);
334 *ierror = 1;
335 sid_r2h = NULL;
336 }
337 return sid_r2h;
338}
339
340static int surf_ID_hecmw2rcap(struct hecmwST_local_mesh *mesh) {
341 struct hecmwST_surf_grp *grp = mesh->surf_group;
342 int i;
343
344 if (grp->n_grp == 0) return HECMW_SUCCESS;
345
346 for (i = 0; i < grp->grp_index[grp->n_grp]; i++) {
347 int eid = grp->grp_item[2 * i];
348 int *sid_p = &(grp->grp_item[2 * i + 1]);
349 const int *sid_h2r;
350 int ierror;
351
352 HECMW_assert(0 < eid && eid <= mesh->n_elem_gross);
353
354 sid_h2r = get_sid_h2r(mesh->elem_type[eid - 1], &ierror);
355 if (ierror) return HECMW_ERROR;
356 if (sid_h2r) *sid_p = sid_h2r[*sid_p];
357 }
358 return HECMW_SUCCESS;
359}
360
361static int surf_ID_rcap2hecmw(struct hecmwST_local_mesh *mesh) {
362 struct hecmwST_surf_grp *grp = mesh->surf_group;
363 int i;
364
365 if (grp->n_grp == 0) return HECMW_SUCCESS;
366
367 for (i = 0; i < grp->grp_index[grp->n_grp]; i++) {
368 int eid = grp->grp_item[2 * i];
369 int *sid_p = &(grp->grp_item[2 * i + 1]);
370 const int *sid_r2h;
371 int ierror;
372
373 HECMW_assert(0 < eid && eid <= mesh->n_elem_gross);
374
375 sid_r2h = get_sid_r2h(mesh->elem_type[eid - 1], &ierror);
376 if (ierror) return HECMW_ERROR;
377 if (sid_r2h) *sid_p = sid_r2h[*sid_p];
378 }
379 return HECMW_SUCCESS;
380}
381
382/*============================================================================*/
383/* */
384/* for models with multiple element types */
385/* */
386/*============================================================================*/
387#define NDIV_SEG 2
388#define NDIV_SURF 4
389#define NDIV_BODY 8
390#define NDIV_OTHER 1
391#define NDIV_MAX 8
392
393static int get_elem_ndiv(int etype, int *ierror) {
394 int ndiv;
395 *ierror = 0;
396 if (HECMW_is_etype_truss(etype)) {
397 ndiv = NDIV_OTHER;
398 } else if (HECMW_is_etype_rod(etype) || HECMW_is_etype_beam(etype)) {
399 ndiv = NDIV_SEG;
400 } else if (HECMW_is_etype_surface(etype) || HECMW_is_etype_shell(etype)
401 || HECMW_is_etype_patch(etype)) {
402 ndiv = NDIV_SURF;
403 } else if (HECMW_is_etype_solid(etype)) {
404 ndiv = NDIV_BODY;
405 } else if (HECMW_is_etype_link(etype)) {
406 ndiv = NDIV_OTHER;
407 } else {
409 "Element type %d not supported for refinement.\n", etype);
410 *ierror = 1;
411 ndiv = NDIV_OTHER;
412 }
413 return ndiv;
414}
415
416/*============================================================================*/
417/* */
418/* prepare, clear and terminate REVOCAP_Refiner */
419/* */
420/*============================================================================*/
421static void register_node(struct hecmwST_local_mesh *mesh) {
422 HECMW_log(HECMW_LOG_DEBUG, "Original Node Count = %d\n", mesh->n_node_gross);
423
424#ifdef DEBUG_REFINER
425 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapSetNode64(num=%d, coords, globalIds, NULL)\n", mesh->n_node_gross);
426 {
427 char fname[256];
428 FILE *fp_dump;
429 sprintf(fname, "dump_rcapSetNode64.%d", HECMW_comm_get_rank());
430 fp_dump = fopen(fname, "w");
431 if (fp_dump == NULL) {
432 fprintf(stderr, "Error: cannot open file %s\n", fname);
433 } else {
434 HECMW_log(HECMW_LOG_DEBUG, "RCAP: writing num, coords[0..3*num-1], globalIds[0..num-1] to %s\n", fname);
435 fwrite(&(mesh->n_node_gross), sizeof(int), 1, fp_dump);
436 fwrite(mesh->node, sizeof(double), 3*mesh->n_node_gross, fp_dump);
437 fwrite(mesh->global_node_ID, sizeof(int), mesh->n_node_gross, fp_dump);
438 fclose(fp_dump);
439 }
440 }
441#endif
442 rcapSetNode64(mesh->n_node_gross, mesh->node, mesh->global_node_ID, NULL);
443
444 HECMW_assert(rcapGetNodeCount() == mesh->n_node_gross);
445}
446
447static void register_node_groups(struct hecmwST_node_grp *grp) {
448 int i;
449 char rcap_name[80];
450
451 HECMW_assert(strcmp(grp->grp_name[0], "ALL") == 0);
452
453 for (i = 1; i < grp->n_grp; i++) {
454 int start = grp->grp_index[i];
455 int num = grp->grp_index[i + 1] - start;
456 int *array = grp->grp_item + start;
457 sprintf(rcap_name, "NG_%s", grp->grp_name[i]);
458
459#ifdef DEBUG_REFINER
460 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapAppendNodeGroup(dataname=%s, num=%d, nodeArray)\n", rcap_name, num);
461 {
462 char fname[256];
463 FILE *fp_dump;
464 sprintf(fname, "dump_rcapAppendNodeGroup_%s.%d", grp->grp_name[i], HECMW_comm_get_rank());
465 fp_dump = fopen(fname, "w");
466 if (fp_dump == NULL) {
467 fprintf(stderr, "Error: cannot open file %s\n", fname);
468 } else {
469 int len = strlen(rcap_name);
470 HECMW_log(HECMW_LOG_DEBUG, "RCAP: writing strlen(dataname), dataname[], num, nodeArray[0..num-1] to %s\n", fname);
471 fwrite(&len, sizeof(int), 1, fp_dump);
472 fwrite(rcap_name, sizeof(char), len, fp_dump);
473 fwrite(&num, sizeof(int), 1, fp_dump);
474 fwrite(array, sizeof(int), num, fp_dump);
475 fclose(fp_dump);
476 }
477 }
478#endif
479 rcapAppendNodeGroup(rcap_name, num, array);
480 }
481}
482
483static int register_surf_groups(struct hecmwST_local_mesh *mesh) {
484 struct hecmwST_surf_grp *grp = mesh->surf_group;
485 int i, j;
486 char rcap_name[80];
487
488 struct hecmw_varray_int other;
489
490 for (i = 0; i < grp->n_grp; i++) {
491 int start = grp->grp_index[i];
492 int num = grp->grp_index[i + 1] - start;
493 int *array = grp->grp_item + start * 2;
494
495 if (HECMW_varray_int_init(&other) != HECMW_SUCCESS) return HECMW_ERROR;
496
497 for (j = 0; j < num; j++) {
498 int elem = array[j * 2];
499 int surf = array[j * 2 + 1];
500 int etype = mesh->elem_type[elem - 1];
501
502 /* ignore shell/patch surface */
503 if (HECMW_is_etype_shell(etype) || HECMW_is_etype_patch(etype)) continue;
504
505 if (HECMW_varray_int_append(&other, elem) != HECMW_SUCCESS)
506 return HECMW_ERROR;
507 if (HECMW_varray_int_append(&other, surf) != HECMW_SUCCESS)
508 return HECMW_ERROR;
509 }
510
511 sprintf(rcap_name, "SG_%s", grp->grp_name[i]);
512 num = HECMW_varray_int_nval(&other) / 2;
513 array = HECMW_varray_int_get_v(&other);
514
515#ifdef DEBUG_REFINER
516 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapAppendFaceGroup(dataname=%s, num=%d, faceArray)\n", rcap_name, num);
517 {
518 char fname[256];
519 FILE *fp_dump;
520 sprintf(fname, "dump_rcapAppendFaceGroup_%s.%d", grp->grp_name[i], HECMW_comm_get_rank());
521 fp_dump = fopen(fname, "w");
522 if (fp_dump == NULL) {
523 fprintf(stderr, "Error: cannot open file %s\n", fname);
524 } else {
525 int len = strlen(rcap_name);
526 HECMW_log(HECMW_LOG_DEBUG, "RCAP: writing strlen(dataname), dataname[], num, faceArray[0..2*num-1] to %s\n", fname);
527 fwrite(&len, sizeof(int), 1, fp_dump);
528 fwrite(rcap_name, sizeof(char), len, fp_dump);
529 fwrite(&num, sizeof(int), 1, fp_dump);
530 fwrite(array, sizeof(int), num*2, fp_dump);
531 fclose(fp_dump);
532 }
533 }
534#endif
535 rcapAppendFaceGroup(rcap_name, num, array);
536
538 }
539 return HECMW_SUCCESS;
540}
541
542static int prepare_refiner(struct hecmwST_local_mesh *mesh,
543 const char *cad_filename,
544 const char *part_filename) {
545 HECMW_log(HECMW_LOG_DEBUG, "Started preparing refiner...\n");
546
547 if (elem_node_item_hecmw2rcap(mesh) != HECMW_SUCCESS) return HECMW_ERROR;
548 if (surf_ID_hecmw2rcap(mesh) != HECMW_SUCCESS) return HECMW_ERROR;
549
550#ifdef DEBUG_REFINER
551 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapInitRefiner(nodeOffset=1, elementOffset=1)\n");
552#endif
553 rcapInitRefiner(1, 1);
554 if (cad_filename != NULL) {
555#ifdef DEBUG_REFINER
556 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapSetCADFilename(filename=%s)\n", cad_filename);
557#endif
558 rcapSetCADFilename(cad_filename);
559 }
560 if (part_filename != NULL) {
561#ifdef DEBUG_REFINER
562 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapSetPartitionFilename(filename=%s)\n", part_filename);
563#endif
564 rcapSetPartitionFilename(part_filename);
565 }
566 register_node(mesh);
567
568 /* node groups are refined by REVOCAP_Refiner */
569 register_node_groups(mesh->node_group);
570
571 /* element groups are refined by myself */
572
573 /* surface groups are refined by REVOCAP_Refiner except for shell/patch surface */
574 if (register_surf_groups(mesh) != HECMW_SUCCESS) return HECMW_ERROR;
575
576 HECMW_log(HECMW_LOG_DEBUG, "Finished preparing refiner.\n");
577 return HECMW_SUCCESS;
578}
579
580static void clear_refiner(void) {
581#ifdef DEBUG_REFINER
582 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapClearRefiner()\n");
583#endif
584 rcapClearRefiner();
585}
586
587static int terminate_refiner(struct hecmwST_local_mesh *mesh) {
588#ifdef DEBUG_REFINER
589 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapTermRefiner()\n");
590#endif
591 rcapTermRefiner();
592 if (surf_ID_rcap2hecmw(mesh) != HECMW_SUCCESS) return HECMW_ERROR;
593 if (elem_node_item_rcap2hecmw(mesh) != HECMW_SUCCESS) return HECMW_ERROR;
594 return HECMW_SUCCESS;
595}
596
597/*============================================================================*/
598/* */
599/* call REVOCAP_Refiner */
600/* */
601/*============================================================================*/
602static int elem_type_hecmw2rcap(int etype) {
603 switch (etype) {
604 case HECMW_ETYPE_ROD1:
605 return RCAP_SEGMENT;
606 case HECMW_ETYPE_ROD2:
607 return RCAP_SEGMENT2;
608 case HECMW_ETYPE_TRI1:
609 return RCAP_TRIANGLE;
610 case HECMW_ETYPE_TRI2:
611 return RCAP_TRIANGLE2;
612 case HECMW_ETYPE_QUA1:
613 return RCAP_QUAD;
614 case HECMW_ETYPE_QUA2:
615 return RCAP_QUAD2;
617 return RCAP_SEGMENT;
618 case HECMW_ETYPE_TET1:
619 return RCAP_TETRAHEDRON;
621 return RCAP_TETRAHEDRON;
622 case HECMW_ETYPE_TET2:
623 return RCAP_TETRAHEDRON2;
624 case HECMW_ETYPE_PRI1:
625 return RCAP_WEDGE;
626 case HECMW_ETYPE_PRI2:
627 return RCAP_WEDGE2;
628 case HECMW_ETYPE_PYR1:
629 return RCAP_PYRAMID;
630 case HECMW_ETYPE_PYR2:
631 return RCAP_PYRAMID2;
632 case HECMW_ETYPE_HEX1:
633 return RCAP_HEXAHEDRON;
635 return RCAP_HEXAHEDRON;
636 case HECMW_ETYPE_HEX2:
637 return RCAP_HEXAHEDRON2;
638 case HECMW_ETYPE_BEM1:
639 return RCAP_SEGMENT;
640 case HECMW_ETYPE_BEM2:
641 return RCAP_SEGMENT2;
642 case HECMW_ETYPE_BEM3:
643 return RCAP_SEGMENT;
644 case HECMW_ETYPE_SHT1:
645 return RCAP_TRIANGLE;
646 case HECMW_ETYPE_SHT2:
647 return RCAP_TRIANGLE2;
648 case HECMW_ETYPE_SHT6:
649 return RCAP_TRIANGLE;
650 case HECMW_ETYPE_SHQ1:
651 return RCAP_QUAD;
652 case HECMW_ETYPE_SHQ2:
653 return RCAP_QUAD2;
654 case HECMW_ETYPE_SHQ8:
655 return RCAP_QUAD;
656 case HECMW_ETYPE_PTT1:
657 return RCAP_TRIANGLE;
658 case HECMW_ETYPE_PTT2:
659 return RCAP_TRIANGLE2;
660 case HECMW_ETYPE_PTQ1:
661 return RCAP_QUAD;
662 case HECMW_ETYPE_PTQ2:
663 return RCAP_QUAD2;
664 }
665 return RCAP_UNKNOWNTYPE;
666}
667
668static int reorder_enode_33struct(int etype, int n_elem,
669 int *elem_node_item_ref) {
670 int nn, nn_2, ndiv, nn_tot, ierror;
671 int j, k, l;
672 int tmp[32];
673 int *enode1, *enode2;
674 nn = HECMW_get_max_node(etype);
675 nn_2 = nn / 2;
676 ndiv = get_elem_ndiv(etype, &ierror);
677 if (ierror) return HECMW_ERROR;
678 nn_tot = ndiv * nn;
679 enode1 = elem_node_item_ref;
680 enode2 = enode1 + nn_2 * ndiv;
681 for (j = 0; j < n_elem; j++) {
682 int *tmp_p = tmp;
683 int *en1 = enode1;
684 for (k = 0; k < ndiv; k++) {
685 for (l = 0; l < nn_2; l++) {
686 tmp_p[l] = enode1[l];
687 tmp_p[nn_2 + l] = enode2[l];
688 }
689 tmp_p += nn;
690 enode1 += nn_2;
691 enode2 += nn_2;
692 }
693 for (k = 0; k < nn_tot; k++) {
694 en1[k] = tmp[k];
695 }
696 enode1 += nn_2 * ndiv;
697 enode2 += nn_2 * ndiv;
698 }
699 return HECMW_SUCCESS;
700}
701
702static int refine_element(struct hecmwST_local_mesh *mesh,
703 struct hecmwST_local_mesh *ref_mesh) {
704 int i, j;
705 int n_elem_ref_tot;
706 size_t n_enode_ref_tot;
707 int *elem_node_index_ref;
708 int *elem_node_item_ref;
709
710 HECMW_log(HECMW_LOG_DEBUG, "Original Element Count = %d\n", mesh->n_elem_gross);
711
712 /*
713 * count num of elems after refinement
714 */
715 n_elem_ref_tot = 0;
716 n_enode_ref_tot = 0;
717 for (i = 0; i < mesh->n_elem_type; i++) {
718 int etype, istart, iend, n_elem, etype_rcap, nn;
719 size_t n_elem_ref;
720 int *elem_node_item;
721 etype = mesh->elem_type_item[i];
722 istart = mesh->elem_type_index[i];
723 iend = mesh->elem_type_index[i + 1];
724 n_elem = iend - istart;
725 if (HECMW_is_etype_link(etype) || HECMW_is_etype_truss(etype)) {
726 n_elem_ref = n_elem;
727 } else {
728 int is33 = HECMW_is_etype_33struct(etype);
729 if (is33) n_elem *= 2;
730 etype_rcap = elem_type_hecmw2rcap(etype);
731 HECMW_assert(etype_rcap != RCAP_UNKNOWNTYPE);
732 elem_node_item = mesh->elem_node_item + mesh->elem_node_index[istart];
733
734#ifdef DEBUG_REFINER
735 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapRefineElement(num=%d, etype=%d, nodeArray, NULL)\n", n_elem, etype_rcap);
736 {
737 char fname[256];
738 FILE *fp_dump;
739 sprintf(fname, "dump_rcapRefineElement_pre_%d.%d", etype, HECMW_comm_get_rank());
740 fp_dump = fopen(fname, "w");
741 if (fp_dump == NULL) {
742 fprintf(stderr, "Error: cannot open file %s\n", fname);
743 } else {
744 int nn = HECMW_get_max_node(etype);
745 HECMW_log(HECMW_LOG_DEBUG, "RCAP: writing num, etype, nodeArray[0..%d*num-1] to %s\n", nn, fname);
746 fwrite(&n_elem, sizeof(int), 1, fp_dump);
747 fwrite(&etype_rcap, sizeof(int), 1, fp_dump);
748 fwrite(elem_node_item, sizeof(int), nn * n_elem, fp_dump);
749 fclose(fp_dump);
750 }
751 }
752#endif
753 n_elem_ref = rcapRefineElement(n_elem, etype_rcap, elem_node_item, NULL);
754#ifdef DEBUG_REFINER
755 HECMW_log(HECMW_LOG_DEBUG, "RCAP: return value was %d\n", n_elem_ref);
756#endif
757 if (is33) {
758 n_elem_ref /= 2;
759 n_elem /= 2;
760 }
761 }
762#ifdef DEBUG
763 {
764 int ndiv, ierror;
765 ndiv = get_elem_ndiv(etype, &ierror);
766 if (ierror) return HECMW_ERROR;
767 HECMW_assert(n_elem_ref == n_elem * ndiv);
768 }
769#endif
770 n_elem_ref_tot += n_elem_ref;
771 nn = HECMW_get_max_node(etype);
772 n_enode_ref_tot += n_elem_ref * nn;
773 }
774 ref_mesh->n_elem_gross = n_elem_ref_tot;
775
776 HECMW_log(HECMW_LOG_DEBUG, "Refined Element Count will be %d\n", ref_mesh->n_elem_gross);
777
778 /*
779 * allocate memory
780 */
781 ref_mesh->elem_node_index =
782 (int *)HECMW_malloc(sizeof(int) * (n_elem_ref_tot + 1));
783 if (ref_mesh->elem_node_index == NULL) {
784 HECMW_set_error(errno, "");
785 return HECMW_ERROR;
786 }
787 ref_mesh->elem_node_item = (int *)HECMW_calloc(n_enode_ref_tot, sizeof(int));
788 if (ref_mesh->elem_node_item == NULL) {
789 HECMW_set_error(errno, "");
790 return HECMW_ERROR;
791 }
792
793 HECMW_assert(rcapGetNodeCount() == mesh->n_node_gross);
794
795 /*
796 * perform actual refinement
797 */
798 ref_mesh->elem_node_index[0] = 0;
799 elem_node_index_ref = ref_mesh->elem_node_index;
800 elem_node_item_ref = ref_mesh->elem_node_item;
801 n_elem_ref_tot = 0;
802 for (i = 0; i < mesh->n_elem_type; i++) {
803 int etype, istart, iend, n_elem, etype_rcap, n_elem_ref, nn, jstart;
804 int *elem_node_item;
805 etype = mesh->elem_type_item[i];
806 istart = mesh->elem_type_index[i];
807 iend = mesh->elem_type_index[i + 1];
808 n_elem = iend - istart;
809 nn = HECMW_get_max_node(etype);
810 if (HECMW_is_etype_link(etype) || HECMW_is_etype_truss(etype)) {
811 n_elem_ref = n_elem;
812 for (j = 0; j < n_elem * nn; j++) {
813 /* elem_node_item_ref[j] = elem_node_item[j]; */
814 jstart = mesh->elem_node_index[istart];
815 /* printf("aa %d %d\n",jstart,j); */
816 elem_node_item_ref[j] = mesh->elem_node_item[jstart + j];
817 }
818 } else {
819 int is33 = HECMW_is_etype_33struct(etype);
820 if (is33) n_elem *= 2;
821 etype_rcap = elem_type_hecmw2rcap(etype);
822 elem_node_item = mesh->elem_node_item + mesh->elem_node_index[istart];
823
824#ifdef DEBUG_REFINER
825 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapRefineElement(num=%d, etype=%d, nodeArray, resultNodeArray)\n", n_elem, etype_rcap);
826#endif
827 n_elem_ref = rcapRefineElement(n_elem, etype_rcap, elem_node_item,
828 elem_node_item_ref);
829#ifdef DEBUG_REFINER
830 {
831 char fname[256];
832 FILE *fp_dump;
833 sprintf(fname, "dump_rcapRefineElement_%d.%d", etype, HECMW_comm_get_rank());
834 fp_dump = fopen(fname, "w");
835 if (fp_dump == NULL) {
836 fprintf(stderr, "Error: cannot open file %s\n", fname);
837 } else {
838 int nn = HECMW_get_max_node(etype);
840 "RCAP: writing num, etype, nodeArray[0..%d*num-1], num_ref, resultNodeArray[0..%d*num_ref-1] to %s\n",
841 nn, nn, fname);
842 fwrite(&n_elem, sizeof(int), 1, fp_dump);
843 fwrite(&etype_rcap, sizeof(int), 1, fp_dump);
844 fwrite(elem_node_item, sizeof(int), nn * n_elem, fp_dump);
845 fwrite(&n_elem_ref, sizeof(int), 1, fp_dump);
846 fwrite(elem_node_item_ref, sizeof(int), nn * n_elem_ref, fp_dump);
847 fclose(fp_dump);
848 }
849 }
850#endif
851
852 if (is33) {
853 n_elem /= 2;
854 n_elem_ref /= 2;
855 reorder_enode_33struct(etype, n_elem, elem_node_item_ref);
856 }
857 }
858 for (j = 0; j < n_elem_ref; j++) {
859 elem_node_index_ref[j + 1] = elem_node_index_ref[j] + nn;
860 if (elem_node_index_ref[j + 1] < elem_node_index_ref[j]) {
862 "Integer overflow detected while creating elem_node_index\n");
863 return HECMW_ERROR;
864 }
865 }
866 elem_node_index_ref += n_elem_ref;
867 elem_node_item_ref += nn * n_elem_ref;
868 n_elem_ref_tot += n_elem_ref;
869 }
870
871 HECMW_log(HECMW_LOG_DEBUG, "Refined Element Count = %d\n", n_elem_ref_tot);
872
873 HECMW_assert(n_elem_ref_tot == ref_mesh->n_elem_gross);
874
875#ifdef DEBUG_REFINER
876 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapCommit()\n");
877#endif
878 rcapCommit();
879 return HECMW_SUCCESS;
880}
881
882static int refine_node(struct hecmwST_local_mesh *mesh,
883 struct hecmwST_local_mesh *ref_mesh) {
884 /* n_node_gross */
885#ifdef DEBUG_REFINER
886 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapGetNodeCount()\n");
887#endif
888 ref_mesh->n_node_gross = rcapGetNodeCount();
889#ifdef DEBUG_REFINER
890 HECMW_log(HECMW_LOG_DEBUG, "RCAP: return value was %d\n", ref_mesh->n_node_gross);
891#endif
892
893 HECMW_log(HECMW_LOG_DEBUG, "Refined Node Count = %d\n", ref_mesh->n_node_gross);
894
896
897 /* node */
898 ref_mesh->node =
899 (double *)HECMW_malloc(sizeof(double) * ref_mesh->n_node_gross * 3);
900 if (ref_mesh->node == NULL) {
901 HECMW_set_error(errno, "");
902 return HECMW_ERROR;
903 }
904#ifdef DEBUG_REFINER
905 HECMW_log(HECMW_LOG_DEBUG, "RCAP: calling rcapGetNodeSeq64(num=%d, initId=1, coords)\n", ref_mesh->n_node_gross);
906#endif
907 rcapGetNodeSeq64(ref_mesh->n_node_gross, 1, ref_mesh->node);
908#ifdef DEBUG_REFINER
909 {
910 char fname[256];
911 FILE *fp_dump;
912 sprintf(fname, "dump_rcapGetNodeSeq64.%d", HECMW_comm_get_rank());
913 fp_dump = fopen(fname, "w");
914 if (fp_dump == NULL) {
915 fprintf(stderr, "Error: cannot open file %s\n", fname);
916 } else {
917 int one = 1;
918 HECMW_log(HECMW_LOG_DEBUG, "RCAP: writing num, initId, coords[0..3*num-1] to %s\n", fname);
919 fwrite(&(ref_mesh->n_node_gross), sizeof(int), 1, fp_dump);
920 fwrite(&one, sizeof(int), 1, fp_dump);
921 fwrite(ref_mesh->node, sizeof(int), 3 * ref_mesh->n_node_gross, fp_dump);
922 fclose(fp_dump);
923 }
924 }
925#endif
926
927 return HECMW_SUCCESS;
928}
929
930static int refine_node_group_info(struct hecmwST_node_grp *grp,
931 struct hecmwST_node_grp *ref_grp,
932 int ref_n_node_gross) {
933 int i;
934 char rcap_name[80];
935
936 ref_grp->n_grp = grp->n_grp;
937
938 /* grp_name (COPY) */
939 ref_grp->grp_name = (char **)HECMW_malloc(sizeof(char *) * grp->n_grp);
940 if (ref_grp->grp_name == NULL) {
941 HECMW_set_error(errno, "");
942 return HECMW_ERROR;
943 }
944 /* ALL */
945 ref_grp->grp_name[0] = HECMW_strdup("ALL");
946 if (ref_grp->grp_name[0] == NULL) {
947 HECMW_set_error(errno, "");
948 return HECMW_ERROR;
949 }
950 /* other groups */
951 for (i = 1; i < grp->n_grp; i++) {
952 ref_grp->grp_name[i] = HECMW_strdup(grp->grp_name[i]);
953 if (ref_grp->grp_name[i] == NULL) {
954 HECMW_set_error(errno, "");
955 return HECMW_ERROR;
956 }
957 }
958
959 /* grp_index */
960 ref_grp->grp_index = (int *)HECMW_malloc(sizeof(int) * (grp->n_grp + 1));
961 if (ref_grp->grp_index == NULL) {
962 HECMW_set_error(errno, "");
963 return HECMW_ERROR;
964 }
965 ref_grp->grp_index[0] = 0;
966 /* ALL */
967 ref_grp->grp_index[1] = ref_n_node_gross;
968 /* other groups */
969 for (i = 1; i < grp->n_grp; i++) {
970 sprintf(rcap_name, "NG_%s", grp->grp_name[i]);
971 ref_grp->grp_index[i + 1] =
972 ref_grp->grp_index[i] + rcapGetNodeGroupCount(rcap_name);
973 }
974
975 /* grp_item */
976 if (ref_grp->grp_index[grp->n_grp] > 0) {
977 ref_grp->grp_item =
978 (int *)HECMW_calloc(ref_grp->grp_index[grp->n_grp], sizeof(int));
979 if (ref_grp->grp_item == NULL) {
980 HECMW_set_error(errno, "");
981 return HECMW_ERROR;
982 }
983 /* ALL */
984 for (i = 0; i < ref_n_node_gross; i++) {
985 ref_grp->grp_item[i] = i + 1;
986 }
987 /* other groups */
988 for (i = 1; i < grp->n_grp; i++) {
989 int start = ref_grp->grp_index[i];
990 int num = ref_grp->grp_index[i + 1] - start;
991 if (num > 0) {
992 int *array = ref_grp->grp_item + start;
993 sprintf(rcap_name, "NG_%s", grp->grp_name[i]);
994 rcapGetNodeGroup(rcap_name, num, array);
995 }
996 }
997 } else {
998 ref_grp->grp_item = NULL;
999 }
1000 return HECMW_SUCCESS;
1001}
1002
1003/*
1004 * static functions for refining element groups
1005 */
1006
1007static int get_refined_element_count(const struct hecmwST_local_mesh *mesh,
1008 int n_elem, int *elem_list) {
1009 int i, eid, etype, ierror, ndiv;
1010 int n_elem_ref = 0;
1011 for (i = 0; i < n_elem; i++) {
1012 eid = elem_list[i];
1013 HECMW_assert(0 < eid && eid <= mesh->n_elem_gross);
1014 etype = mesh->elem_type[eid - 1];
1015 ndiv = get_elem_ndiv(etype, &ierror);
1016 if (ierror) return -1;
1017 n_elem_ref += ndiv;
1018 }
1019 return n_elem_ref;
1020}
1021
1022static int get_refined_elements(const struct hecmwST_local_mesh *mesh,
1023 struct hecmwST_local_mesh *ref_mesh, int eid,
1024 int *eid_ref) {
1025 int i, etype, ierror, ndiv;
1026 int etype_idx, eid_type;
1027 int offset;
1028
1029 HECMW_assert(0 < eid);
1030 HECMW_assert(eid <= mesh->n_elem_gross);
1031 HECMW_assert(eid_ref);
1032
1033 etype = mesh->elem_type[eid - 1];
1034 ndiv = get_elem_ndiv(etype, &ierror);
1035 if (ierror) return -1;
1036 etype_idx = -1;
1037 eid_type = -1;
1038 for (i = 0; i < mesh->n_elem_type; i++) {
1039 if (mesh->elem_type_item[i] == etype) {
1040 etype_idx = i;
1042 HECMW_assert(eid <= mesh->elem_type_index[i + 1]);
1043 eid_type = eid - mesh->elem_type_index[i];
1044 break;
1045 }
1046 }
1047 HECMW_assert(etype_idx >= 0);
1048 HECMW_assert(eid_type > 0);
1049
1050 offset = ref_mesh->elem_type_index[etype_idx] + (eid_type - 1) * ndiv + 1;
1051 for (i = 0; i < ndiv; i++) {
1052 eid_ref[i] = offset + i;
1053 }
1054 return ndiv;
1055}
1056
1057static int get_refined_elem_list(const struct hecmwST_local_mesh *mesh,
1058 struct hecmwST_local_mesh *ref_mesh,
1059 int n_elem, int *elem_list,
1060 int *elem_list_ref) {
1061 int i, j;
1062 int eid_ref[NDIV_MAX];
1063 int *el_ref = elem_list_ref;
1064 int n_elem_ref = 0;
1065 for (i = 0; i < n_elem; i++) {
1066 int ndiv = get_refined_elements(mesh, ref_mesh, elem_list[i], eid_ref);
1067 if (ndiv < 0) return -1;
1068 for (j = 0; j < ndiv; j++) {
1069 el_ref[j] = eid_ref[j];
1070 }
1071 el_ref += ndiv;
1072 n_elem_ref += ndiv;
1073 }
1074 return n_elem_ref;
1075}
1076
1077/*
1078 * refinement of element groups is done by myself
1079 */
1080
1081static int refine_elem_group_info(struct hecmwST_local_mesh *mesh,
1082 struct hecmwST_local_mesh *ref_mesh) {
1083 int i;
1084 struct hecmwST_elem_grp *grp = mesh->elem_group;
1085 struct hecmwST_elem_grp *ref_grp = ref_mesh->elem_group;
1086 int ref_n_elem_gross = ref_mesh->n_elem_gross;
1087
1088 ref_grp->n_grp = grp->n_grp;
1089
1090 /* grp_name (COPY) */
1091 ref_grp->grp_name = (char **)HECMW_malloc(sizeof(char *) * grp->n_grp);
1092 if (ref_grp->grp_name == NULL) {
1093 HECMW_set_error(errno, "");
1094 return HECMW_ERROR;
1095 }
1096 /* ALL */
1097 HECMW_assert(strcmp(grp->grp_name[0], "ALL") == 0);
1098 ref_grp->grp_name[0] = HECMW_strdup("ALL");
1099 if (ref_grp->grp_name[0] == NULL) {
1100 HECMW_set_error(errno, "");
1101 return HECMW_ERROR;
1102 }
1103 /* other groups */
1104 for (i = 1; i < grp->n_grp; i++) {
1105 ref_grp->grp_name[i] = HECMW_strdup(grp->grp_name[i]);
1106 if (ref_grp->grp_name[i] == NULL) {
1107 HECMW_set_error(errno, "");
1108 return HECMW_ERROR;
1109 }
1110 }
1111
1112 /* grp_index */
1113 ref_grp->grp_index = (int *)HECMW_malloc(sizeof(int) * (grp->n_grp + 1));
1114 if (ref_grp->grp_index == NULL) {
1115 HECMW_set_error(errno, "");
1116 return HECMW_ERROR;
1117 }
1118 ref_grp->grp_index[0] = 0;
1119 /* ALL */
1120 ref_grp->grp_index[1] = ref_n_elem_gross;
1121 /* other groups */
1122 for (i = 1; i < grp->n_grp; i++) {
1123 int start, n_elem, n_elem_ref;
1124 int *elem_list;
1125 start = grp->grp_index[i];
1126 n_elem = grp->grp_index[i + 1] - start;
1127 elem_list = grp->grp_item + start;
1128 n_elem_ref = get_refined_element_count(mesh, n_elem, elem_list);
1129 if (n_elem_ref < 0) return HECMW_ERROR;
1130 ref_grp->grp_index[i + 1] = ref_grp->grp_index[i] + n_elem_ref;
1131 }
1132
1133 /* grp_item */
1134 if (ref_grp->grp_index[grp->n_grp] > 0) {
1135 ref_grp->grp_item =
1136 (int *)HECMW_calloc(ref_grp->grp_index[grp->n_grp], sizeof(int));
1137 if (ref_grp->grp_item == NULL) {
1138 HECMW_set_error(errno, "");
1139 return HECMW_ERROR;
1140 }
1141 /* ALL */
1142 for (i = 0; i < ref_n_elem_gross; i++) {
1143 ref_grp->grp_item[i] = i + 1;
1144 }
1145 /* other groups */
1146 for (i = 1; i < grp->n_grp; i++) {
1147 int start, start_ref, n_elem, n_elem_ref, ret;
1148 int *elem_list, *elem_list_ref;
1149 start_ref = ref_grp->grp_index[i];
1150 n_elem_ref = ref_grp->grp_index[i + 1] - start_ref;
1151 HECMW_assert(n_elem_ref >= 0);
1152 if (n_elem_ref == 0) continue;
1153 start = grp->grp_index[i];
1154 n_elem = grp->grp_index[i + 1] - start;
1155 elem_list = grp->grp_item + start;
1156 elem_list_ref = ref_grp->grp_item + start_ref;
1157 ret = get_refined_elem_list(mesh, ref_mesh, n_elem, elem_list,
1158 elem_list_ref);
1159 HECMW_assert(ret == n_elem_ref);
1160 }
1161 } else {
1162 ref_grp->grp_item = NULL;
1163 }
1164 return HECMW_SUCCESS;
1165}
1166
1167/*
1168 * refinement of surface groups except for shell/patch surface is done
1169 * by REVOCAP_Refiner
1170 * refinement of shell/patch surfaces is done by myself
1171 */
1172
1173static int refine_surf_group_info(struct hecmwST_local_mesh *mesh,
1174 struct hecmwST_local_mesh *ref_mesh) {
1175 struct hecmwST_surf_grp *grp = mesh->surf_group;
1176 struct hecmwST_surf_grp *ref_grp = ref_mesh->surf_group;
1177 int i, j;
1178 char rcap_name[80];
1179
1180 ref_grp->n_grp = grp->n_grp;
1181
1182 /* grp_name (COPY) */
1183 ref_grp->grp_name = (char **)HECMW_malloc(sizeof(char *) * grp->n_grp);
1184 if (ref_grp->grp_name == NULL) {
1185 HECMW_set_error(errno, "");
1186 return HECMW_ERROR;
1187 }
1188 for (i = 0; i < grp->n_grp; i++) {
1189 ref_grp->grp_name[i] = HECMW_strdup(grp->grp_name[i]);
1190 if (ref_grp->grp_name[i] == NULL) {
1191 HECMW_set_error(errno, "");
1192 return HECMW_ERROR;
1193 }
1194 }
1195
1196 /* grp_index */
1197 ref_grp->grp_index = (int *)HECMW_malloc(sizeof(int) * (grp->n_grp + 1));
1198 if (ref_grp->grp_index == NULL) {
1199 HECMW_set_error(errno, "");
1200 return HECMW_ERROR;
1201 }
1202 ref_grp->grp_index[0] = 0;
1203 for (i = 0; i < grp->n_grp; i++) {
1204 int start = grp->grp_index[i];
1205 int num = grp->grp_index[i + 1] - start;
1206 int *array = grp->grp_item + start * 2;
1207 int num_sh, num_other;
1208
1209 /* count surfaces except for shell/patch */
1210 sprintf(rcap_name, "SG_%s", grp->grp_name[i]);
1211 num_other = rcapGetFaceGroupCount(rcap_name);
1212
1213 /* count shell/patch surfaces */
1214 num_sh = 0;
1215 for (j = 0; j < num; j++) {
1216 int elem = array[j * 2];
1217 int etype = mesh->elem_type[elem - 1];
1218 if (!(HECMW_is_etype_shell(etype) || HECMW_is_etype_patch(etype))) continue;
1219 num_sh += 1;
1220 }
1221 num_sh *= NDIV_SURF;
1222
1223 ref_grp->grp_index[i + 1] = ref_grp->grp_index[i] + num_other + num_sh;
1224 }
1225 HECMW_assert(ref_grp->grp_index[grp->n_grp] >= 0);
1226
1227 /* grp_item */
1228 if (ref_grp->grp_index[grp->n_grp] == 0) {
1229 ref_grp->grp_item = NULL;
1230 return HECMW_SUCCESS;
1231 }
1232 ref_grp->grp_item =
1233 (int *)HECMW_calloc(ref_grp->grp_index[grp->n_grp] * 2, sizeof(int));
1234 if (ref_grp->grp_item == NULL) {
1235 HECMW_set_error(errno, "");
1236 return HECMW_ERROR;
1237 }
1238 for (i = 0; i < grp->n_grp; i++) {
1239 struct hecmw_varray_int sh1, sh2;
1240 int start, num_tot, n_elem, ret;
1241 int start_ref, num_ref, num_tot_ref, n_elem_ref;
1242 int *array, *array_ref, *elem_list, *elem_list_ref;
1243
1244 start_ref = ref_grp->grp_index[i];
1245 num_tot_ref = ref_grp->grp_index[i + 1] - start_ref;
1246 if (num_tot_ref == 0) continue;
1247 array_ref = ref_grp->grp_item + start_ref * 2;
1248
1249 /* get surfaces from REVOCAP_Refiner */
1250 sprintf(rcap_name, "SG_%s", grp->grp_name[i]);
1251 num_ref = rcapGetFaceGroupCount(rcap_name);
1252 rcapGetFaceGroup(rcap_name, num_ref, array_ref);
1253
1254 num_tot_ref -= num_ref;
1255 if (num_tot_ref == 0) continue;
1256 array_ref += num_ref * 2;
1257
1258 /*
1259 * make refined shell/patch surfaces by myself
1260 */
1261
1262 if (HECMW_varray_int_init(&sh1) != HECMW_SUCCESS) return HECMW_ERROR;
1263 if (HECMW_varray_int_init(&sh2) != HECMW_SUCCESS) return HECMW_ERROR;
1264
1265 /* collect shell/patch surface */
1266 start = grp->grp_index[i];
1267 num_tot = grp->grp_index[i + 1] - start;
1268 array = grp->grp_item + start * 2;
1269 for (j = 0; j < num_tot; j++) {
1270 int elem, etype, surf;
1271 elem = array[j * 2];
1272 surf = array[j * 2 + 1];
1273 etype = mesh->elem_type[elem - 1];
1274 if (!(HECMW_is_etype_shell(etype) || HECMW_is_etype_patch(etype))) continue;
1275 if (surf == 1) {
1276 if (HECMW_varray_int_append(&sh1, elem) != HECMW_SUCCESS)
1277 return HECMW_ERROR;
1278 } else {
1279 HECMW_assert(surf == 2);
1280 if (HECMW_varray_int_append(&sh2, elem) != HECMW_SUCCESS)
1281 return HECMW_ERROR;
1282 }
1283 }
1284
1285 /* 1st surface of shells/patches */
1286 n_elem = HECMW_varray_int_nval(&sh1);
1287 if (n_elem > 0) {
1289 n_elem_ref = n_elem * NDIV_SURF;
1290 elem_list_ref = (int *)HECMW_calloc(n_elem_ref, sizeof(int));
1291 if (elem_list_ref == NULL) {
1292 HECMW_set_error(errno, "");
1293 return HECMW_ERROR;
1294 }
1295 ret = get_refined_elem_list(mesh, ref_mesh, n_elem, elem_list, elem_list_ref);
1296 HECMW_assert(ret == n_elem_ref);
1297 for (j = 0; j < n_elem_ref; j++) {
1298 array_ref[j * 2] = elem_list_ref[j];
1299 array_ref[j * 2 + 1] = 1;
1300 }
1301 HECMW_free(elem_list_ref);
1302
1303 num_tot_ref -= n_elem_ref;
1304 array_ref += n_elem_ref * 2;
1305 }
1306
1307 /* 2nd surface of shells */
1308 n_elem = HECMW_varray_int_nval(&sh2);
1309 if (n_elem > 0) {
1311 n_elem_ref = n_elem * NDIV_SURF;
1312 elem_list_ref = (int *)HECMW_calloc(n_elem_ref, sizeof(int));
1313 if (elem_list_ref == NULL) {
1314 HECMW_set_error(errno, "");
1315 return HECMW_ERROR;
1316 }
1317 ret = get_refined_elem_list(mesh, ref_mesh, n_elem, elem_list, elem_list_ref);
1318 HECMW_assert(ret == n_elem_ref);
1319 for (j = 0; j < n_elem_ref; j++) {
1320 array_ref[j * 2] = elem_list_ref[j];
1321 array_ref[j * 2 + 1] = 2;
1322 }
1323 HECMW_free(elem_list_ref);
1324
1325 num_tot_ref -= n_elem_ref;
1326 }
1327
1330
1331 HECMW_assert(num_tot_ref == 0);
1332 }
1333 return HECMW_SUCCESS;
1334}
1335
1336static int call_refiner(struct hecmwST_local_mesh *mesh,
1337 struct hecmwST_local_mesh *ref_mesh) {
1338 HECMW_log(HECMW_LOG_DEBUG, "Started calling refiner...\n");
1339
1340 if (refine_element(mesh, ref_mesh) != HECMW_SUCCESS) {
1341 return HECMW_ERROR;
1342 }
1343 if (refine_node(mesh, ref_mesh) != HECMW_SUCCESS) {
1344 return HECMW_ERROR;
1345 }
1346 if (refine_node_group_info(mesh->node_group, ref_mesh->node_group,
1347 ref_mesh->n_node_gross) != HECMW_SUCCESS) {
1348 return HECMW_ERROR;
1349 }
1350 if (refine_elem_group_info(mesh, ref_mesh) != HECMW_SUCCESS) {
1351 return HECMW_ERROR;
1352 }
1353 if (refine_surf_group_info(mesh, ref_mesh) != HECMW_SUCCESS) {
1354 return HECMW_ERROR;
1355 }
1356 HECMW_log(HECMW_LOG_DEBUG, "Finished calling refiner.\n");
1357 return HECMW_SUCCESS;
1358}
1359
1360/*============================================================================*/
1361/* */
1362/* rebuild information */
1363/* */
1364/*============================================================================*/
1365static int rebuild_elem_ID(const struct hecmwST_local_mesh *mesh,
1366 struct hecmwST_local_mesh *ref_mesh) {
1367 int i, j, k;
1368 int count;
1369 int *elem_ID_org;
1370 int *elem_ID_ref;
1371
1372 ref_mesh->elem_ID =
1373 (int *)HECMW_malloc(sizeof(int) * ref_mesh->n_elem_gross * 2);
1374 if (ref_mesh->elem_ID == NULL) {
1375 HECMW_set_error(errno, "");
1376 return HECMW_ERROR;
1377 }
1378 count = 0;
1379 elem_ID_org = mesh->elem_ID;
1380 elem_ID_ref = ref_mesh->elem_ID;
1381 for (i = 0; i < mesh->n_elem_type; i++) {
1382 int etype, istart, iend, n_elem, ierror, ndiv;
1383 etype = mesh->elem_type_item[i];
1384 istart = mesh->elem_type_index[i];
1385 iend = mesh->elem_type_index[i + 1];
1386 n_elem = iend - istart;
1387 ndiv = get_elem_ndiv(etype, &ierror);
1388 if (ierror) return HECMW_ERROR;
1389
1390 for (j = 0; j < n_elem; j++) {
1391 int rank = elem_ID_org[1];
1392 for (k = 0; k < ndiv; k++) {
1393 count++;
1394 elem_ID_ref[0] = -count; /* TEMPORARY */
1395 elem_ID_ref[1] = rank; /* to be corrected later */
1396 elem_ID_ref += 2;
1397 }
1398 elem_ID_org += 2;
1399 }
1400 }
1401 return HECMW_SUCCESS;
1402}
1403
1404static int rebuild_global_elem_ID(const struct hecmwST_local_mesh *mesh,
1405 struct hecmwST_local_mesh *ref_mesh) {
1406 int i, j, k;
1407 int min_ID;
1408 int *global_elem_ID_org;
1409 int *global_elem_ID_ref;
1410
1411 ref_mesh->global_elem_ID =
1412 (int *)HECMW_malloc(sizeof(int) * ref_mesh->n_elem_gross);
1413 if (ref_mesh->global_elem_ID == NULL) {
1414 HECMW_set_error(errno, "");
1415 return HECMW_ERROR;
1416 }
1417 min_ID = 0;
1418 global_elem_ID_org = mesh->global_elem_ID;
1419 for (i = 0; i < mesh->n_elem_gross; i++) {
1420 if (global_elem_ID_org[i] < min_ID)
1421 min_ID = global_elem_ID_org[i];
1422 }
1423 global_elem_ID_ref = ref_mesh->global_elem_ID;
1424 for (i = 0; i < mesh->n_elem_type; i++) {
1425 int etype, istart, iend, n_elem, ierror, ndiv;
1426 etype = mesh->elem_type_item[i];
1427 istart = mesh->elem_type_index[i];
1428 iend = mesh->elem_type_index[i + 1];
1429 n_elem = iend - istart;
1430 ndiv = get_elem_ndiv(etype, &ierror);
1431 if (ierror) return HECMW_ERROR;
1432
1433 for (j = 0; j < n_elem; j++) {
1434 global_elem_ID_ref[0] = global_elem_ID_org[0];
1435 for (k = 1; k < ndiv; k++) {
1436 min_ID--;
1437 global_elem_ID_ref[k] = min_ID;
1438 }
1439 global_elem_ID_org += 1;
1440 global_elem_ID_ref += ndiv;
1441 }
1442 }
1443 return HECMW_SUCCESS;
1444}
1445
1446static int rebuild_elem_type(const struct hecmwST_local_mesh *mesh,
1447 struct hecmwST_local_mesh *ref_mesh) {
1448 int i, j, k;
1449 int *elem_type_ref;
1450
1451 ref_mesh->elem_type =
1452 (int *)HECMW_malloc(sizeof(int) * ref_mesh->n_elem_gross);
1453 if (ref_mesh->elem_type == NULL) {
1454 HECMW_set_error(errno, "");
1455 return HECMW_ERROR;
1456 }
1457 elem_type_ref = ref_mesh->elem_type;
1458 for (i = 0; i < mesh->n_elem_type; i++) {
1459 int etype, istart, iend, n_elem, ierror, ndiv;
1460 etype = mesh->elem_type_item[i];
1461 istart = mesh->elem_type_index[i];
1462 iend = mesh->elem_type_index[i + 1];
1463 n_elem = iend - istart;
1464 ndiv = get_elem_ndiv(etype, &ierror);
1465 if (ierror) return HECMW_ERROR;
1466
1467 for (j = 0; j < n_elem; j++) {
1468 for (k = 0; k < ndiv; k++) {
1469 elem_type_ref[k] = etype;
1470 }
1471 elem_type_ref += ndiv;
1472 }
1473 }
1474 return HECMW_SUCCESS;
1475}
1476
1477static int rebuild_section_ID(const struct hecmwST_local_mesh *mesh,
1478 struct hecmwST_local_mesh *ref_mesh) {
1479 int i, j, k;
1480 int *section_ID_org;
1481 int *section_ID_ref;
1482 ref_mesh->section_ID =
1483 (int *)HECMW_malloc(sizeof(int) * ref_mesh->n_elem_gross);
1484 if (ref_mesh->section_ID == NULL) {
1485 HECMW_set_error(errno, "");
1486 return HECMW_ERROR;
1487 }
1488 section_ID_org = mesh->section_ID;
1489 section_ID_ref = ref_mesh->section_ID;
1490 for (i = 0; i < mesh->n_elem_type; i++) {
1491 int etype, istart, iend, n_elem, ierror, ndiv;
1492 etype = mesh->elem_type_item[i];
1493 istart = mesh->elem_type_index[i];
1494 iend = mesh->elem_type_index[i + 1];
1495 n_elem = iend - istart;
1496 ndiv = get_elem_ndiv(etype, &ierror);
1497 if (ierror) return HECMW_ERROR;
1498
1499 for (j = 0; j < n_elem; j++) {
1500 int sect_id = section_ID_org[0];
1501 for (k = 0; k < ndiv; k++) {
1502 section_ID_ref[k] = sect_id;
1503 }
1504 section_ID_org += 1;
1505 section_ID_ref += ndiv;
1506 }
1507 }
1508 return HECMW_SUCCESS;
1509}
1510
1511static int rebuild_elem_mat_ID_index(const struct hecmwST_local_mesh *mesh,
1512 struct hecmwST_local_mesh *ref_mesh) {
1513 int i, j, k;
1514 int *elem_mat_ID_index_org;
1515 int *elem_mat_ID_index_ref;
1516
1517 ref_mesh->elem_mat_ID_index =
1518 (int *)HECMW_malloc(sizeof(int) * (ref_mesh->n_elem_gross + 1));
1519 if (ref_mesh->elem_mat_ID_index == NULL) {
1520 HECMW_set_error(errno, "");
1521 return HECMW_ERROR;
1522 }
1523 ref_mesh->elem_mat_ID_index[0] = 0;
1524 elem_mat_ID_index_org = mesh->elem_mat_ID_index;
1525 elem_mat_ID_index_ref = ref_mesh->elem_mat_ID_index;
1526 for (i = 0; i < mesh->n_elem_type; i++) {
1527 int etype, istart, iend, n_elem, ierror, ndiv;
1528 etype = mesh->elem_type_item[i];
1529 istart = mesh->elem_type_index[i];
1530 iend = mesh->elem_type_index[i + 1];
1531 n_elem = iend - istart;
1532 ndiv = get_elem_ndiv(etype, &ierror);
1533 if (ierror) return HECMW_ERROR;
1534
1535 for (j = 0; j < n_elem; j++) {
1536 int n = elem_mat_ID_index_org[1] - elem_mat_ID_index_org[0];
1537 for (k = 0; k < ndiv; k++) {
1538 elem_mat_ID_index_ref[k + 1] = elem_mat_ID_index_ref[k] + n;
1539 }
1540 elem_mat_ID_index_org += 1;
1541 elem_mat_ID_index_ref += ndiv;
1542 }
1543 }
1544 return HECMW_SUCCESS;
1545}
1546
1547static int rebuild_elem_mat_ID_item(const struct hecmwST_local_mesh *mesh,
1548 struct hecmwST_local_mesh *ref_mesh) {
1549 int i, j, k, l;
1550 int *elem_mat_ID_index_org;
1551 int *elem_mat_ID_item_org;
1552 int *elem_mat_ID_item_ref;
1553
1554 ref_mesh->elem_mat_ID_item = (int *)HECMW_malloc(
1555 sizeof(int) * ref_mesh->elem_mat_ID_index[ref_mesh->n_elem_gross]);
1556 if (ref_mesh->elem_mat_ID_item == NULL) {
1557 HECMW_set_error(errno, "");
1558 return HECMW_ERROR;
1559 }
1560 elem_mat_ID_index_org = mesh->elem_mat_ID_index;
1561 elem_mat_ID_item_org = mesh->elem_mat_ID_item;
1562 elem_mat_ID_item_ref = ref_mesh->elem_mat_ID_item;
1563 for (i = 0; i < mesh->n_elem_type; i++) {
1564 int etype, istart, iend, n_elem, ierror, ndiv;
1565 etype = mesh->elem_type_item[i];
1566 istart = mesh->elem_type_index[i];
1567 iend = mesh->elem_type_index[i + 1];
1568 n_elem = iend - istart;
1569 ndiv = get_elem_ndiv(etype, &ierror);
1570 if (ierror) return HECMW_ERROR;
1571
1572 for (j = 0; j < n_elem; j++) {
1573 int n = elem_mat_ID_index_org[1] - elem_mat_ID_index_org[0];
1574 for (k = 0; k < ndiv; k++) {
1575 for (l = 0; l < n; l++) {
1576 elem_mat_ID_item_ref[l] = elem_mat_ID_item_org[l];
1577 }
1578 elem_mat_ID_item_ref += n;
1579 }
1580 elem_mat_ID_index_org += 1;
1581 elem_mat_ID_item_org += n;
1582 }
1583 }
1584 return HECMW_SUCCESS;
1585}
1586
1587static int rebuild_elem_info(const struct hecmwST_local_mesh *mesh,
1588 struct hecmwST_local_mesh *ref_mesh) {
1589 HECMW_log(HECMW_LOG_DEBUG, "Started rebuilding element info...\n");
1590
1591 if (rebuild_elem_ID(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
1592 if (rebuild_global_elem_ID(mesh, ref_mesh) != HECMW_SUCCESS)
1593 return HECMW_ERROR;
1594 if (rebuild_elem_type(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
1595 if (rebuild_section_ID(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
1596 if (rebuild_elem_mat_ID_index(mesh, ref_mesh) != HECMW_SUCCESS)
1597 return HECMW_ERROR;
1598 if (rebuild_elem_mat_ID_item(mesh, ref_mesh) != HECMW_SUCCESS)
1599 return HECMW_ERROR;
1600
1601 ref_mesh->n_elem_mat_ID = ref_mesh->elem_mat_ID_index[ref_mesh->n_elem_gross];
1602
1603 HECMW_log(HECMW_LOG_DEBUG, "Finished rebuilding element info.\n");
1604 return HECMW_SUCCESS;
1605}
1606
1607static int rebuild_node_info(const struct hecmwST_local_mesh *mesh,
1608 struct hecmwST_local_mesh *ref_mesh) {
1609 int i;
1610 int count, min_ID;
1611
1612 HECMW_log(HECMW_LOG_DEBUG, "Started rebuilding node info...\n");
1613
1614 /* allocate node_ID */
1615 ref_mesh->node_ID =
1616 (int *)HECMW_malloc(sizeof(int) * ref_mesh->n_node_gross * 2);
1617 if (ref_mesh->node_ID == NULL) {
1618 HECMW_set_error(errno, "");
1619 return HECMW_ERROR;
1620 }
1621
1622 /* node_ID */
1623 for (i = 0; i < mesh->n_node_gross; i++) {
1624 ref_mesh->node_ID[2 * i] = mesh->node_ID[2 * i];
1625 ref_mesh->node_ID[2 * i + 1] = mesh->node_ID[2 * i + 1];
1626 }
1627 count = 0;
1628 for (i = mesh->n_node_gross; i < ref_mesh->n_node_gross; i++) {
1629 count++;
1630 ref_mesh->node_ID[2 * i] = -count; /* TEMPORARY */
1631 ref_mesh->node_ID[2 * i + 1] = mesh->my_rank; /* to be corrected later */
1632 }
1633
1634 /* global_node_ID */
1635 ref_mesh->global_node_ID =
1636 (int *)HECMW_malloc(sizeof(int) * ref_mesh->n_node_gross);
1637 if (ref_mesh->global_node_ID == NULL) {
1638 HECMW_set_error(errno, "");
1639 return HECMW_ERROR;
1640 }
1641 min_ID = 0;
1642 for (i = 0; i < mesh->n_node_gross; i++) {
1643 ref_mesh->global_node_ID[i] = mesh->global_node_ID[i];
1644 if (mesh->global_node_ID[i] < min_ID)
1645 min_ID = mesh->global_node_ID[i];
1646 }
1647 for (i = mesh->n_node_gross; i < ref_mesh->n_node_gross; i++) {
1648 min_ID--;
1649 ref_mesh->global_node_ID[i] = min_ID;
1650 }
1651
1653 /* node_init_val_index */
1654 ref_mesh->node_init_val_index =
1655 (int *)HECMW_malloc(sizeof(int) * (ref_mesh->n_node_gross + 1));
1656 if (ref_mesh->node_init_val_index == NULL) {
1657 HECMW_set_error(errno, "");
1658 return HECMW_ERROR;
1659 }
1660 for (i = 0; i < mesh->n_node_gross + 1; i++) {
1661 ref_mesh->node_init_val_index[i] = mesh->node_init_val_index[i];
1662 }
1663 for (i = mesh->n_node_gross + 1; i < ref_mesh->n_node_gross + 1; i++) {
1664 ref_mesh->node_init_val_index[i] =
1666 }
1667
1668 /* node_init_val_item */
1669 if (ref_mesh->node_init_val_index[ref_mesh->n_node_gross] == 0) {
1670 ref_mesh->node_init_val_item = NULL;
1671 } else {
1672 ref_mesh->node_init_val_item = (double *)HECMW_malloc(
1673 sizeof(double) * mesh->node_init_val_index[mesh->n_node_gross]);
1674 if (ref_mesh->node_init_val_item == NULL) {
1675 HECMW_set_error(errno, "");
1676 return HECMW_ERROR;
1677 }
1678 for (i = 0; i < mesh->node_init_val_index[mesh->n_node_gross]; i++) {
1679 ref_mesh->node_init_val_item[i] = mesh->node_init_val_item[i];
1680 }
1681 }
1682 }
1683
1684 HECMW_log(HECMW_LOG_DEBUG, "Finished rebuilding node info.\n");
1685 return HECMW_SUCCESS;
1686}
1687
1688/*
1689 * static functions called by rebuild_comm_tables
1690 */
1691
1692static int get_refined_shared(const struct hecmwST_local_mesh *mesh,
1693 struct hecmwST_local_mesh *ref_mesh,
1694 struct hecmw_varray_int *shared) {
1695 int i;
1696
1697 for (i = 0; i < mesh->n_neighbor_pe; i++) {
1698 int start, n_elem, n_elem_ref, ret;
1699 int *elem_list, *elem_list_ref;
1700 start = mesh->shared_index[i];
1701 n_elem = mesh->shared_index[i + 1] - start;
1702 elem_list = mesh->shared_item + start;
1703 n_elem_ref = get_refined_element_count(mesh, n_elem, elem_list);
1704 if (n_elem_ref < 0) return HECMW_ERROR;
1705 HECMW_varray_int_resize(shared + i, n_elem_ref);
1706 elem_list_ref = HECMW_varray_int_get_v(shared + i);
1707 ret =
1708 get_refined_elem_list(mesh, ref_mesh, n_elem, elem_list, elem_list_ref);
1709 HECMW_assert(ret == n_elem_ref);
1710 }
1711 return HECMW_SUCCESS;
1712}
1713
1714static int elem_type_rcap2hecmw(int etype) {
1715 switch (etype) {
1716 case RCAP_SEGMENT:
1717 return HECMW_ETYPE_ROD1;
1718 case RCAP_SEGMENT2:
1719 return HECMW_ETYPE_ROD2;
1720 case RCAP_TRIANGLE:
1721 return HECMW_ETYPE_TRI1;
1722 case RCAP_TRIANGLE2:
1723 return HECMW_ETYPE_TRI2;
1724 case RCAP_QUAD:
1725 return HECMW_ETYPE_QUA1;
1726 case RCAP_QUAD2:
1727 return HECMW_ETYPE_QUA2;
1728 case RCAP_TETRAHEDRON:
1729 return HECMW_ETYPE_TET1;
1730 case RCAP_TETRAHEDRON2:
1731 return HECMW_ETYPE_TET2;
1732 case RCAP_WEDGE:
1733 return HECMW_ETYPE_PRI1;
1734 case RCAP_WEDGE2:
1735 return HECMW_ETYPE_PRI2;
1736 case RCAP_PYRAMID:
1737 return HECMW_ETYPE_PYR1;
1738 case RCAP_PYRAMID2:
1739 return HECMW_ETYPE_PYR2;
1740 case RCAP_HEXAHEDRON:
1741 return HECMW_ETYPE_HEX1;
1742 case RCAP_HEXAHEDRON2:
1743 return HECMW_ETYPE_HEX2;
1744 }
1745 return HECMW_ERROR;
1746}
1747
1748static int get_rcap_elem_max_node(int etype_rcap) {
1749 int etype_hecmw = elem_type_rcap2hecmw(etype_rcap);
1750 return HECMW_get_max_node(etype_hecmw);
1751}
1752
1753static int determine_node_rank_of_purely_external_elems(
1754 const struct hecmwST_local_mesh *mesh,
1755 struct hecmwST_local_mesh *ref_mesh) {
1756 int eid, rank, nn, k, nid;
1757 int *elem_nodes;
1758
1759 for (eid = 1; eid <= ref_mesh->n_elem_gross; eid++) {
1760 rank = ref_mesh->elem_ID[2 * (eid - 1) + 1];
1761
1762 if (rank >= 0) continue; /* purely external elems have rank < 0 */
1763
1764 nn = ref_mesh->elem_node_index[eid] - ref_mesh->elem_node_index[eid - 1];
1765 elem_nodes = ref_mesh->elem_node_item + ref_mesh->elem_node_index[eid - 1];
1766 for (k = 0; k < nn; k++) {
1767 nid = elem_nodes[k];
1768 if (nid <= mesh->n_node_gross) continue; /* skip old nodes */
1769
1770 /* new nodes that have not been assigned new rank have default rank = my_rank */
1771 if (ref_mesh->node_ID[2 * (nid - 1) + 1] == ref_mesh->my_rank) {
1772 /* assign out-of-range rank = n_subdomain as purely external number = -rank-1 */
1773 ref_mesh->node_ID[2 * (nid - 1) + 1] = - ref_mesh->n_subdomain - 1;
1774 }
1775 }
1776 }
1777 return HECMW_SUCCESS;
1778}
1779
1780static int determine_node_rank(const struct hecmwST_local_mesh *mesh,
1781 const struct hecmw_varray_int *shared,
1782 struct hecmwST_local_mesh *ref_mesh) {
1783 int i, j, k;
1784 struct hecmw_set_int boundary_nodes;
1785 int bnode;
1786
1787 HECMW_log(HECMW_LOG_DEBUG, "Started determine_node_rank...\n");
1788
1789 HECMW_set_int_init(&boundary_nodes);
1790
1791 /* Collect nodes of the shared elements, except for the old nodes */
1792 for (i = 0; i < mesh->n_neighbor_pe; i++) {
1793 for (j = 0; j < HECMW_varray_int_nval(shared + i); j++) {
1794 int eid = HECMW_varray_int_get(shared + i, j);
1795 int nn =
1796 ref_mesh->elem_node_index[eid] - ref_mesh->elem_node_index[eid - 1];
1797 int *elem_nodes =
1798 ref_mesh->elem_node_item + ref_mesh->elem_node_index[eid - 1];
1799
1800 for (k = 0; k < nn; k++) {
1801 if (elem_nodes[k] <= mesh->n_node_gross) continue;
1802
1803 HECMW_set_int_add(&boundary_nodes, elem_nodes[k]);
1804 }
1805 }
1806 }
1807 HECMW_set_int_check_dup(&boundary_nodes);
1808
1809 /* for each node: determine which domain it belongs to; save in node_ID */
1810 HECMW_set_int_iter_init(&boundary_nodes);
1811 while (HECMW_set_int_iter_next(&boundary_nodes, &bnode)) {
1812 int original[HECMW_MAX_NODE_MAX];
1813 int orig_type_rcap = rcapGetOriginal(bnode, original);
1814 int min_rank = mesh->n_subdomain;
1815 int max_rank = -1;
1816 int nn;
1817
1818 HECMW_assert(orig_type_rcap != RCAP_UNKNOWNTYPE);
1819
1820 nn = get_rcap_elem_max_node(orig_type_rcap);
1821
1822 for (k = 0; k < nn; k++) {
1823 int rank;
1824
1825 HECMW_assert(1 <= original[k] && original[k] < bnode);
1826
1827 rank = ref_mesh->node_ID[2 * (original[k] - 1) + 1];
1828
1829 /* for purely external nodes, rank (0, 1, ...) is stored as -rank-1 (-1,
1830 * -2, ...) */
1831 if (rank < 0) rank = -rank - 1;
1832
1833 if (rank < min_rank) min_rank = rank;
1834 if (max_rank < rank) max_rank = rank;
1835 }
1836 HECMW_assert(0 <= min_rank && min_rank < mesh->n_subdomain);
1837 HECMW_assert(0 <= max_rank && max_rank < mesh->n_subdomain);
1838
1839 /* min-rule for the first step, max-rule for the following steps */
1840 if (mesh->n_refine == 0) {
1841 ref_mesh->node_ID[2 * (bnode - 1) + 1] = min_rank;
1842 } else {
1843 ref_mesh->node_ID[2 * (bnode - 1) + 1] = max_rank;
1844 }
1845 }
1846
1847 HECMW_set_int_finalize(&boundary_nodes);
1848
1849 /* new nodes in purely external elems */
1850 if (determine_node_rank_of_purely_external_elems(mesh, ref_mesh) != HECMW_SUCCESS)
1851 return HECMW_ERROR;
1852
1853 HECMW_log(HECMW_LOG_DEBUG, "Finished determine_node_rank.\n");
1854 return HECMW_SUCCESS;
1855}
1856
1857static int *conv_index_ucd2hec(int etype)
1858{
1859 static int conv_index_ucd2hec_rod1[] = {0, 1};
1860 static int conv_index_ucd2hec_rod2[] = {0, -1, 2};
1861 static int conv_index_ucd2hec_tri1[] = {0, 1, 2};
1862 static int conv_index_ucd2hec_tri2[] = {0, 1, 2, -1, -1, -1};
1863 static int conv_index_ucd2hec_qua1[] = {0, 1, 2, 3};
1864 static int conv_index_ucd2hec_qua2[] = {0, 1, 2, 3, -1, -1, -1, -1};
1865 static int conv_index_ucd2hec_tet1[] = {0, 3, 2, 1};
1866 static int conv_index_ucd2hec_tet2[] = {0, 3, 2, 1, -1, -1, -1, -1, -1, -1};
1867 static int conv_index_ucd2hec_pri1[] = {3, 4, 5, 0, 1, 2};
1868 static int conv_index_ucd2hec_pri2[] = {3, 4, 5, 0, 1, 2, -1, -1,
1869 -1, -1, -1, -1, -1, -1, -1};
1870 static int conv_index_ucd2hec_hex1[] = {4, 5, 6, 7, 0, 1, 2, 3};
1871 static int conv_index_ucd2hec_hex2[] = {4, 5, 6, 7, 0, 1, 2, 3, -1, -1,
1872 -1, -1, -1, -1, -1, -1, -1, -1, -1, -1};
1873 static int conv_index_ucd2hec_pyr1[] = {4, 0, 1, 2, 3};
1874 static int conv_index_ucd2hec_pyr2[] = {4, 0, 1, 2, 3, -1, -1,
1875 -1, -1, -1, -1, -1, -1};
1876 static int conv_index_ucd2hec_mst1[] = {0, 1, 2, 3};
1877 static int conv_index_ucd2hec_mst2[] = {0, 1, 2, 3, -1, -1, -1, -1};
1878 static int conv_index_ucd2hec_msq1[] = {0, 1, 2, 3, 4};
1879 static int conv_index_ucd2hec_msq2[] = {0, 1, 2, 3, 4, -1, -1, -1, -1};
1880 static int conv_index_ucd2hec_jtt1[] = {3, 4, 5, 0, 1, 2};
1881 static int conv_index_ucd2hec_jtt2[] = {3, 4, 5, 0, 1, 2,
1882 -1, -1, -1, -1, -1, -1};
1883 static int conv_index_ucd2hec_jtq1[] = {4, 5, 6, 7, 0, 1, 2, 3};
1884 static int conv_index_ucd2hec_jtq2[] = {4, 5, 6, 7, 0, 1, 2, 3,
1885 -1, -1, -1, -1, -1, -1, -1, -1};
1886#if 0
1887 static int conv_index_ucd2hec_bem1[] = {0, 1};
1888 static int conv_index_ucd2hec_bem2[] = {0, -1, 1};
1889#endif
1890 static int conv_index_ucd2hec_sht1[] = {0, 1, 2};
1891 static int conv_index_ucd2hec_sht2[] = {0, 1, 2, -1, -1, -1};
1892 static int conv_index_ucd2hec_shq1[] = {0, 1, 2, 3};
1893 static int conv_index_ucd2hec_shq2[] = {0, 1, 2, 3, -1, -1, -1, -1};
1894 static int conv_index_ucd2hec_ln[] = {0, 1};
1895 switch (etype) {
1896 case HECMW_ETYPE_ROD1:
1897 return conv_index_ucd2hec_rod1;
1898 case HECMW_ETYPE_ROD2:
1899 return conv_index_ucd2hec_rod2;
1900 case HECMW_ETYPE_TRI1:
1901 return conv_index_ucd2hec_tri1;
1902 case HECMW_ETYPE_TRI2:
1903 return conv_index_ucd2hec_tri2;
1904 case HECMW_ETYPE_QUA1:
1905 return conv_index_ucd2hec_qua1;
1906 case HECMW_ETYPE_QUA2:
1907 return conv_index_ucd2hec_qua2;
1908 case HECMW_ETYPE_TET1:
1909 return conv_index_ucd2hec_tet1;
1910 case HECMW_ETYPE_TET2:
1911 return conv_index_ucd2hec_tet2;
1912 case HECMW_ETYPE_PRI1:
1913 return conv_index_ucd2hec_pri1;
1914 case HECMW_ETYPE_PRI2:
1915 return conv_index_ucd2hec_pri2;
1916 case HECMW_ETYPE_HEX1:
1917 return conv_index_ucd2hec_hex1;
1918 case HECMW_ETYPE_HEX2:
1919 return conv_index_ucd2hec_hex2;
1920 case HECMW_ETYPE_PYR1:
1921 return conv_index_ucd2hec_pyr1;
1922 case HECMW_ETYPE_PYR2:
1923 return conv_index_ucd2hec_pyr2;
1924 case HECMW_ETYPE_MST1:
1925 return conv_index_ucd2hec_mst1;
1926 case HECMW_ETYPE_MST2:
1927 return conv_index_ucd2hec_mst2;
1928 case HECMW_ETYPE_MSQ1:
1929 return conv_index_ucd2hec_msq1;
1930 case HECMW_ETYPE_MSQ2:
1931 return conv_index_ucd2hec_msq2;
1932 case HECMW_ETYPE_JTT1:
1933 return conv_index_ucd2hec_jtt1;
1934 case HECMW_ETYPE_JTT2:
1935 return conv_index_ucd2hec_jtt2;
1936 case HECMW_ETYPE_JTQ1:
1937 return conv_index_ucd2hec_jtq1;
1938 case HECMW_ETYPE_JTQ2:
1939 return conv_index_ucd2hec_jtq2;
1940 case HECMW_ETYPE_SHT1:
1941 return conv_index_ucd2hec_sht1;
1942 case HECMW_ETYPE_SHT2:
1943 return conv_index_ucd2hec_sht2;
1944 case HECMW_ETYPE_SHQ1:
1945 return conv_index_ucd2hec_shq1;
1946 case HECMW_ETYPE_SHQ2:
1947 return conv_index_ucd2hec_shq2;
1948 case HECMW_ETYPE_LN11:
1949 case HECMW_ETYPE_LN12:
1950 case HECMW_ETYPE_LN13:
1951 case HECMW_ETYPE_LN14:
1952 case HECMW_ETYPE_LN15:
1953 case HECMW_ETYPE_LN16:
1954 case HECMW_ETYPE_LN21:
1955 case HECMW_ETYPE_LN22:
1956 case HECMW_ETYPE_LN23:
1957 case HECMW_ETYPE_LN24:
1958 case HECMW_ETYPE_LN25:
1959 case HECMW_ETYPE_LN26:
1960 case HECMW_ETYPE_LN31:
1961 case HECMW_ETYPE_LN32:
1962 case HECMW_ETYPE_LN33:
1963 case HECMW_ETYPE_LN34:
1964 case HECMW_ETYPE_LN35:
1965 case HECMW_ETYPE_LN36:
1966 case HECMW_ETYPE_LN41:
1967 case HECMW_ETYPE_LN42:
1968 case HECMW_ETYPE_LN43:
1969 case HECMW_ETYPE_LN44:
1970 case HECMW_ETYPE_LN45:
1971 case HECMW_ETYPE_LN46:
1972 case HECMW_ETYPE_LN51:
1973 case HECMW_ETYPE_LN52:
1974 case HECMW_ETYPE_LN53:
1975 case HECMW_ETYPE_LN54:
1976 case HECMW_ETYPE_LN55:
1977 case HECMW_ETYPE_LN56:
1978 case HECMW_ETYPE_LN61:
1979 case HECMW_ETYPE_LN62:
1980 case HECMW_ETYPE_LN63:
1981 case HECMW_ETYPE_LN64:
1982 case HECMW_ETYPE_LN65:
1983 case HECMW_ETYPE_LN66:
1984 return conv_index_ucd2hec_ln;
1985 default:
1986 return NULL;
1987 }
1988}
1989
1990static void write_elem_with_mismatching_refinement(const struct hecmwST_local_mesh *mesh,
1991 const struct hecmwST_local_mesh *ref_mesh,
1992 const struct hecmw_varray_int *shared,
1993 const struct hecmw_varray_int *shared_nodes,
1994 int i)
1995{
1996 int irank = mesh->neighbor_pe[i];
1997 int j, k, l, nid;
1998 int n_shared_elem = HECMW_varray_int_nval(shared + i);
1999 int n_shared_node = HECMW_varray_int_nval(shared_nodes + i);
2000 struct hecmw_varray_int shared_elem_node_index;
2001 struct hecmw_varray_int shared_elem_node_item;
2002 int n_shared_elem_node;
2003 int n_shared_elem_r;
2004 int n_shared_elem_node_r;
2005 int *shared_elem_node_index_l = NULL;
2006 int *shared_elem_node_item_l = NULL;
2007 int *shared_elem_node_index_r = NULL;
2008 int *shared_elem_node_item_r = NULL;
2009 { /* create elem_node_index/item for shared elements */
2010 if (HECMW_varray_int_init(&shared_elem_node_index) != HECMW_SUCCESS) goto exit;
2011 if (HECMW_varray_int_init(&shared_elem_node_item) != HECMW_SUCCESS) goto exit;
2012 if (HECMW_varray_int_append(&shared_elem_node_index, 0) != HECMW_SUCCESS) goto exit;
2013 n_shared_elem_node = 0;
2014 for (j = 0; j < n_shared_elem; j++) {
2015 int eid = HECMW_varray_int_get(shared + i, j);
2016 int *elem_nodes = ref_mesh->elem_node_item + ref_mesh->elem_node_index[eid-1];
2017 int nn = ref_mesh->elem_node_index[eid] - ref_mesh->elem_node_index[eid-1];
2018 for (k = 0; k < nn; k++) {
2019 nid = elem_nodes[k];
2020 /* find local ID within shared_nodes */
2021 for (l = 0; l < n_shared_node; l++) {
2022 if (HECMW_varray_int_get(shared_nodes + i, l) == nid) break;
2023 }
2024 if (HECMW_varray_int_append(&shared_elem_node_item, l+1) != HECMW_SUCCESS) goto exit;
2025 }
2026 n_shared_elem_node += nn;
2027 if (HECMW_varray_int_append(&shared_elem_node_index, n_shared_elem_node) != HECMW_SUCCESS) goto exit;
2028 }
2029 }
2030 shared_elem_node_index_l = HECMW_varray_int_get_v(&shared_elem_node_index);
2031 shared_elem_node_item_l = HECMW_varray_int_get_v(&shared_elem_node_item);
2032 { /* send/receive shared_elem_node_index/item */
2033 HECMW_Comm comm;
2034 HECMW_Request requests[3];
2035 HECMW_Status statuses[3];
2036 int tag;
2037 { /* send elem_node_index, elem_node_item */
2038 int send_buf[2];
2039 comm = HECMW_comm_get_comm();
2040 send_buf[0] = n_shared_elem;
2041 send_buf[1] = n_shared_elem_node;
2042 tag = 9801;
2043 HECMW_Isend(send_buf, 2, HECMW_INT, irank, tag, comm, &requests[0]);
2044 tag = 9802;
2045 HECMW_Isend(shared_elem_node_index_l, n_shared_elem + 1, HECMW_INT, irank, tag, comm, &requests[1]);
2046 tag = 9803;
2047 HECMW_Isend(shared_elem_node_item_l, n_shared_elem_node, HECMW_INT, irank, tag, comm, &requests[2]);
2048 }
2049 { /* receiv elem_node_index, elem_node_item */
2050 int recv_buf[2];
2051 /* receive n_shared_elem/n_shared_elem_node */
2052 tag = 9801;
2053 HECMW_Recv(recv_buf, 2, HECMW_INT, irank, tag, comm, &statuses[0]);
2054 n_shared_elem_r = recv_buf[0];
2055 n_shared_elem_node_r = recv_buf[1];
2056 /* allocate */
2057 shared_elem_node_index_r = (int *)HECMW_malloc(sizeof(int) * (n_shared_elem_r + 1));
2058 if (shared_elem_node_index_r == NULL) {
2059 HECMW_set_error(errno, "");
2060 goto exit;
2061 }
2062 shared_elem_node_item_r = (int *)HECMW_malloc(sizeof(int) * n_shared_elem_node_r);
2063 if (shared_elem_node_item_r == NULL) {
2064 HECMW_set_error(errno, "");
2065 goto exit;
2066 }
2067 /* receive index/item */
2068 tag = 9802;
2069 HECMW_Recv(shared_elem_node_index_r, n_shared_elem_r + 1, HECMW_INT, irank, tag, comm, &statuses[1]);
2070 tag = 9803;
2071 HECMW_Recv(shared_elem_node_item_r, n_shared_elem_node_r, HECMW_INT, irank, tag, comm, &statuses[2]);
2072 }
2073 HECMW_Waitall(3, requests, statuses);
2074 }
2075 { /* compare */
2076 if (n_shared_elem != n_shared_elem_r) {
2077 HECMW_log(HECMW_LOG_ERROR, "rank %d has different number of shared elements\n", irank);
2078 goto exit;
2079 }
2080 for (j = 0; j < n_shared_elem; j++) {
2081 int is_error = 0;
2082 if (shared_elem_node_index_l[j+1] != shared_elem_node_index_r[j+1]) {
2083 HECMW_log(HECMW_LOG_ERROR, "rank %d has different number of nodes for %dth shared element\n", irank, j+1);
2084 goto exit;
2085 }
2086 for (k = shared_elem_node_index_l[j]; k < shared_elem_node_index_l[j+1]; k++) {
2087 if (shared_elem_node_item_l[k] != shared_elem_node_item_r[k]) {
2088 HECMW_log(HECMW_LOG_ERROR, "rank %d has different node for %dth node of %dth shared element\n", irank,
2089 k - shared_elem_node_index_l[j] + 1, j+1);
2090 is_error = 1;
2091 break;
2092 }
2093 }
2094 if (is_error) {
2095 int eid = HECMW_varray_int_get(shared + i, j);
2096 int *elem_nodes = ref_mesh->elem_node_item + ref_mesh->elem_node_index[eid-1];
2097 int nn = ref_mesh->elem_node_index[eid] - ref_mesh->elem_node_index[eid-1];
2098 HECMW_log(HECMW_LOG_ERROR, "elem local ID: %d\n", eid);
2099 for (k = 0; k < nn; k++) {
2100 nid = elem_nodes[k];
2101 HECMW_log(HECMW_LOG_ERROR, "elem node (%2d) local ID: %10d, coordinates: %14.6e, %14.6e, %14.6e\n", k+1, nid,
2102 ref_mesh->node[3*(nid-1)], ref_mesh->node[3*(nid-1)+1], ref_mesh->node[3*(nid-1)+2]);
2103 }
2104 /* TODO: original element */
2105 /* TODO: origin of the node */
2106 goto exit;
2107 }
2108 }
2109 }
2110
2111exit:
2112 if (shared_elem_node_index_r != NULL) HECMW_free(shared_elem_node_index_r);
2113 if (shared_elem_node_item_r != NULL) HECMW_free(shared_elem_node_item_r);
2114 HECMW_varray_int_finalize(&shared_elem_node_index);
2115 HECMW_varray_int_finalize(&shared_elem_node_item);
2116}
2117
2118static void write_refined_shared_in_ucd(const struct hecmwST_local_mesh *mesh,
2119 const struct hecmwST_local_mesh *ref_mesh,
2120 const struct hecmw_varray_int *shared,
2121 const struct hecmw_varray_int *shared_nodes,
2122 int i)
2123{
2124 char fname[64];
2125 FILE *fp;
2126 int irank = mesh->neighbor_pe[i];
2127 int j, k, l, nid;
2128 int n_shared_elem = HECMW_varray_int_nval(shared + i);
2129 int n_shared_node = HECMW_varray_int_nval(shared_nodes + i);
2130 sprintf(fname, "shared_%d_%d.inp", mesh->my_rank, irank);
2131 HECMW_log(HECMW_LOG_DEBUG, "writing refined shared elements to %s\n", fname);
2132 fp = fopen(fname, "w");
2133 if (fp == NULL) {
2134 HECMW_log(HECMW_LOG_ERROR, "failed to open file %s\n", fname);
2135 } else {
2136 fprintf(fp, "%d %d 1 1 0\n", n_shared_node, n_shared_elem);
2137 for (j = 0; j < n_shared_node; j++) {
2138 nid = HECMW_varray_int_get(shared_nodes + i, j);
2139 fprintf(fp, "%d %e %e %e\n", j+1, ref_mesh->node[3*(nid-1)], ref_mesh->node[3*(nid-1)+1], ref_mesh->node[3*(nid-1)+2]);
2140 }
2141 for (j = 0; j < n_shared_elem; j++) {
2142 int eid = HECMW_varray_int_get(shared + i, j);
2143 int etype = ref_mesh->elem_type[eid-1];
2144 int *elem_nodes = ref_mesh->elem_node_item + ref_mesh->elem_node_index[eid-1];
2145 int nn = ref_mesh->elem_node_index[eid] - ref_mesh->elem_node_index[eid-1];
2146 int *index_ucd = conv_index_ucd2hec(etype);
2147 HECMW_assert(index_ucd);
2148 fprintf(fp, "%d 0 %s", j+1, HECMW_get_ucd_label(etype));
2149 for (k = 0; k < nn; k++) {
2150 if (index_ucd[k] < 0) continue;
2151 nid = elem_nodes[index_ucd[k]];
2152 for (l = 0; l < n_shared_node; l++) {
2153 if (HECMW_varray_int_get(shared_nodes + i, l) == nid) break;
2154 }
2155 fprintf(fp, " %d", l+1);
2156 }
2157 fprintf(fp, "\n");
2158 }
2159 fclose(fp);
2160 }
2161}
2162
2163static int check_node_rank(const struct hecmwST_local_mesh *mesh,
2164 const struct hecmw_varray_int *shared,
2165 const struct hecmwST_local_mesh *ref_mesh)
2166{
2167 int *n_shared_nodes;
2168 struct hecmw_varray_int *shared_nodes, *shared_ranks;
2169 HECMW_Request *requests;
2170 HECMW_Status *statuses;
2171 HECMW_Comm comm;
2172 int *send_buf;
2173 int recv_buf[4];
2174 int n_error = 0, n_error_g;
2175 int i, j, k, irank, tag, nid, rank;
2176
2177 HECMW_log(HECMW_LOG_DEBUG, "Started check_node_rank...\n");
2178
2179 /* alloc memory */
2180 n_shared_nodes = (int *)HECMW_malloc(sizeof(int) * mesh->n_neighbor_pe);
2181 if (n_shared_nodes == NULL) {
2182 HECMW_set_error(errno, "");
2183 return HECMW_ERROR;
2184 }
2185
2186 shared_nodes = (struct hecmw_varray_int *)HECMW_malloc(
2187 sizeof(struct hecmw_varray_int) * mesh->n_neighbor_pe);
2188 if (shared_nodes == NULL) {
2189 HECMW_set_error(errno, "");
2190 return HECMW_ERROR;
2191 }
2192
2193 shared_ranks = (struct hecmw_varray_int *)HECMW_malloc(
2194 sizeof(struct hecmw_varray_int) * mesh->n_neighbor_pe);
2195 if (shared_ranks == NULL) {
2196 HECMW_set_error(errno, "");
2197 return HECMW_ERROR;
2198 }
2199
2200 send_buf = (int *)HECMW_malloc(sizeof(int) * 4 * mesh->n_neighbor_pe);
2201 if (send_buf == NULL) {
2202 HECMW_set_error(errno, "");
2203 return HECMW_ERROR;
2204 }
2205
2206 comm = HECMW_comm_get_comm();
2207
2208 requests = (HECMW_Request *)HECMW_malloc(sizeof(HECMW_Request) * mesh->n_neighbor_pe);
2209 if (requests == NULL) {
2210 HECMW_set_error(errno, "");
2211 return HECMW_ERROR;
2212 }
2213
2214 statuses = (HECMW_Status *)HECMW_malloc(sizeof(HECMW_Status) * mesh->n_neighbor_pe);
2215 if (statuses == NULL) {
2216 HECMW_set_error(errno, "");
2217 return HECMW_ERROR;
2218 }
2219
2220 /*
2221 * check number of refined shared nodes
2222 */
2223
2224 /* Collect nodes of the shared elements, and send size */
2225 for (i = 0; i < mesh->n_neighbor_pe; i++) {
2226 /* numver of shared elements in original mesh */
2227 send_buf[4*i] = mesh->shared_index[i+1] - mesh->shared_index[i];
2228
2229 /* number of refined shared elements */
2230 send_buf[4*i+1] = HECMW_varray_int_nval(shared + i);
2231
2232 HECMW_varray_int_init(shared_nodes + i);
2233 for (j = 0; j < HECMW_varray_int_nval(shared + i); j++) {
2234 int eid = HECMW_varray_int_get(shared + i, j);
2235 int nn = ref_mesh->elem_node_index[eid] - ref_mesh->elem_node_index[eid - 1];
2236 int *elem_nodes = ref_mesh->elem_node_item + ref_mesh->elem_node_index[eid - 1];
2237
2238 for (k = 0; k < nn; k++) {
2239#if 0 /* keep old nodes for plotting */
2240 if (elem_nodes[k] <= mesh->n_node_gross) continue; /* skip old nodes */
2241#endif
2242
2243 HECMW_varray_int_append(shared_nodes + i, elem_nodes[k]);
2244 }
2245 }
2246 /* number of shared nodes before removing duplication */
2247 send_buf[4*i+2] = HECMW_varray_int_nval(shared_nodes + i);
2248
2249#if 1
2250 HECMW_varray_int_rmdup(shared_nodes + i);
2251#endif
2252 n_shared_nodes[i] = HECMW_varray_int_nval(shared_nodes + i);
2253
2254 /* number of shared nodes after removing duplication */
2255 send_buf[4*i+3] = n_shared_nodes[i];
2256
2257 /* send sizes */
2258 irank = mesh->neighbor_pe[i];
2259 tag = 9901;
2260 HECMW_Isend(send_buf + 4*i, 4, HECMW_INT, irank, tag, comm, requests + i);
2261 }
2262
2263 /* recieve size and check */
2264 for (i = 0; i < mesh->n_neighbor_pe; i++) {
2265 irank = mesh->neighbor_pe[i];
2266 tag = 9901;
2267 HECMW_Recv(recv_buf, 4, HECMW_INT, irank, tag, comm, statuses + i);
2268
2269 if (recv_buf[0] != send_buf[4*i]) {
2271 "number of shared elements %d not match with %d in rank %d\n",
2272 send_buf[4*i], recv_buf[0], irank);
2273 n_error++;
2274 }
2275 if (recv_buf[1] != send_buf[4*i+1]) {
2277 "number of refined shared elements %d not match with %d in rank %d\n",
2278 send_buf[4*i+1], recv_buf[1], irank);
2279 }
2280 if (recv_buf[2] != send_buf[4*i+2]) {
2282 "number of new shared nodes before rmdup %d not match with %d in rank %d\n",
2283 send_buf[4*i+2], recv_buf[2], irank);
2284 n_error++;
2285 }
2286 if (recv_buf[3] != send_buf[4*i+3]) {
2288 "number of new shared nodes %d not match with %d in rank %d\n",
2289 send_buf[4*i+3], recv_buf[3], irank);
2290 n_error++;
2291 write_elem_with_mismatching_refinement(mesh, ref_mesh, shared, shared_nodes, i);
2292 write_refined_shared_in_ucd(mesh, ref_mesh, shared, shared_nodes, i);
2293 }
2294 }
2295 HECMW_Waitall(mesh->n_neighbor_pe, requests, statuses);
2296
2297 HECMW_Allreduce(&n_error, &n_error_g, 1, HECMW_INT, HECMW_SUM, comm);
2298 if (n_error_g > 0) return HECMW_ERROR;
2299
2300 /*
2301 * check rank of new shared nodes
2302 */
2303
2304 n_error = 0;
2305
2306 /* send rank of new shared nodes */
2307 for (i = 0; i < mesh->n_neighbor_pe; i++) {
2308 HECMW_varray_int_init(shared_ranks + i);
2309 /* HECMW_varray_int_resize(shared_ranks + i, n_shared_nodes[i]); */
2310 for (j = 0; j < n_shared_nodes[i]; j++) {
2311 nid = HECMW_varray_int_get(shared_nodes + i, j);
2312 rank = ref_mesh->node_ID[2 * (nid - 1) + 1];
2313 /* HECMW_varray_int_set(shared_ranks + i, j, rank); */
2314 HECMW_varray_int_append(shared_ranks + i, rank);
2315 }
2316
2317 irank = mesh->neighbor_pe[i];
2318 tag = 9902;
2319 int *send_ranks = HECMW_varray_int_get_v(shared_ranks + i);
2320 HECMW_Isend(send_ranks, n_shared_nodes[i], HECMW_INT, irank, tag, comm, requests + i);
2321 }
2322
2323 /* recieve rank of new shared nodes */
2324 for (i = 0; i < mesh->n_neighbor_pe; i++) {
2325 int *recv_ranks = (int *)HECMW_malloc(sizeof(int) * n_shared_nodes[i]);
2326 if (recv_ranks == NULL) {
2327 HECMW_set_error(errno, "");
2328 return HECMW_ERROR;
2329 }
2330
2331 irank = mesh->neighbor_pe[i];
2332 tag = 9902;
2333 HECMW_Recv(recv_ranks, n_shared_nodes[i], HECMW_INT, irank, tag, comm, statuses + i);
2334
2335 for (j = 0; j < n_shared_nodes[i]; j++) {
2336 if (recv_ranks[j] != HECMW_varray_int_get(shared_ranks + i, j)) {
2338 "rank of new shared node %d: %d not match with %d in rank %d\n",
2339 HECMW_varray_int_get(shared_nodes + i, j),
2340 HECMW_varray_int_get(shared_ranks + i, j), recv_ranks[j], irank);
2341 n_error++;
2342 }
2343 }
2344 HECMW_free(recv_ranks);
2345 }
2346 HECMW_Waitall(mesh->n_neighbor_pe, requests, statuses);
2347
2348 HECMW_Allreduce(&n_error, &n_error_g, 1, HECMW_INT, HECMW_SUM, comm);
2349 if (n_error_g > 0) return HECMW_ERROR;
2350
2351 /* free memory */
2352 for (i = 0; i < mesh->n_neighbor_pe; i++) {
2353 HECMW_varray_int_finalize(shared_nodes + i);
2354 HECMW_varray_int_finalize(shared_ranks + i);
2355 }
2356 HECMW_free(shared_nodes);
2357 HECMW_free(shared_ranks);
2358 HECMW_free(n_shared_nodes);
2359 HECMW_free(send_buf);
2360 HECMW_free(requests);
2361 HECMW_free(statuses);
2362
2363 HECMW_log(HECMW_LOG_DEBUG, "Finished check_node_rank.\n");
2364
2365 return HECMW_SUCCESS;
2366}
2367
2368static int count_internal(int my_rank, int n, int *ID) {
2369 int i, count = 0;
2370 for (i = 0; i < n; i++) {
2371 if (ID[2 * i + 1] == my_rank) {
2372 ID[2 * i] = ++count;
2373 }
2374 }
2375 return count;
2376}
2377
2378static int count_purely_external(int n, int *ID) {
2379 int i, count = 0;
2380 for (i = 0; i < n; i++) {
2381 if (ID[2 * i + 1] < 0) {
2382 count++;
2383 }
2384 }
2385 return count;
2386}
2387
2388static int *new_internal_list(int my_rank, int n, int n_internal, int *ID) {
2389 int i, k;
2390 int *internal_list;
2391
2392 internal_list = (int *)HECMW_malloc(sizeof(int) * n_internal);
2393 if (internal_list == NULL) {
2394 HECMW_set_error(errno, "");
2395 return NULL;
2396 }
2397
2398 for (i = 0, k = 0; i < n; i++) {
2399 if (ID[2 * i + 1] == my_rank) {
2400 internal_list[k] = i + 1;
2401 k++;
2402 }
2403 }
2404 return internal_list;
2405}
2406
2407static int create_index_item(int n, const struct hecmw_varray_int *arrays,
2408 int **index, int **item) {
2409 int i, j, k;
2410
2411 /* index */
2412 *index = (int *)HECMW_malloc(sizeof(int) * (n + 1));
2413 if (*index == NULL) {
2414 HECMW_set_error(errno, "");
2415 return HECMW_ERROR;
2416 }
2417 (*index)[0] = 0;
2418 for (i = 0; i < n; i++) {
2419 (*index)[i + 1] = (*index)[i] + HECMW_varray_int_nval(arrays + i);
2420 if ((*index)[i + 1] < (*index)[i]) {
2422 "Integer overflow detected while creating index array\n");
2423 return HECMW_ERROR;
2424 }
2425 }
2426
2427 /* item */
2428 if ((*index)[n] == 0) {
2429 (*item) = NULL;
2430 return HECMW_SUCCESS;
2431 }
2432 *item = (int *)HECMW_malloc(sizeof(int) * (*index)[n]);
2433 if (*item == NULL) {
2434 HECMW_set_error(errno, "");
2435 return HECMW_ERROR;
2436 }
2437 k = 0;
2438 for (i = 0; i < n; i++) {
2439 for (j = 0; j < HECMW_varray_int_nval(arrays + i); j++) {
2440 (*item)[k] = HECMW_varray_int_get(arrays + i, j);
2441 k++;
2442 }
2443 }
2444 return HECMW_SUCCESS;
2445}
2446
2447static int new_shared_export_import(const struct hecmwST_local_mesh *mesh,
2448 const struct hecmw_varray_int *shared,
2449 struct hecmwST_local_mesh *ref_mesh) {
2450 struct hecmw_varray_int *new_shared, *new_export, *new_import;
2451 int i;
2452
2453 /* prepare new shared, import and export lists */
2454 new_shared = (struct hecmw_varray_int *)HECMW_malloc(
2455 sizeof(struct hecmw_varray_int) * mesh->n_neighbor_pe);
2456 if (new_shared == NULL) {
2457 HECMW_set_error(errno, "");
2458 return HECMW_ERROR;
2459 }
2460 new_export = (struct hecmw_varray_int *)HECMW_malloc(
2461 sizeof(struct hecmw_varray_int) * mesh->n_neighbor_pe);
2462 if (new_export == NULL) {
2463 HECMW_set_error(errno, "");
2464 return HECMW_ERROR;
2465 }
2466 new_import = (struct hecmw_varray_int *)HECMW_malloc(
2467 sizeof(struct hecmw_varray_int) * mesh->n_neighbor_pe);
2468 if (new_import == NULL) {
2469 HECMW_set_error(errno, "");
2470 return HECMW_ERROR;
2471 }
2472
2473 /* for each refined shared list */
2474 for (i = 0; i < mesh->n_neighbor_pe; i++) {
2475 int j;
2476 int rank_i = mesh->neighbor_pe[i];
2477
2478 HECMW_varray_int_init(new_shared + i);
2479 HECMW_varray_int_init(new_export + i);
2480 HECMW_varray_int_init(new_import + i);
2481
2482 /* for each element in the list */
2483 for (j = 0; j < HECMW_varray_int_nval(shared + i); j++) {
2484 struct hecmw_varray_int exp_cand, imp_cand;
2485
2486 int eid = HECMW_varray_int_get(shared + i, j);
2487 int nn =
2488 ref_mesh->elem_node_index[eid] - ref_mesh->elem_node_index[eid - 1];
2489 int *elem_nodes =
2490 ref_mesh->elem_node_item + ref_mesh->elem_node_index[eid - 1];
2491 int etype = ref_mesh->elem_type[eid - 1];
2492
2493 int min_rank = mesh->n_subdomain;
2494 int myrank_included = 0;
2495 int rank_i_included = 0;
2496 int k;
2497
2498 HECMW_varray_int_init(&exp_cand);
2499 HECMW_varray_int_init(&imp_cand);
2500
2501 /* for each node in the element */
2502 for (k = 0; k < nn; k++) {
2503 int rank = ref_mesh->node_ID[2 * (elem_nodes[k] - 1) + 1];
2504
2505 /* for purely external nodes, rank (0, 1, ...) is stored as -rank-1 (-1,
2506 * -2, ...) */
2507 if (rank < 0) rank = -rank - 1;
2508
2509 if (rank < min_rank) min_rank = rank;
2510
2511 if (rank == mesh->my_rank) {
2512 myrank_included = 1;
2513 HECMW_varray_int_append(&exp_cand, elem_nodes[k]);
2514 } else if (rank == rank_i) {
2515 rank_i_included = 1;
2516 HECMW_varray_int_append(&imp_cand, elem_nodes[k]);
2517 }
2518 }
2519 HECMW_assert(0 <= min_rank && min_rank < mesh->n_subdomain);
2520
2521 ref_mesh->elem_ID[2 * (eid - 1) + 1] = min_rank;
2522
2523 if (HECMW_is_etype_patch(etype)) {
2524 if (myrank_included || rank_i_included) HECMW_varray_int_append(new_shared + i, eid);
2525 if (myrank_included) HECMW_varray_int_cat(new_export + i, &exp_cand);
2526 if (rank_i_included) HECMW_varray_int_cat(new_import + i, &imp_cand);
2527 } else {
2528 if (myrank_included && rank_i_included) { /* the element is shared */
2529 HECMW_varray_int_append(new_shared + i, eid);
2530 HECMW_varray_int_cat(new_export + i, &exp_cand);
2531 HECMW_varray_int_cat(new_import + i, &imp_cand);
2532 } else if (myrank_included == 0) {
2533 /* for purely external elements, rank (0, 1, ...) is stored as -rank-1
2534 * (-1, -2, ...) */
2535 ref_mesh->elem_ID[2 * (eid - 1) + 1] = -min_rank - 1;
2536 }
2537 }
2538
2539 HECMW_varray_int_finalize(&exp_cand);
2540 HECMW_varray_int_finalize(&imp_cand);
2541 }
2542 /* remove duplication */
2543 HECMW_varray_int_rmdup(new_export + i);
2544 HECMW_varray_int_rmdup(new_import + i);
2545 }
2546
2547 /* new shared */
2548 if (create_index_item(mesh->n_neighbor_pe, new_shared, &(ref_mesh->shared_index),
2549 &(ref_mesh->shared_item)) != HECMW_SUCCESS) {
2550 HECMW_log(HECMW_LOG_ERROR, "Create shared_index and shared_item failed\n");
2551 return HECMW_ERROR;
2552 }
2553 HECMW_log(HECMW_LOG_DEBUG, "Total number of shared elements = %d\n",
2554 ref_mesh->shared_index[mesh->n_neighbor_pe]);
2555
2556 /* new export */
2557 if (create_index_item(mesh->n_neighbor_pe, new_export, &(ref_mesh->export_index),
2558 &(ref_mesh->export_item)) != HECMW_SUCCESS) {
2559 HECMW_log(HECMW_LOG_ERROR, "Create export_index and export_item failed\n");
2560 return HECMW_ERROR;
2561 }
2562 HECMW_log(HECMW_LOG_DEBUG, "Total number of export nodes = %d\n",
2563 ref_mesh->export_index[mesh->n_neighbor_pe]);
2564
2565 /* new import */
2566 if (create_index_item(mesh->n_neighbor_pe, new_import, &(ref_mesh->import_index),
2567 &(ref_mesh->import_item)) != HECMW_SUCCESS) {
2568 HECMW_log(HECMW_LOG_ERROR, "Create import_index and import_item failed\n");
2569 return HECMW_ERROR;
2570 }
2571 HECMW_log(HECMW_LOG_DEBUG, "Total number of import nodes = %d\n",
2572 ref_mesh->import_index[mesh->n_neighbor_pe]);
2573
2574 /* deallocate new shared, import and export lists */
2575 for (i = 0; i < mesh->n_neighbor_pe; i++) {
2576 HECMW_varray_int_finalize(new_shared + i);
2577 HECMW_varray_int_finalize(new_export + i);
2578 HECMW_varray_int_finalize(new_import + i);
2579 }
2580 HECMW_free(new_shared);
2581 HECMW_free(new_export);
2582 HECMW_free(new_import);
2583
2584 return HECMW_SUCCESS;
2585}
2586
2587static int rebuild_comm_tables_serial(const struct hecmwST_local_mesh *mesh,
2588 struct hecmwST_local_mesh *ref_mesh) {
2589 int i;
2590
2593
2594 /* nn_internal, n_node */
2595 ref_mesh->nn_internal = count_internal(mesh->my_rank, ref_mesh->n_node_gross, ref_mesh->node_ID);
2596 ref_mesh->n_node = ref_mesh->n_node_gross;
2597 ref_mesh->nn_middle = ref_mesh->n_node;
2598
2599 HECMW_log(HECMW_LOG_DEBUG, "nn_internal = %d, n_node = %d, nn_middle = %d\n",
2600 ref_mesh->nn_internal, ref_mesh->n_node, ref_mesh->nn_middle);
2601
2602 /* node_internal_list */
2603 ref_mesh->node_internal_list =
2604 (int *)HECMW_malloc(sizeof(int) * ref_mesh->nn_internal);
2605 if (ref_mesh->node_internal_list == NULL) {
2606 HECMW_set_error(errno, "");
2607 return HECMW_ERROR;
2608 }
2609 for (i = 0; i < ref_mesh->nn_internal; i++) {
2610 ref_mesh->node_internal_list[i] = i + 1;
2611 }
2612
2613 /* ne_internal, n_elem */
2614 ref_mesh->ne_internal = count_internal(mesh->my_rank, ref_mesh->n_elem_gross, ref_mesh->elem_ID);
2615 ref_mesh->n_elem = ref_mesh->n_elem_gross;
2616
2617 HECMW_log(HECMW_LOG_DEBUG, "ne_internal = %d, n_elem = %d\n",
2618 ref_mesh->ne_internal, ref_mesh->n_elem);
2619
2620 /* elem_internal_list */
2621 ref_mesh->elem_internal_list =
2622 (int *)HECMW_malloc(sizeof(int) * ref_mesh->ne_internal);
2623 if (ref_mesh->elem_internal_list == NULL) {
2624 HECMW_set_error(errno, "");
2625 return HECMW_ERROR;
2626 }
2627 for (i = 0; i < ref_mesh->ne_internal; i++) {
2628 ref_mesh->elem_internal_list[i] = i + 1;
2629 }
2630 return HECMW_SUCCESS;
2631}
2632
2633static int mark_purely_external_nodes(int my_rank,
2634 struct hecmwST_local_mesh *ref_mesh) {
2635 struct hecmw_bit_array mark;
2636 int i;
2637
2638 if (HECMW_bit_array_init(&mark, ref_mesh->n_node_gross) != HECMW_SUCCESS) {
2639 return HECMW_ERROR;
2640 }
2641 /* mark all nodes in the new import lists */
2642 for (i = 0; i < ref_mesh->import_index[ref_mesh->n_neighbor_pe]; i++) {
2643 HECMW_bit_array_set(&mark, ref_mesh->import_item[i] - 1);
2644 }
2645 /* scan for all external, unmarked nodes and set their rank to -rank-1 */
2646 for (i = 0; i < ref_mesh->n_node_gross; i++) {
2647 int rank = ref_mesh->node_ID[2 * i + 1];
2648 if (rank == my_rank || rank < 0) continue;
2649 if (!HECMW_bit_array_get(&mark, i)) {
2650 /* for purely external nodes, rank (0, 1, ...) is stored as -rank-1 (-1,
2651 * -2, ...) */
2652 ref_mesh->node_ID[2 * i + 1] = -rank - 1;
2653 }
2654 }
2656 return HECMW_SUCCESS;
2657}
2658
2659static int rebuild_comm_tables(const struct hecmwST_local_mesh *mesh,
2660 struct hecmwST_local_mesh *ref_mesh) {
2661 int i;
2662 struct hecmw_varray_int *shared;
2663
2664 if (mesh->n_neighbor_pe == 0) { /* SERIAL */
2665 return rebuild_comm_tables_serial(mesh, ref_mesh);
2666 }
2667
2668 /* PARALLEL */
2669 HECMW_log(HECMW_LOG_DEBUG, "Started rebuilding communication tables...\n");
2670
2671 /* check part-type */
2674 "Partitioning-type must be node-based for refinement.\n");
2675 return HECMW_ERROR;
2676 }
2677
2678 /* Get refined shared element list */
2679 shared = (struct hecmw_varray_int *)HECMW_malloc(
2680 sizeof(struct hecmw_varray_int) * mesh->n_neighbor_pe);
2681 if (shared == NULL) {
2682 HECMW_set_error(errno, "");
2683 return HECMW_ERROR;
2684 }
2685 for (i = 0; i < mesh->n_neighbor_pe; i++) {
2686 HECMW_varray_int_init(shared + i);
2687 }
2688 if (get_refined_shared(mesh, ref_mesh, shared) != HECMW_SUCCESS) {
2689 return HECMW_ERROR;
2690 }
2691
2692 /* Determine rank of new nodes */
2693 if (determine_node_rank(mesh, shared, ref_mesh) != HECMW_SUCCESS) {
2694 return HECMW_ERROR;
2695 }
2696
2697 /* Check rank of new nodes with neighbor processes */
2698 if (check_node_rank(mesh, shared, ref_mesh) != HECMW_SUCCESS) {
2699 return HECMW_ERROR;
2700 }
2701
2702 /* nn_internal */
2703 ref_mesh->nn_internal =
2704 count_internal(mesh->my_rank, ref_mesh->n_node_gross, ref_mesh->node_ID);
2705
2707 NULL); /* because PARTTYPE is NODEBASED. */
2708
2709 /* new shared, import and export lists */
2710 if (new_shared_export_import(mesh, shared, ref_mesh) != HECMW_SUCCESS) {
2711 return HECMW_ERROR;
2712 }
2713
2714 /* deallocate refined shared list */
2715 for (i = 0; i < mesh->n_neighbor_pe; i++) {
2716 HECMW_varray_int_finalize(shared + i);
2717 }
2718 HECMW_free(shared);
2719
2720 /* ne_internal */
2721 ref_mesh->ne_internal =
2722 count_internal(mesh->my_rank, ref_mesh->n_elem_gross, ref_mesh->elem_ID);
2723
2724 /* n_elem */
2725 ref_mesh->n_elem =
2726 ref_mesh->n_elem_gross -
2727 count_purely_external(ref_mesh->n_elem_gross, ref_mesh->elem_ID);
2728
2729 /* elem_internal_list */
2730 if (ref_mesh->ne_internal > 0) {
2731 ref_mesh->elem_internal_list =
2732 new_internal_list(mesh->my_rank, ref_mesh->n_elem_gross,
2733 ref_mesh->ne_internal, ref_mesh->elem_ID);
2734 if (ref_mesh->elem_internal_list == NULL) return HECMW_ERROR;
2735 }
2736
2737 /* mark purely external nodes */
2738 if (mark_purely_external_nodes(mesh->my_rank, ref_mesh) != HECMW_SUCCESS) {
2739 return HECMW_ERROR;
2740 }
2741
2742 /* n_node */
2743 ref_mesh->n_node =
2744 ref_mesh->n_node_gross -
2745 count_purely_external(ref_mesh->n_node_gross, ref_mesh->node_ID);
2746 ref_mesh->nn_middle = ref_mesh->n_node; /* not set properly (since it's never used) */
2747
2748 HECMW_log(HECMW_LOG_DEBUG, "nn_internal = %d, n_node = %d nn_middle = %d, \n",
2749 ref_mesh->nn_internal, ref_mesh->n_node, ref_mesh->nn_middle);
2750 HECMW_log(HECMW_LOG_DEBUG, "ne_internal = %d, nelem = %d\n",
2751 ref_mesh->ne_internal, ref_mesh->n_elem);
2752
2753 HECMW_log(HECMW_LOG_DEBUG, "Finished rebuilding communication tables.\n");
2754 return HECMW_SUCCESS;
2755}
2756
2757static int check_comm_table_len(const struct hecmwST_local_mesh *ref_mesh) {
2758 int i, irank, js, je, len, tag, nsend;
2759 HECMW_Request *requests;
2760 HECMW_Status *statuses;
2761 HECMW_Comm comm;
2762 int *send_buf;
2763 int n_error = 0, n_error_g;
2764
2765 if (ref_mesh->n_neighbor_pe == 0) return HECMW_SUCCESS;
2766
2767 HECMW_log(HECMW_LOG_DEBUG, "Started checking communication tables...\n");
2768
2769 comm = HECMW_comm_get_comm();
2770
2771 send_buf = (int *)HECMW_malloc(sizeof(int) * ref_mesh->n_neighbor_pe);
2772 if (send_buf == 0) {
2773 HECMW_set_error(errno, "");
2774 HECMW_abort(comm);
2775 }
2776
2777 requests = (HECMW_Request *)HECMW_malloc(sizeof(HECMW_Request) *
2778 ref_mesh->n_neighbor_pe);
2779 if (requests == NULL) {
2780 HECMW_set_error(errno, "");
2781 HECMW_abort(comm);
2782 }
2783 statuses = (HECMW_Status *)HECMW_malloc(sizeof(HECMW_Status) *
2784 ref_mesh->n_neighbor_pe);
2785 if (statuses == NULL) {
2786 HECMW_set_error(errno, "");
2787 HECMW_abort(comm);
2788 }
2789
2790 /* export and import */
2791
2792 for (i = 0; i < ref_mesh->n_neighbor_pe; i++) {
2793 irank = ref_mesh->neighbor_pe[i];
2794 js = ref_mesh->export_index[i];
2795 je = ref_mesh->export_index[i + 1];
2796 send_buf[i] = je - js;
2797 /* isend number of export nodes */
2798 tag = 0;
2799 HECMW_Isend(send_buf + i, 1, HECMW_INT, irank, tag, comm, requests + i);
2800 }
2801 for (i = 0; i < ref_mesh->n_neighbor_pe; i++) {
2802 irank = ref_mesh->neighbor_pe[i];
2803 js = ref_mesh->import_index[i];
2804 je = ref_mesh->import_index[i + 1];
2805 /* recv number of import nodes */
2806 tag = 0;
2807 HECMW_Recv(&len, 1, HECMW_INT, irank, tag, comm, statuses + i);
2808 if (len != je - js) {
2810 "inconsistent length of import (%d) with export on rank %d (%d)\n",
2811 je-js, irank, len);
2812 n_error++;
2813 }
2814 }
2815 HECMW_Waitall(ref_mesh->n_neighbor_pe, requests, statuses);
2816
2817 /* shared */
2818
2819 nsend = 0;
2820 for (i = 0; i < ref_mesh->n_neighbor_pe; i++) {
2821 irank = ref_mesh->neighbor_pe[i];
2822 if (irank < ref_mesh->my_rank) continue;
2823 js = ref_mesh->shared_index[i];
2824 je = ref_mesh->shared_index[i + 1];
2825 send_buf[i] = je - js;
2826 /* isend number of shared elems */
2827 tag = 1;
2828 HECMW_Isend(send_buf + i, 1, HECMW_INT, irank, tag, comm, requests + nsend);
2829 nsend++;
2830 }
2831 for (i = 0; i < ref_mesh->n_neighbor_pe; i++) {
2832 irank = ref_mesh->neighbor_pe[i];
2833 if (ref_mesh->my_rank < irank) continue;
2834 js = ref_mesh->shared_index[i];
2835 je = ref_mesh->shared_index[i + 1];
2836 /* recv number of shared elems */
2837 tag = 1;
2838 HECMW_Recv(&len, 1, HECMW_INT, irank, tag, comm, statuses);
2839 if (len != je - js) {
2841 "inconsistent length of shared (%d) with rank %d (%d)\n",
2842 je-js, irank, len);
2843 n_error++;
2844 }
2845 }
2846 HECMW_Waitall(nsend, requests, statuses);
2847
2848 HECMW_free(send_buf);
2849 HECMW_free(requests);
2850 HECMW_free(statuses);
2851
2852 HECMW_Allreduce(&n_error, &n_error_g, 1, HECMW_INT, HECMW_SUM, comm);
2853 if (n_error_g > 0) return HECMW_ERROR;
2854
2855 HECMW_log(HECMW_LOG_DEBUG, "Finished checking communication tables.\n");
2856 return HECMW_SUCCESS;
2857}
2858
2859static int rebuild_ID_external(struct hecmwST_local_mesh *ref_mesh) {
2860 int len_tot;
2861 int *sendbuf, *recvbuf, *srbuf;
2862 int i, j, irank, js, je, len, nid, tag, nsend, rank, lid;
2863 int *item_p, *sendbuf_p, *recvbuf_p, *srbuf_p;
2864 HECMW_Request *requests;
2865 HECMW_Status *statuses;
2866 HECMW_Comm comm;
2867
2868 if (ref_mesh->n_neighbor_pe == 0) return HECMW_SUCCESS;
2869
2870 HECMW_log(HECMW_LOG_DEBUG, "Started rebuilding external IDs...\n");
2871
2872 comm = HECMW_comm_get_comm();
2873
2874 requests = (HECMW_Request *)HECMW_malloc(sizeof(HECMW_Request) *
2875 ref_mesh->n_neighbor_pe);
2876 if (requests == NULL) {
2877 HECMW_set_error(errno, "");
2878 return HECMW_ERROR;
2879 }
2880 statuses = (HECMW_Status *)HECMW_malloc(sizeof(HECMW_Status) *
2881 ref_mesh->n_neighbor_pe);
2882 if (statuses == NULL) {
2883 HECMW_set_error(errno, "");
2884 return HECMW_ERROR;
2885 }
2886
2887 /* node_ID (for external nodes) */
2888
2889 /* send local IDs of export nodes */
2890 len_tot = ref_mesh->export_index[ref_mesh->n_neighbor_pe];
2891 sendbuf = (int *)HECMW_malloc(sizeof(int) * len_tot);
2892 if (sendbuf == NULL) {
2893 HECMW_set_error(errno, "");
2894 return HECMW_ERROR;
2895 }
2896 for (i = 0; i < ref_mesh->n_neighbor_pe; i++) {
2897 irank = ref_mesh->neighbor_pe[i];
2898 js = ref_mesh->export_index[i];
2899 je = ref_mesh->export_index[i + 1];
2900 len = je - js;
2901 item_p = ref_mesh->export_item + js;
2902 sendbuf_p = sendbuf + js;
2903 for (j = 0; j < len; j++) {
2904 nid = item_p[j] - 1;
2905 lid = ref_mesh->node_ID[2 * nid];
2906 HECMW_assert(0 < lid);
2907 HECMW_assert(lid <= ref_mesh->nn_internal);
2908 sendbuf_p[j] = lid;
2909 }
2910 /* isend local id of export nodes */
2911 tag = 0;
2912 HECMW_Isend(sendbuf_p, len, HECMW_INT, irank, tag, comm, requests + i);
2913 }
2914
2915 /* recv local IDs of import nodes */
2916 len_tot = ref_mesh->import_index[ref_mesh->n_neighbor_pe];
2917 recvbuf = (int *)HECMW_malloc(sizeof(int) * len_tot);
2918 if (recvbuf == NULL) {
2919 HECMW_set_error(errno, "");
2920 return HECMW_ERROR;
2921 }
2922 for (i = 0; i < ref_mesh->n_neighbor_pe; i++) {
2923 irank = ref_mesh->neighbor_pe[i];
2924 js = ref_mesh->import_index[i];
2925 je = ref_mesh->import_index[i + 1];
2926 len = je - js;
2927 item_p = ref_mesh->import_item + js;
2928 recvbuf_p = recvbuf + js;
2929 /* recv local id of import nodes */
2930 tag = 0;
2931 HECMW_Recv(recvbuf_p, len, HECMW_INT, irank, tag, comm, statuses + i);
2932 /* set node_ID[2*j] */
2933 for (j = 0; j < len; j++) {
2934 nid = item_p[j] - 1;
2935 lid = recvbuf_p[j];
2936 HECMW_assert(0 < lid);
2937 ref_mesh->node_ID[2 * nid] = lid;
2938 }
2939 }
2940
2941 HECMW_Waitall(ref_mesh->n_neighbor_pe, requests, statuses);
2942 HECMW_free(sendbuf);
2943 HECMW_free(recvbuf);
2944
2945 /* elem_ID (for external elements) */
2946
2947 len_tot = ref_mesh->shared_index[ref_mesh->n_neighbor_pe];
2948 srbuf = (int *)HECMW_malloc(sizeof(int) * len_tot);
2949 if (srbuf == NULL) {
2950 HECMW_set_error(errno, "");
2951 return HECMW_ERROR;
2952 }
2953 nsend = 0;
2954 for (i = 0; i < ref_mesh->n_neighbor_pe; i++) {
2955 irank = ref_mesh->neighbor_pe[i];
2956 if (irank < ref_mesh->my_rank) continue;
2957 js = ref_mesh->shared_index[i];
2958 je = ref_mesh->shared_index[i + 1];
2959 len = je - js;
2960 srbuf_p = srbuf + js;
2961 item_p = ref_mesh->shared_item + js;
2962 for (j = 0; j < len; j++) {
2963 nid = item_p[j] - 1;
2964 rank = ref_mesh->elem_ID[2 * nid + 1];
2965 lid = ref_mesh->elem_ID[2 * nid];
2966 if (rank == ref_mesh->my_rank && lid > 0) {
2967 HECMW_assert(lid <= ref_mesh->ne_internal);
2968 srbuf_p[j] = lid;
2969 } else {
2970 srbuf_p[j] = -1;
2971 }
2972 }
2973 /* isend list of local_ID of those elems */
2974 tag = 1;
2975 HECMW_Isend(srbuf_p, len, HECMW_INT, irank, tag, comm, requests + nsend);
2976 nsend++;
2977 }
2978 for (i = 0; i < ref_mesh->n_neighbor_pe; i++) {
2979 irank = ref_mesh->neighbor_pe[i];
2980 if (ref_mesh->my_rank < irank) continue;
2981 js = ref_mesh->shared_index[i];
2982 je = ref_mesh->shared_index[i + 1];
2983 len = je - js;
2984 srbuf_p = srbuf + js;
2985 item_p = ref_mesh->shared_item + js;
2986 /* recv list of local_ID of those elems */
2987 tag = 1;
2988 HECMW_Recv(srbuf_p, len, HECMW_INT, irank, tag, comm, statuses);
2989 /* set elem_ID[2*j] */
2990 for (j = 0; j < len; j++) {
2991 lid = srbuf_p[j];
2992 if (lid < 0) continue;
2993 nid = item_p[j] - 1;
2994 rank = ref_mesh->elem_ID[2 * nid + 1];
2995 HECMW_assert(rank == irank || rank == -irank - 1);
2996 if (rank < 0) continue;
2997 ref_mesh->elem_ID[2 * nid] = lid;
2998 }
2999 }
3000 HECMW_Waitall(nsend, requests, statuses);
3001 HECMW_free(srbuf);
3002
3003 HECMW_free(requests);
3004 HECMW_free(statuses);
3005
3006 HECMW_log(HECMW_LOG_DEBUG, "Finished rebuilding external IDs.\n");
3007 return HECMW_SUCCESS;
3008}
3009
3010static int rebuild_refine_origin(const struct hecmwST_local_mesh *mesh,
3011 struct hecmwST_local_mesh *ref_mesh) {
3012 const struct hecmwST_refine_origin *reforg;
3013 struct hecmwST_refine_origin *ref_reforg;
3014 struct hecmw_varray_int ro_item;
3015 const int *ro_item_v;
3016 int n_refine, offset, cnt, cnt2, i, j;
3017
3018 reforg = mesh->refine_origin;
3019 n_refine = mesh->n_refine + 1;
3020
3021 ref_reforg = (struct hecmwST_refine_origin *)HECMW_malloc(
3022 sizeof(struct hecmwST_refine_origin));
3023 if (ref_reforg == NULL) {
3024 HECMW_set_error(errno, "");
3025 return HECMW_ERROR;
3026 }
3027
3028 ref_reforg->index = (int *)HECMW_malloc(sizeof(int) * (n_refine + 1));
3029 if (ref_reforg->index == NULL) {
3030 HECMW_set_error(errno, "");
3031 return HECMW_ERROR;
3032 }
3033 ref_reforg->index[0] = 0;
3034 for (i = 1; i < n_refine; i++) {
3035 ref_reforg->index[i] = reforg->index[i];
3036 }
3037 ref_reforg->index[n_refine] =
3038 ref_reforg->index[n_refine - 1] + ref_mesh->n_node;
3039
3040 ref_reforg->item_index =
3041 (int *)HECMW_malloc(sizeof(int) * (ref_reforg->index[n_refine] + 1));
3042 if (ref_reforg->item_index == NULL) {
3043 HECMW_set_error(errno, "");
3044 return HECMW_ERROR;
3045 }
3046
3047 /* copy original-node infomation up to previous refinement */
3048 ref_reforg->item_index[0] = 0;
3049 for (i = 1; i <= ref_reforg->index[n_refine - 1]; i++) {
3050 ref_reforg->item_index[i] = reforg->item_index[i];
3051 }
3052 offset = ref_reforg->index[n_refine - 1];
3053 cnt = ref_reforg->item_index[offset];
3054 /* fprintf(stderr, "offset=%d, cnt=%d\n", offset, cnt); */
3055
3056 /* get original nodes of newly generated nodes at current refinement */
3057 HECMW_varray_int_init(&ro_item);
3058 for (i = 0; i < ref_mesh->n_node; i++) {
3059 int iold;
3060 iold = ref_mesh->node_new2old ? ref_mesh->node_new2old[i] : i;
3061 if (iold < mesh->n_node_gross) {
3062 int org;
3063 org = mesh->node_old2new ? mesh->node_old2new[iold] + 1 : iold + 1;
3064 HECMW_varray_int_append(&ro_item, org);
3065 cnt++;
3066 } else {
3067 int original[HECMW_MAX_NODE_MAX];
3068 int orig_type_rcap = rcapGetOriginal(iold + 1, original);
3069 int nn;
3070 HECMW_assert(orig_type_rcap != RCAP_UNKNOWNTYPE);
3071 nn = get_rcap_elem_max_node(orig_type_rcap);
3072 for (j = 0; j < nn; j++) {
3073 int org = original[j];
3074 if (mesh->node_old2new) org = mesh->node_old2new[org - 1] + 1;
3075 HECMW_varray_int_append(&ro_item, org);
3076 }
3077 cnt += nn;
3078 }
3079 ref_reforg->item_index[offset + i + 1] = cnt;
3080 }
3081
3082 ref_reforg->item_item = (int *)HECMW_malloc(sizeof(int) * cnt);
3083 if (ref_reforg->item_item == NULL) {
3084 HECMW_set_error(errno, "");
3085 return HECMW_ERROR;
3086 }
3087
3088 /* copy original nodes generated until previous refinements */
3089 cnt = ref_reforg->item_index[ref_reforg->index[n_refine - 1]];
3090 for (i = 0; i < cnt; i++) {
3091 ref_reforg->item_item[i] = reforg->item_item[i];
3092 }
3093 /* copy original nodes generated at current refinement */
3094 cnt2 = HECMW_varray_int_nval(&ro_item);
3095 ro_item_v = HECMW_varray_int_get_cv(&ro_item);
3096 for (i = 0; i < cnt2; i++) {
3097 ref_reforg->item_item[cnt + i] = ro_item_v[i];
3098 }
3099 HECMW_varray_int_finalize(&ro_item);
3100#if 0
3101 fprintf(stderr, "refine_origin->index:\n");
3102 for( i=0; i <= n_refine; i++ ) {
3103 fprintf(stderr, " %d", ref_reforg->index[i]);
3104 if( i % 10 == 9 ) fprintf(stderr, "\n");
3105 }
3106 if( i % 10 != 0 ) fprintf(stderr, "\n");
3107 fprintf(stderr, "refine_origin->item_index:\n");
3108 for( i=0; i <= ref_reforg->index[n_refine]; i++) {
3109 fprintf(stderr, " %d", ref_reforg->item_index[i]);
3110 if( i % 10 == 9 ) fprintf(stderr, "\n");
3111 }
3112 if( i % 10 != 0 ) fprintf(stderr, "\n");
3113 fprintf(stderr, "refine_origin->item_item:\n");
3114 for( i=0; i < ref_reforg->item_index[ref_reforg->index[n_refine]]; i++ ) {
3115 fprintf(stderr, " %d", ref_reforg->item_item[i]);
3116 if( i % 10 == 9 ) fprintf(stderr, "\n");
3117 }
3118 if( i % 10 != 0 ) fprintf(stderr, "\n");
3119#endif
3120 ref_mesh->refine_origin = ref_reforg;
3121 return HECMW_SUCCESS;
3122}
3123
3124static int rebuild_n_node_refine_hist(const struct hecmwST_local_mesh *mesh,
3125 struct hecmwST_local_mesh *ref_mesh) {
3126 int i;
3127
3128 ref_mesh->n_node_refine_hist =
3129 (int *)HECMW_malloc((mesh->n_refine + 1) * sizeof(int));
3130 if (ref_mesh->n_node_refine_hist == NULL) {
3131 HECMW_set_error(errno, "");
3132 return HECMW_ERROR;
3133 }
3134
3135 for (i = 0; i < mesh->n_refine; i++) {
3136 ref_mesh->n_node_refine_hist[i] = mesh->n_node_refine_hist[i];
3137 }
3138 ref_mesh->n_node_refine_hist[mesh->n_refine] = ref_mesh->n_node_gross;
3139
3140 return HECMW_SUCCESS;
3141}
3142
3143static int rebuild_refine_info(const struct hecmwST_local_mesh *mesh,
3144 struct hecmwST_local_mesh *ref_mesh) {
3145 ref_mesh->n_refine = mesh->n_refine + 1;
3146
3147 if (rebuild_refine_origin(mesh, ref_mesh) != HECMW_SUCCESS)
3148 return HECMW_ERROR;
3149 if (rebuild_n_node_refine_hist(mesh, ref_mesh) != HECMW_SUCCESS)
3150 return HECMW_ERROR;
3151
3152 return HECMW_SUCCESS;
3153}
3154
3155static int renumber_nodes_generate_tables(struct hecmwST_local_mesh *mesh);
3156static int renumber_elements_generate_tables(struct hecmwST_local_mesh *mesh);
3157
3158static int rebuild_info(const struct hecmwST_local_mesh *mesh,
3159 struct hecmwST_local_mesh *ref_mesh) {
3160 HECMW_log(HECMW_LOG_DEBUG, "Started rebuilding info...\n");
3161
3162 if (rebuild_elem_info(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3163 if (rebuild_node_info(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3164 if (rebuild_comm_tables(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3165 if (check_comm_table_len(ref_mesh) != HECMW_SUCCESS) {
3166 HECMW_log(HECMW_LOG_ERROR, "Check communication table failed. Checking original mesh...\n");
3167 if (check_comm_table_len(mesh) != HECMW_SUCCESS) {
3168 HECMW_log(HECMW_LOG_ERROR, "Original mesh has error\n");
3169 } else {
3170 HECMW_log(HECMW_LOG_ERROR, "Original mesh is OK\n");
3171 }
3172 return HECMW_ERROR;
3173 }
3174 if (rebuild_ID_external(ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3175
3176 if (renumber_nodes_generate_tables(ref_mesh) != HECMW_SUCCESS)
3177 return HECMW_ERROR;
3178 if (renumber_elements_generate_tables(ref_mesh) != HECMW_SUCCESS)
3179 return HECMW_ERROR;
3180
3181 if (rebuild_refine_info(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3182
3183 HECMW_log(HECMW_LOG_DEBUG, "Finished rebuilding info.\n");
3184 return HECMW_SUCCESS;
3185}
3186
3187/*============================================================================*/
3188/* */
3189/* copy unchanging infomation */
3190/* */
3191/*============================================================================*/
3192static int copy_global_info(const struct hecmwST_local_mesh *mesh,
3193 struct hecmwST_local_mesh *ref_mesh) {
3194 int i;
3195
3202 strcpy(ref_mesh->gridfile, mesh->gridfile);
3203 ref_mesh->hecmw_n_file = mesh->hecmw_n_file;
3204
3205 /* files */
3206 if (mesh->hecmw_n_file > 0) {
3207 ref_mesh->files = (char **)HECMW_calloc(mesh->hecmw_n_file, sizeof(char *));
3208 if (ref_mesh->files == NULL) {
3209 HECMW_set_error(errno, "");
3210 return HECMW_ERROR;
3211 }
3212 for (i = 0; i < mesh->hecmw_n_file; i++) {
3213 ref_mesh->files[i] = HECMW_strdup(mesh->files[i]);
3214 if (ref_mesh->files[i] == NULL) {
3215 HECMW_set_error(errno, "");
3216 return HECMW_ERROR;
3217 }
3218 }
3219 } else {
3220 ref_mesh->files = NULL;
3221 }
3222
3223 strcpy(ref_mesh->header, mesh->header);
3224 ref_mesh->zero_temp = mesh->zero_temp;
3225
3226 return HECMW_SUCCESS;
3227}
3228
3229static int copy_node_ndof(const struct hecmwST_local_mesh *mesh,
3230 struct hecmwST_local_mesh *ref_mesh) {
3231 int i;
3232
3233 ref_mesh->n_dof = mesh->n_dof;
3234 ref_mesh->n_dof_grp = mesh->n_dof_grp;
3235 ref_mesh->n_dof_tot = mesh->n_dof_tot;
3236
3237 if (mesh->n_dof_grp > 0) {
3238 /* node_dof_index */
3239 ref_mesh->node_dof_index =
3240 (int *)HECMW_malloc(sizeof(int) * (mesh->n_dof_grp + 1));
3241 if (ref_mesh->node_dof_index == NULL) {
3242 HECMW_set_error(errno, "");
3243 return HECMW_ERROR;
3244 }
3245 for (i = 0; i < mesh->n_dof_grp + 1; i++) {
3246 ref_mesh->node_dof_index[i] = mesh->node_dof_index[i];
3247 }
3248
3249 /* node_dof_item */
3250 ref_mesh->node_dof_item =
3251 (int *)HECMW_malloc(sizeof(int) * mesh->n_dof_grp);
3252 if (ref_mesh->node_dof_item == NULL) {
3253 HECMW_set_error(errno, "");
3254 return HECMW_ERROR;
3255 }
3256 for (i = 0; i < mesh->n_dof_grp; i++) {
3257 ref_mesh->node_dof_item[i] = mesh->node_dof_item[i];
3258 }
3259 }
3260 return HECMW_SUCCESS;
3261}
3262
3263static int copy_elem_etype(const struct hecmwST_local_mesh *mesh,
3264 struct hecmwST_local_mesh *ref_mesh) {
3265 int i;
3266
3267 ref_mesh->n_elem_type = mesh->n_elem_type;
3268
3269 if (mesh->n_elem_type > 0) {
3270 /* elem_type_index */
3271 ref_mesh->elem_type_index =
3272 (int *)HECMW_malloc(sizeof(int) * (mesh->n_elem_type + 1));
3273 if (ref_mesh->elem_type_index == NULL) {
3274 HECMW_set_error(errno, "");
3275 return HECMW_ERROR;
3276 }
3277
3278 /* elem_type_item */
3279 ref_mesh->elem_type_item =
3280 (int *)HECMW_malloc(sizeof(int) * mesh->n_elem_type);
3281 if (ref_mesh->elem_type_item == NULL) {
3282 HECMW_set_error(errno, "");
3283 return HECMW_ERROR;
3284 }
3285
3286 ref_mesh->elem_type_index[0] = 0;
3287 for (i = 0; i < mesh->n_elem_type; i++) {
3288 int etype, ierror, ndiv, nelem;
3289 etype = mesh->elem_type_item[i];
3290 ndiv = get_elem_ndiv(etype, &ierror);
3291 if (ierror) return HECMW_ERROR;
3293 ref_mesh->elem_type_index[i + 1] =
3294 ref_mesh->elem_type_index[i] + nelem * ndiv;
3295 ref_mesh->elem_type_item[i] = etype;
3296 }
3297 }
3298 return HECMW_SUCCESS;
3299}
3300
3301static int copy_comm_info(const struct hecmwST_local_mesh *mesh,
3302 struct hecmwST_local_mesh *ref_mesh) {
3303 int i;
3304
3305 ref_mesh->zero = mesh->zero;
3306 ref_mesh->HECMW_COMM = mesh->HECMW_COMM;
3307 ref_mesh->PETOT = mesh->PETOT;
3308 ref_mesh->PEsmpTOT = mesh->PEsmpTOT;
3309 ref_mesh->my_rank = mesh->my_rank;
3310 ref_mesh->errnof = mesh->errnof;
3311 ref_mesh->n_subdomain = mesh->n_subdomain;
3312 ref_mesh->n_neighbor_pe = mesh->n_neighbor_pe;
3313
3314 if (mesh->n_neighbor_pe == 0) {
3315 ref_mesh->neighbor_pe = NULL;
3316 ref_mesh->import_item = NULL;
3317 ref_mesh->export_item = NULL;
3318 ref_mesh->shared_item = NULL;
3319
3320 ref_mesh->import_index = (int *)HECMW_malloc(sizeof(int));
3321 if (ref_mesh->import_index == NULL) {
3322 HECMW_set_error(errno, "");
3323 return HECMW_ERROR;
3324 }
3325 ref_mesh->import_index[0] = 0;
3326
3327 ref_mesh->export_index = (int *)HECMW_malloc(sizeof(int));
3328 if (ref_mesh->export_index == NULL) {
3329 HECMW_set_error(errno, "");
3330 return HECMW_ERROR;
3331 }
3332 ref_mesh->export_index[0] = 0;
3333
3334 ref_mesh->shared_index = (int *)HECMW_malloc(sizeof(int));
3335 if (ref_mesh->shared_index == NULL) {
3336 HECMW_set_error(errno, "");
3337 return HECMW_ERROR;
3338 }
3339 ref_mesh->shared_index[0] = 0;
3340
3341 return HECMW_SUCCESS;
3342 }
3343
3344 /* neighbor_pe */
3345 ref_mesh->neighbor_pe =
3346 (int *)HECMW_malloc(sizeof(int) * mesh->n_neighbor_pe);
3347 if (ref_mesh->neighbor_pe == NULL) {
3348 HECMW_set_error(errno, "");
3349 return HECMW_ERROR;
3350 }
3351 for (i = 0; i < mesh->n_neighbor_pe; i++) {
3352 ref_mesh->neighbor_pe[i] = mesh->neighbor_pe[i];
3353 }
3354
3355 return HECMW_SUCCESS;
3356}
3357
3358static int copy_adapt_info(const struct hecmwST_local_mesh *mesh,
3359 struct hecmwST_local_mesh *ref_mesh) {
3360 if (mesh->hecmw_flag_adapt != 0) {
3361 HECMW_log(HECMW_LOG_ERROR, "Refinement of adaptation not supported\n");
3362 return HECMW_ERROR;
3363 }
3364 ref_mesh->coarse_grid_level = 0;
3365 ref_mesh->n_adapt = 0;
3366 ref_mesh->when_i_was_refined_node = NULL;
3367 ref_mesh->when_i_was_refined_elem = NULL;
3368 ref_mesh->adapt_parent_type = NULL;
3369 ref_mesh->adapt_type = NULL;
3370 ref_mesh->adapt_level = NULL;
3371 ref_mesh->adapt_parent = NULL;
3372 ref_mesh->adapt_children_index = NULL;
3373 ref_mesh->adapt_children_item = NULL;
3374 return HECMW_SUCCESS;
3375}
3376
3377static int copy_section_info(const struct hecmwST_section *sect,
3378 struct hecmwST_section *ref_sect) {
3379 int i;
3380
3381 ref_sect->n_sect = sect->n_sect;
3382
3383 if (sect->n_sect == 0) {
3384 ref_sect->sect_type = NULL;
3385 ref_sect->sect_opt = NULL;
3386 ref_sect->sect_mat_ID_index = NULL;
3387 ref_sect->sect_mat_ID_item = NULL;
3388 ref_sect->sect_I_index = NULL;
3389 ref_sect->sect_I_item = NULL;
3390 ref_sect->sect_R_index = NULL;
3391 ref_sect->sect_R_item = NULL;
3392 return HECMW_SUCCESS;
3393 }
3394
3395 /* sect_type */
3396 ref_sect->sect_type = (int *)HECMW_malloc(sizeof(int) * sect->n_sect);
3397 if (ref_sect->sect_type == NULL) {
3398 HECMW_set_error(errno, "");
3399 return HECMW_ERROR;
3400 }
3401 for (i = 0; i < sect->n_sect; i++) {
3402 ref_sect->sect_type[i] = sect->sect_type[i];
3403 }
3404
3405 /* sect_opt */
3406 ref_sect->sect_opt = (int *)HECMW_malloc(sizeof(int) * sect->n_sect);
3407 if (ref_sect->sect_opt == NULL) {
3408 HECMW_set_error(errno, "");
3409 return HECMW_ERROR;
3410 }
3411 for (i = 0; i < sect->n_sect; i++) {
3412 ref_sect->sect_opt[i] = sect->sect_opt[i];
3413 }
3414
3415 /* sect_mat_ID_index */
3416 ref_sect->sect_mat_ID_index =
3417 (int *)HECMW_malloc(sizeof(int) * (sect->n_sect + 1));
3418 if (ref_sect->sect_mat_ID_index == NULL) {
3419 HECMW_set_error(errno, "");
3420 return HECMW_ERROR;
3421 }
3422 for (i = 0; i < sect->n_sect + 1; i++) {
3423 ref_sect->sect_mat_ID_index[i] = sect->sect_mat_ID_index[i];
3424 }
3425
3426 /* sect_mat_ID_item */
3427 if (sect->sect_mat_ID_index[sect->n_sect] > 0) {
3428 ref_sect->sect_mat_ID_item = (int *)HECMW_malloc(
3429 sizeof(int) * sect->sect_mat_ID_index[sect->n_sect]);
3430 if (ref_sect->sect_mat_ID_item == NULL) {
3431 HECMW_set_error(errno, "");
3432 return HECMW_ERROR;
3433 }
3434 for (i = 0; i < sect->sect_mat_ID_index[sect->n_sect]; i++) {
3435 ref_sect->sect_mat_ID_item[i] = sect->sect_mat_ID_item[i];
3436 }
3437 }
3438
3439 /* sect_I_index */
3440 ref_sect->sect_I_index =
3441 (int *)HECMW_malloc(sizeof(int) * (sect->n_sect + 1));
3442 if (ref_sect->sect_I_index == NULL) {
3443 HECMW_set_error(errno, "");
3444 return HECMW_ERROR;
3445 }
3446 for (i = 0; i < sect->n_sect + 1; i++) {
3447 ref_sect->sect_I_index[i] = sect->sect_I_index[i];
3448 }
3449
3450 /* sect_I_item */
3451 if (sect->sect_I_index[sect->n_sect] > 0) {
3452 ref_sect->sect_I_item =
3453 (int *)HECMW_malloc(sizeof(int) * sect->sect_I_index[sect->n_sect]);
3454 if (ref_sect->sect_I_item == NULL) {
3455 HECMW_set_error(errno, "");
3456 return HECMW_ERROR;
3457 }
3458 for (i = 0; i < sect->sect_I_index[sect->n_sect]; i++) {
3459 ref_sect->sect_I_item[i] = sect->sect_I_item[i];
3460 }
3461 }
3462
3463 /* sect_R_index */
3464 ref_sect->sect_R_index =
3465 (int *)HECMW_malloc(sizeof(int) * (sect->n_sect + 1));
3466 if (ref_sect->sect_R_index == NULL) {
3467 HECMW_set_error(errno, "");
3468 return HECMW_ERROR;
3469 }
3470 for (i = 0; i < sect->n_sect + 1; i++) {
3471 ref_sect->sect_R_index[i] = sect->sect_R_index[i];
3472 }
3473
3474 /* sect_R_item */
3475 if (sect->sect_R_index[sect->n_sect] > 0) {
3476 ref_sect->sect_R_item = (double *)HECMW_malloc(
3477 sizeof(double) * sect->sect_R_index[sect->n_sect]);
3478 if (ref_sect->sect_R_item == NULL) {
3479 HECMW_set_error(errno, "");
3480 return HECMW_ERROR;
3481 }
3482 for (i = 0; i < sect->sect_R_index[sect->n_sect]; i++) {
3483 ref_sect->sect_R_item[i] = sect->sect_R_item[i];
3484 }
3485 }
3486
3487 return HECMW_SUCCESS;
3488}
3489
3490static int copy_material_info(const struct hecmwST_material *mat,
3491 struct hecmwST_material *ref_mat) {
3492 int i;
3493
3494 ref_mat->n_mat = mat->n_mat;
3495
3496 if (mat->n_mat == 0) {
3497 ref_mat->n_mat_item = 0;
3498 ref_mat->n_mat_subitem = 0;
3499 ref_mat->n_mat_table = 0;
3500 ref_mat->mat_name = NULL;
3501 ref_mat->mat_item_index = NULL;
3502 ref_mat->mat_subitem_index = NULL;
3503 ref_mat->mat_table_index = NULL;
3504 ref_mat->mat_val = NULL;
3505 ref_mat->mat_temp = NULL;
3506 return HECMW_SUCCESS;
3507 }
3508
3509 ref_mat->n_mat_item = mat->n_mat_item;
3510 ref_mat->n_mat_subitem = mat->n_mat_subitem;
3511 ref_mat->n_mat_table = mat->n_mat_table;
3512
3513 /* mat_name */
3514 ref_mat->mat_name = (char **)HECMW_malloc(sizeof(char *) * mat->n_mat);
3515 if (ref_mat->mat_name == NULL) {
3516 HECMW_set_error(errno, "");
3517 return HECMW_ERROR;
3518 }
3519 for (i = 0; i < mat->n_mat; i++) {
3520 ref_mat->mat_name[i] = HECMW_strdup(mat->mat_name[i]);
3521 if (ref_mat->mat_name[i] == NULL) {
3522 HECMW_set_error(errno, "");
3523 return HECMW_ERROR;
3524 }
3525 }
3526
3527 /* mat_item_index */
3528 ref_mat->mat_item_index = (int *)HECMW_malloc(sizeof(int) * (mat->n_mat + 1));
3529 if (ref_mat->mat_item_index == NULL) {
3530 HECMW_set_error(errno, "");
3531 return HECMW_ERROR;
3532 }
3533 for (i = 0; i < mat->n_mat + 1; i++) {
3534 ref_mat->mat_item_index[i] = mat->mat_item_index[i];
3535 }
3536
3537 /* mat_subitem_index */
3538 ref_mat->mat_subitem_index =
3539 (int *)HECMW_malloc(sizeof(int) * (mat->n_mat_item + 1));
3540 if (ref_mat->mat_subitem_index == NULL) {
3541 HECMW_set_error(errno, "");
3542 return HECMW_ERROR;
3543 }
3544 for (i = 0; i < mat->n_mat_item + 1; i++) {
3545 ref_mat->mat_subitem_index[i] = mat->mat_subitem_index[i];
3546 }
3547
3548 /* mat_table_index */
3549 ref_mat->mat_table_index =
3550 (int *)HECMW_malloc(sizeof(int) * (mat->n_mat_subitem + 1));
3551 if (ref_mat->mat_table_index == NULL) {
3552 HECMW_set_error(errno, "");
3553 return HECMW_ERROR;
3554 }
3555 for (i = 0; i < mat->n_mat_subitem + 1; i++) {
3556 ref_mat->mat_table_index[i] = mat->mat_table_index[i];
3557 }
3558
3559 /* mat_val */
3560 ref_mat->mat_val = (double *)HECMW_malloc(sizeof(double) * mat->n_mat_table);
3561 if (ref_mat->mat_val == NULL) {
3562 HECMW_set_error(errno, "");
3563 return HECMW_ERROR;
3564 }
3565 for (i = 0; i < mat->n_mat_table; i++) {
3566 ref_mat->mat_val[i] = mat->mat_val[i];
3567 }
3568
3569 /* mat_temp */
3570 ref_mat->mat_temp = (double *)HECMW_malloc(sizeof(double) * mat->n_mat_table);
3571 if (ref_mat->mat_temp == NULL) {
3572 HECMW_set_error(errno, "");
3573 return HECMW_ERROR;
3574 }
3575 for (i = 0; i < mat->n_mat_table; i++) {
3576 ref_mat->mat_temp[i] = mat->mat_temp[i];
3577 }
3578
3579 return HECMW_SUCCESS;
3580}
3581
3582static int copy_mpc_info(const struct hecmwST_mpc *mpc,
3583 struct hecmwST_mpc *ref_mpc) {
3584 int i;
3585
3586 ref_mpc->n_mpc = mpc->n_mpc;
3587
3588 if (mpc->n_mpc == 0) {
3589 ref_mpc->mpc_item = NULL;
3590 ref_mpc->mpc_dof = NULL;
3591 ref_mpc->mpc_val = NULL;
3592 ref_mpc->mpc_const = NULL;
3593 ref_mpc->mpc_index = HECMW_malloc(sizeof(int));
3594 if (ref_mpc->mpc_index == NULL) {
3595 HECMW_set_error(errno, "");
3596 return HECMW_ERROR;
3597 }
3598 ref_mpc->mpc_index[0] = 0;
3599 return HECMW_SUCCESS;
3600 }
3601
3602#if 0
3603 HECMW_log(HECMW_LOG_ERROR, "refinement of MPC information is not supported\n");
3604 return HECMW_ERROR;
3605
3606#else
3607 HECMW_log(HECMW_LOG_WARN, "MPC information is not refined\n");
3608
3609 /* mpc_index */
3610 ref_mpc->mpc_index = (int *) HECMW_malloc( sizeof(int) * (mpc->n_mpc+1) );
3611 if( ref_mpc->mpc_index == NULL ) {
3612 HECMW_set_error(errno, "");
3613 return HECMW_ERROR;
3614 }
3615 for ( i = 0; i < mpc->n_mpc + 1; i++ ) {
3616 ref_mpc->mpc_index[i] = mpc->mpc_index[i];
3617 }
3618
3619 /* mpc_item */
3620 ref_mpc->mpc_item = (int *) HECMW_malloc( sizeof(int) * mpc->mpc_index[mpc->n_mpc] );
3621 if( ref_mpc->mpc_item == NULL ) {
3622 HECMW_set_error(errno, "");
3623 return HECMW_ERROR;
3624 }
3625 for ( i = 0; i < mpc->mpc_index[mpc->n_mpc]; i++ ) {
3626 ref_mpc->mpc_item[i] = mpc->mpc_item[i];
3627 }
3628
3629 /* mpc_dof */
3630 ref_mpc->mpc_dof = (int *) HECMW_malloc( sizeof(int) * mpc->mpc_index[mpc->n_mpc] );
3631 if( ref_mpc->mpc_dof == NULL ) {
3632 HECMW_set_error(errno, "");
3633 return HECMW_ERROR;
3634 }
3635 for ( i = 0; i < mpc->mpc_index[mpc->n_mpc]; i++ ) {
3636 ref_mpc->mpc_dof[i] = mpc->mpc_dof[i];
3637 }
3638
3639 /* mpc_val */
3640 ref_mpc->mpc_val = (double *) HECMW_malloc( sizeof(double) * mpc->mpc_index[mpc->n_mpc] );
3641 if( ref_mpc->mpc_val == NULL ) {
3642 HECMW_set_error(errno, "");
3643 return HECMW_ERROR;
3644 }
3645 for ( i = 0; i < mpc->mpc_index[mpc->n_mpc]; i++ ) {
3646 ref_mpc->mpc_val[i] = mpc->mpc_val[i];
3647 }
3648
3649 /* mpc_const */
3650 ref_mpc->mpc_const = (double *) HECMW_malloc( sizeof(double) * mpc->n_mpc );
3651 if( ref_mpc->mpc_const == NULL ) {
3652 HECMW_set_error(errno, "");
3653 return HECMW_ERROR;
3654 }
3655 for ( i = 0; i < mpc->n_mpc; i++ ) {
3656 ref_mpc->mpc_const[i] = mpc->mpc_const[i];
3657 }
3658
3659 return HECMW_SUCCESS;
3660#endif
3661}
3662
3663static int copy_amp_info(const struct hecmwST_amplitude *amp,
3664 struct hecmwST_amplitude *ref_amp) {
3665 int i;
3666
3667 ref_amp->n_amp = amp->n_amp;
3668
3669 if (amp->n_amp == 0) {
3670 ref_amp->amp_name = NULL;
3671 ref_amp->amp_type_definition = NULL;
3672 ref_amp->amp_type_time = NULL;
3673 ref_amp->amp_type_value = NULL;
3674 ref_amp->amp_val = NULL;
3675 ref_amp->amp_table = NULL;
3676 ref_amp->amp_index = (int *)HECMW_malloc(sizeof(int));
3677 if (ref_amp->amp_index == NULL) {
3678 HECMW_set_error(errno, "");
3679 return HECMW_ERROR;
3680 }
3681 ref_amp->amp_index[0] = 0;
3682 return HECMW_SUCCESS;
3683 }
3684
3685 /* amp_name */
3686 ref_amp->amp_name = (char **)HECMW_malloc(sizeof(char *) * amp->n_amp);
3687 if (ref_amp->amp_name == NULL) {
3688 HECMW_set_error(errno, "");
3689 return HECMW_ERROR;
3690 }
3691 for (i = 0; i < amp->n_amp; i++) {
3692 ref_amp->amp_name[i] = HECMW_strdup(amp->amp_name[i]);
3693 if (ref_amp->amp_name[i] == NULL) {
3694 HECMW_set_error(errno, "");
3695 return HECMW_ERROR;
3696 }
3697 }
3698
3699 /* amp_type_definition */
3700 ref_amp->amp_type_definition = (int *)HECMW_malloc(sizeof(int) * amp->n_amp);
3701 if (ref_amp->amp_type_definition == NULL) {
3702 HECMW_set_error(errno, "");
3703 return HECMW_ERROR;
3704 }
3705 for (i = 0; i < amp->n_amp; i++) {
3706 ref_amp->amp_type_definition[i] = amp->amp_type_definition[i];
3707 }
3708
3709 /* amp_type_time */
3710 ref_amp->amp_type_time = (int *)HECMW_malloc(sizeof(int) * amp->n_amp);
3711 if (ref_amp->amp_type_time == NULL) {
3712 HECMW_set_error(errno, "");
3713 return HECMW_ERROR;
3714 }
3715 for (i = 0; i < amp->n_amp; i++) {
3716 ref_amp->amp_type_time[i] = amp->amp_type_time[i];
3717 }
3718
3719 /* amp_type_value */
3720 ref_amp->amp_type_value = (int *)HECMW_malloc(sizeof(int) * amp->n_amp);
3721 if (ref_amp->amp_type_value == NULL) {
3722 HECMW_set_error(errno, "");
3723 return HECMW_ERROR;
3724 }
3725 for (i = 0; i < amp->n_amp; i++) {
3726 ref_amp->amp_type_value[i] = amp->amp_type_value[i];
3727 }
3728
3729 /* amp_index */
3730 ref_amp->amp_index = (int *)HECMW_malloc(sizeof(int) * (amp->n_amp + 1));
3731 if (ref_amp->amp_index == NULL) {
3732 HECMW_set_error(errno, "");
3733 return HECMW_ERROR;
3734 }
3735 for (i = 0; i < amp->n_amp + 1; i++) {
3736 ref_amp->amp_index[i] = amp->amp_index[i];
3737 }
3738
3739 /* amp_val */
3740 ref_amp->amp_val =
3741 (double *)HECMW_malloc(sizeof(double) * amp->amp_index[amp->n_amp]);
3742 if (ref_amp->amp_val == NULL) {
3743 HECMW_set_error(errno, "");
3744 return HECMW_ERROR;
3745 }
3746 for (i = 0; i < amp->amp_index[amp->n_amp]; i++) {
3747 ref_amp->amp_val[i] = amp->amp_val[i];
3748 }
3749
3750 /* amp_table */
3751 ref_amp->amp_table =
3752 (double *)HECMW_malloc(sizeof(double) * amp->amp_index[amp->n_amp]);
3753 if (ref_amp->amp_table == NULL) {
3754 HECMW_set_error(errno, "");
3755 return HECMW_ERROR;
3756 }
3757 for (i = 0; i < amp->amp_index[amp->n_amp]; i++) {
3758 ref_amp->amp_table[i] = amp->amp_table[i];
3759 }
3760
3761 return HECMW_SUCCESS;
3762}
3763
3764static int copy_contact_pair(const struct hecmwST_contact_pair *cp,
3765 struct hecmwST_contact_pair *ref_cp) {
3766 int i;
3767
3768 ref_cp->n_pair = cp->n_pair;
3769
3770 /* name */
3771 ref_cp->name = (char **)HECMW_malloc(sizeof(char *) * cp->n_pair);
3772 if (ref_cp->name == NULL) {
3773 HECMW_set_error(errno, "");
3774 return HECMW_ERROR;
3775 }
3776 for (i = 0; i < cp->n_pair; i++) {
3777 ref_cp->name[i] = HECMW_strdup(cp->name[i]);
3778 if (ref_cp->name[i] == NULL) {
3779 HECMW_set_error(errno, "");
3780 return HECMW_ERROR;
3781 }
3782 }
3783
3784 /* type */
3785 ref_cp->type = (int *)HECMW_malloc(sizeof(int) * cp->n_pair);
3786 if (ref_cp->type == NULL) {
3787 HECMW_set_error(errno, "");
3788 return HECMW_ERROR;
3789 }
3790 for (i = 0; i < cp->n_pair; i++) {
3791 ref_cp->type[i] = cp->type[i];
3792 }
3793
3794 /* slave_grp_id */
3795 ref_cp->slave_grp_id = (int *)HECMW_malloc(sizeof(int) * cp->n_pair);
3796 if (ref_cp->slave_grp_id == NULL) {
3797 HECMW_set_error(errno, "");
3798 return HECMW_ERROR;
3799 }
3800 for (i = 0; i < cp->n_pair; i++) {
3801 ref_cp->slave_grp_id[i] = cp->slave_grp_id[i];
3802 }
3803
3804 /* slave_orisgrp_id */
3805 ref_cp->slave_orisgrp_id = (int *)HECMW_malloc(sizeof(int) * cp->n_pair);
3806 if (ref_cp->slave_orisgrp_id == NULL) {
3807 HECMW_set_error(errno, "");
3808 return HECMW_ERROR;
3809 }
3810 for (i = 0; i < cp->n_pair; i++) {
3811 ref_cp->slave_orisgrp_id[i] = cp->slave_orisgrp_id[i];
3812 }
3813
3814 /* master_grp_id */
3815 ref_cp->master_grp_id = (int *)HECMW_malloc(sizeof(int) * cp->n_pair);
3816 if (ref_cp->master_grp_id == NULL) {
3817 HECMW_set_error(errno, "");
3818 return HECMW_ERROR;
3819 }
3820 for (i = 0; i < cp->n_pair; i++) {
3821 ref_cp->master_grp_id[i] = cp->master_grp_id[i];
3822 }
3823
3824 return HECMW_SUCCESS;
3825}
3826
3827static int copy_unchanging_info(const struct hecmwST_local_mesh *mesh,
3828 struct hecmwST_local_mesh *ref_mesh) {
3829 HECMW_log(HECMW_LOG_DEBUG, "Started copying unchanging info...\n");
3830
3831 if (copy_global_info(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3832 if (copy_node_ndof(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3833 if (copy_elem_etype(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3834 if (copy_comm_info(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3835 if (copy_adapt_info(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
3836 if (copy_section_info(mesh->section, ref_mesh->section) != HECMW_SUCCESS)
3837 return HECMW_ERROR;
3838 if (copy_material_info(mesh->material, ref_mesh->material) != HECMW_SUCCESS)
3839 return HECMW_ERROR;
3840 if (copy_mpc_info(mesh->mpc, ref_mesh->mpc) != HECMW_SUCCESS)
3841 return HECMW_ERROR;
3842 if (copy_amp_info(mesh->amp, ref_mesh->amp) != HECMW_SUCCESS)
3843 return HECMW_ERROR;
3844 if (copy_contact_pair(mesh->contact_pair, ref_mesh->contact_pair) !=
3846 return HECMW_ERROR;
3847
3848 HECMW_log(HECMW_LOG_DEBUG, "Finished copying unchanging info.\n");
3849 return HECMW_SUCCESS;
3850}
3851
3852/*============================================================================*/
3853/* */
3854/* renumber nodes */
3855/* */
3856/*============================================================================*/
3857static int renumber_nodes_generate_tables(struct hecmwST_local_mesh *mesh) {
3858 int i;
3859 int count_in, count_ex, count_pex;
3860 int *old2new, *new2old;
3861
3863
3864 HECMW_log(HECMW_LOG_DEBUG, "Started generating renumbering tables...\n");
3865
3866 old2new = (int *)HECMW_malloc(sizeof(int) * mesh->n_node_gross);
3867 if (old2new == NULL) {
3868 HECMW_set_error(errno, "");
3869 return HECMW_ERROR;
3870 }
3871 new2old = (int *)HECMW_malloc(sizeof(int) * mesh->n_node_gross);
3872 if (new2old == NULL) {
3873 HECMW_set_error(errno, "");
3874 return HECMW_ERROR;
3875 }
3876 count_in = 0;
3877 count_ex = mesh->nn_internal;
3878 count_pex = mesh->n_node;
3879 for (i = 0; i < mesh->n_node_gross; i++) {
3880 int rank = mesh->node_ID[i * 2 + 1];
3881 if (rank == mesh->my_rank) {
3882 old2new[i] = count_in;
3883 new2old[count_in] = i;
3884 count_in++;
3885 } else if (rank >= 0) {
3886 old2new[i] = count_ex;
3887 new2old[count_ex] = i;
3888 count_ex++;
3889 } else {
3890 old2new[i] = count_pex;
3891 new2old[count_pex] = i;
3892 count_pex++;
3893 }
3894 }
3895 mesh->node_old2new = old2new;
3896 mesh->node_new2old = new2old;
3897
3898 HECMW_log(HECMW_LOG_DEBUG, "Finished generating renumbering tables.\n");
3899 return HECMW_SUCCESS;
3900}
3901
3902static int reorder_nodes(struct hecmwST_local_mesh *mesh, int *old2new,
3903 int *new2old) {
3904 int i;
3905
3906#ifndef NDEBUG
3907 for (i = 0; i < mesh->n_node_gross; i++) {
3908 HECMW_assert(new2old[old2new[i]] == i);
3909 }
3910#endif
3911
3912 /*
3913 * Reorder using new2old
3914 */
3915 /* node_ID */
3916 {
3917 int *new_node_ID =
3918 (int *)HECMW_malloc(sizeof(int) * mesh->n_node_gross * 2);
3919 if (new_node_ID == NULL) {
3920 HECMW_set_error(errno, "");
3921 return HECMW_ERROR;
3922 }
3923 for (i = 0; i < mesh->n_node_gross; i++) {
3924 int old = new2old[i];
3925 new_node_ID[2 * i] = mesh->node_ID[2 * old];
3926 new_node_ID[2 * i + 1] = mesh->node_ID[2 * old + 1];
3927 }
3929 mesh->node_ID = new_node_ID;
3930 }
3931 /* global_node_ID */
3932 {
3933 int *new_global_node_ID =
3934 (int *)HECMW_malloc(sizeof(int) * mesh->n_node_gross);
3935 if (new_global_node_ID == NULL) {
3936 HECMW_set_error(errno, "");
3937 return HECMW_ERROR;
3938 }
3939 for (i = 0; i < mesh->n_node_gross; i++) {
3940 int old = new2old[i];
3941 new_global_node_ID[i] = mesh->global_node_ID[old];
3942 }
3944 mesh->global_node_ID = new_global_node_ID;
3945 }
3946 /* node */
3947 {
3948 double *new_node =
3949 (double *)HECMW_malloc(sizeof(double) * mesh->n_node_gross * 3);
3950 if (new_node == NULL) {
3951 HECMW_set_error(errno, "");
3952 return HECMW_ERROR;
3953 }
3954 for (i = 0; i < mesh->n_node_gross; i++) {
3955 int old = new2old[i];
3956 int j;
3957 for (j = 0; j < 3; j++) new_node[3 * i + j] = mesh->node[3 * old + j];
3958 }
3960 mesh->node = new_node;
3961 }
3962 /* node_init_val_index, node_init_val_item */
3964 /* node_init_val_index */
3965 int *new_node_init_val_index =
3966 (int *)HECMW_malloc(sizeof(int) * (mesh->n_node_gross + 1));
3967 if (new_node_init_val_index == NULL) {
3968 HECMW_set_error(errno, "");
3969 return HECMW_ERROR;
3970 }
3971 new_node_init_val_index[0] = 0;
3972 for (i = 0; i < mesh->n_node_gross; i++) {
3973 int old = new2old[i];
3974 int ninit =
3976 new_node_init_val_index[i + 1] = new_node_init_val_index[i] + ninit;
3977 }
3978
3979 /* node_init_val_item */
3981 double *new_node_init_val_item = (double *)HECMW_malloc(
3982 sizeof(double) * mesh->node_init_val_index[mesh->n_node_gross]);
3983 if (new_node_init_val_item == NULL) {
3984 HECMW_set_error(errno, "");
3985 return HECMW_ERROR;
3986 }
3987 for (i = 0; i < mesh->n_node_gross; i++) {
3988 int old = new2old[i];
3989 int idx_new = new_node_init_val_index[i];
3990 int idx_old = mesh->node_init_val_index[old];
3991 int ninit = mesh->node_init_val_index[old + 1] - idx_old;
3992 int j;
3993 for (j = 0; j < ninit; j++) {
3994 new_node_init_val_item[idx_new + j] =
3995 mesh->node_init_val_item[idx_old + j];
3996 }
3997 }
3999 mesh->node_init_val_item = new_node_init_val_item;
4000 }
4002 mesh->node_init_val_index = new_node_init_val_index;
4003 }
4004
4005 /*
4006 * Update using old2new
4007 */
4008 /* elem_node_item */
4009 for (i = 0; i < mesh->elem_node_index[mesh->n_elem_gross]; i++) {
4010 int old = mesh->elem_node_item[i];
4011 int new = old2new[old - 1] + 1;
4012 mesh->elem_node_item[i] = new;
4013 }
4014 /* import_item */
4015 for (i = 0; i < mesh->import_index[mesh->n_neighbor_pe]; i++) {
4016 int old = mesh->import_item[i];
4017 int new = old2new[old - 1] + 1;
4018 mesh->import_item[i] = new;
4019 }
4020 /* export_item */
4021 for (i = 0; i < mesh->export_index[mesh->n_neighbor_pe]; i++) {
4022 int old = mesh->export_item[i];
4023 int new = old2new[old - 1] + 1;
4024 mesh->export_item[i] = new;
4025 }
4026 /* mpc->mpc_item */
4027 for (i = 0; i < mesh->mpc->mpc_index[mesh->mpc->n_mpc]; i++) {
4028 int old = mesh->mpc->mpc_item[i];
4029 int new = old2new[old - 1] + 1;
4030 mesh->mpc->mpc_item[i] = new;
4031 }
4032 /* node_group->grp_item */
4033 for (i = 0; i < mesh->node_group->grp_index[mesh->node_group->n_grp]; i++) {
4034 int old = mesh->node_group->grp_item[i];
4035 int new = old2new[old - 1] + 1;
4036 mesh->node_group->grp_item[i] = new;
4037 }
4038 return HECMW_SUCCESS;
4039}
4040
4041static int delete_external_items(int n_elem, int n_grp, int *index, int *item,
4042 int n_memb) {
4043 int i, j, k;
4044 int n_del = 0;
4045 int start = 0;
4046 for (i = 0; i < n_grp; i++) {
4047 int end = index[i + 1];
4048 for (j = start; j < end; j++) {
4049 if (item[j * n_memb] > n_elem) {
4050 n_del++;
4051 } else {
4052 for (k = 0; k < n_memb; k++)
4053 item[(j - n_del) * n_memb + k] = item[j * n_memb + k];
4054 }
4055 }
4056 /* start for next group */
4057 start = end;
4058 index[i + 1] = end - n_del;
4059 }
4060 return HECMW_SUCCESS;
4061}
4062
4063static int delete_external_equations(int n_node, struct hecmwST_mpc *mpc) {
4064 int iold, inew, jold, jnew;
4065 int *active_eqn;
4066 int n_mpc_new;
4067 int *mpc_index_new, *mpc_item_new, *mpc_dof_new;
4068 double *mpc_val_new, *mpc_const_new;
4069 active_eqn = (int *) HECMW_malloc(sizeof(int) * mpc->n_mpc);
4070 if (active_eqn == NULL) {
4071 HECMW_set_error(errno, "");
4072 return HECMW_ERROR;
4073 }
4074 n_mpc_new = 0;
4075 for (iold = 0; iold < mpc->n_mpc; iold++) {
4076 int is_ext = 0;
4077 for (jold = mpc->mpc_index[iold]; jold < mpc->mpc_index[iold+1]; jold++) {
4078 if (mpc->mpc_item[jold] > n_node) {
4079 is_ext = 1;
4080 break;
4081 }
4082 }
4083 if (is_ext == 0) {
4084 active_eqn[n_mpc_new] = iold;
4085 n_mpc_new++;
4086 }
4087 }
4088 mpc_index_new = (int *) HECMW_malloc(sizeof(int) * (n_mpc_new + 1));
4089 if (mpc_index_new == NULL) {
4090 HECMW_set_error(errno, "");
4091 return HECMW_ERROR;
4092 }
4093 mpc_index_new[0] = 0;
4094 if (n_mpc_new == 0) {
4095 mpc_item_new = NULL;
4096 mpc_dof_new = NULL;
4097 mpc_val_new = NULL;
4098 mpc_const_new = NULL;
4099 } else {
4100 for (inew = 0; inew < n_mpc_new; inew++) {
4101 iold = active_eqn[inew];
4102 mpc_index_new[inew+1] = mpc_index_new[inew] + (mpc->mpc_index[iold+1] - mpc->mpc_index[iold]);
4103 }
4104 mpc_item_new = (int *) HECMW_malloc(sizeof(int) * mpc_index_new[n_mpc_new]);
4105 if (mpc_item_new == NULL) {
4106 HECMW_set_error(errno, "");
4107 return HECMW_ERROR;
4108 }
4109 mpc_dof_new = (int *) HECMW_malloc(sizeof(int) * mpc_index_new[n_mpc_new]);
4110 if (mpc_dof_new == NULL) {
4111 HECMW_set_error(errno, "");
4112 return HECMW_ERROR;
4113 }
4114 mpc_val_new = (double *) HECMW_malloc(sizeof(double) * mpc_index_new[n_mpc_new]);
4115 if (mpc_val_new == NULL) {
4116 HECMW_set_error(errno, "");
4117 return HECMW_ERROR;
4118 }
4119 mpc_const_new = (double *) HECMW_malloc(sizeof(double) * n_mpc_new);
4120 if (mpc_const_new == NULL) {
4121 HECMW_set_error(errno, "");
4122 return HECMW_ERROR;
4123 }
4124 for (inew = 0; inew < n_mpc_new; inew++) {
4125 iold = active_eqn[inew];
4126 for (jold = mpc->mpc_index[iold], jnew = mpc_index_new[inew]; jold < mpc->mpc_index[iold+1]; jold++, jnew++) {
4127 mpc_item_new[jnew] = mpc->mpc_item[jold];
4128 mpc_dof_new[jnew] = mpc->mpc_dof[jold];
4129 mpc_val_new[jnew] = mpc->mpc_val[jold];
4130 }
4131 mpc_const_new[inew] = mpc->mpc_const[iold];
4132 }
4133 }
4134 HECMW_free(active_eqn);
4135 mpc->n_mpc = n_mpc_new;
4136 HECMW_free(mpc->mpc_index);
4137 HECMW_free(mpc->mpc_item);
4138 HECMW_free(mpc->mpc_dof);
4139 HECMW_free(mpc->mpc_val);
4140 HECMW_free(mpc->mpc_const);
4141 mpc->mpc_index = mpc_index_new;
4142 mpc->mpc_item = mpc_item_new;
4143 mpc->mpc_dof = mpc_dof_new;
4144 mpc->mpc_val = mpc_val_new;
4145 mpc->mpc_const = mpc_const_new;
4146 return HECMW_SUCCESS;
4147}
4148
4149static int renumber_nodes(struct hecmwST_local_mesh *mesh) {
4151
4152 HECMW_log(HECMW_LOG_DEBUG, "Started renumbering nodes...\n");
4153
4154 if (reorder_nodes(mesh, mesh->node_old2new, mesh->node_new2old) !=
4155 HECMW_SUCCESS) {
4156 return HECMW_ERROR;
4157 }
4158
4159 if (delete_external_items(mesh->n_node, mesh->node_group->n_grp,
4161 1) != HECMW_SUCCESS) {
4162 return HECMW_ERROR;
4163 }
4164
4165 if (delete_external_equations(mesh->n_node, mesh->mpc) != HECMW_SUCCESS) {
4166 return HECMW_ERROR;
4167 }
4168
4169 HECMW_log(HECMW_LOG_DEBUG, "Finished renumbering nodes.\n");
4170 return HECMW_SUCCESS;
4171}
4172
4173static int renumber_back_nodes(struct hecmwST_local_mesh *mesh) {
4175
4176 HECMW_log(HECMW_LOG_DEBUG, "Started renumbering back nodes...\n");
4177
4178 if (reorder_nodes(mesh, mesh->node_new2old, mesh->node_old2new) !=
4179 HECMW_SUCCESS) {
4180 return HECMW_ERROR;
4181 }
4182
4183 HECMW_log(HECMW_LOG_DEBUG, "Finished renumbering back nodes.\n");
4184 return HECMW_SUCCESS;
4185}
4186
4187/*============================================================================*/
4188/* */
4189/* renumber elements */
4190/* */
4191/*============================================================================*/
4192static int renumber_elements_generate_tables(struct hecmwST_local_mesh *mesh) {
4193 int i;
4194 int count_rel, count_pex;
4195 int *old2new, *new2old;
4196
4198
4199 HECMW_log(HECMW_LOG_DEBUG, "Started generating renumbering tables...\n");
4200
4201 old2new = (int *)HECMW_malloc(sizeof(int) * mesh->n_elem_gross);
4202 if (old2new == NULL) {
4203 HECMW_set_error(errno, "");
4204 return HECMW_ERROR;
4205 }
4206 new2old = (int *)HECMW_malloc(sizeof(int) * mesh->n_elem_gross);
4207 if (new2old == NULL) {
4208 HECMW_set_error(errno, "");
4209 return HECMW_ERROR;
4210 }
4211 count_rel = 0;
4212 count_pex = mesh->n_elem;
4213 for (i = 0; i < mesh->n_elem_gross; i++) {
4214 int rank = mesh->elem_ID[i * 2 + 1];
4215 if (rank >= 0) {
4216 old2new[i] = count_rel;
4217 new2old[count_rel] = i;
4218 count_rel++;
4219 } else {
4220 old2new[i] = count_pex;
4221 new2old[count_pex] = i;
4222 count_pex++;
4223 }
4224 }
4225 mesh->elem_old2new = old2new;
4226 mesh->elem_new2old = new2old;
4227
4228 HECMW_log(HECMW_LOG_DEBUG, "Finished generating renumbering tables.\n");
4229 return HECMW_SUCCESS;
4230}
4231
4232static int reorder_elems(struct hecmwST_local_mesh *mesh, int *old2new,
4233 int *new2old) {
4234 int i;
4235
4236#ifndef NDEBUG
4237 for (i = 0; i < mesh->n_elem_gross; i++) {
4238 HECMW_assert(new2old[old2new[i]] == i);
4239 }
4240#endif
4241
4242 /*
4243 * Reorder using new2old
4244 */
4245 /* elem_ID */
4246 {
4247 int *new_elem_ID =
4248 (int *)HECMW_malloc(sizeof(int) * mesh->n_elem_gross * 2);
4249 if (new_elem_ID == NULL) {
4250 HECMW_set_error(errno, "");
4251 return HECMW_ERROR;
4252 }
4253 for (i = 0; i < mesh->n_elem_gross; i++) {
4254 int old = new2old[i];
4255 new_elem_ID[2 * i] = mesh->elem_ID[2 * old];
4256 new_elem_ID[2 * i + 1] = mesh->elem_ID[2 * old + 1];
4257 }
4259 mesh->elem_ID = new_elem_ID;
4260 }
4261 /* global_elem_ID */
4262 {
4263 int *new_global_elem_ID =
4264 (int *)HECMW_malloc(sizeof(int) * mesh->n_elem_gross);
4265 if (new_global_elem_ID == NULL) {
4266 HECMW_set_error(errno, "");
4267 return HECMW_ERROR;
4268 }
4269 for (i = 0; i < mesh->n_elem_gross; i++) {
4270 int old = new2old[i];
4271 new_global_elem_ID[i] = mesh->global_elem_ID[old];
4272 }
4274 mesh->global_elem_ID = new_global_elem_ID;
4275 }
4276 /* elem_type */
4277 {
4278 int *new_elem_type = (int *)HECMW_malloc(sizeof(int) * mesh->n_elem_gross);
4279 if (new_elem_type == NULL) {
4280 HECMW_set_error(errno, "");
4281 return HECMW_ERROR;
4282 }
4283 for (i = 0; i < mesh->n_elem_gross; i++) {
4284 int old = new2old[i];
4285 new_elem_type[i] = mesh->elem_type[old];
4286 }
4288 mesh->elem_type = new_elem_type;
4289 }
4290 /* elem_type_index, elem_type_item */
4291 {
4292 for (i = 0; i < mesh->n_elem_type; i++) {
4293 int j;
4294 for (j = mesh->elem_type_index[i]; j < mesh->n_elem; j++) {
4295 if (mesh->elem_type[j] != mesh->elem_type_item[i]) {
4296 break;
4297 }
4298 }
4299 mesh->elem_type_index[i + 1] = j;
4300 }
4301 }
4302 /* elem_node_index, elem_node_item */
4303 {
4304 int *new_elem_node_index =
4305 (int *)HECMW_malloc(sizeof(int) * (mesh->n_elem_gross + 1));
4306 int *new_elem_node_item = (int *)HECMW_malloc(
4307 sizeof(int) * mesh->elem_node_index[mesh->n_elem_gross]);
4308 if (new_elem_node_index == NULL || new_elem_node_item == NULL) {
4309 HECMW_set_error(errno, "");
4310 return HECMW_ERROR;
4311 }
4312 new_elem_node_index[0] = 0;
4313 for (i = 0; i < mesh->n_elem_gross; i++) {
4314 int old = new2old[i];
4315 int nn = mesh->elem_node_index[old + 1] - mesh->elem_node_index[old];
4316 new_elem_node_index[i + 1] = new_elem_node_index[i] + nn;
4317 }
4318 for (i = 0; i < mesh->n_elem_gross; i++) {
4319 int old = new2old[i];
4320 int start_old = mesh->elem_node_index[old];
4321 int start_new = new_elem_node_index[i];
4322 int nn = new_elem_node_index[i + 1] - start_new;
4323 int j;
4324 for (j = 0; j < nn; j++)
4325 new_elem_node_item[start_new + j] = mesh->elem_node_item[start_old + j];
4326 }
4329 mesh->elem_node_index = new_elem_node_index;
4330 mesh->elem_node_item = new_elem_node_item;
4331 }
4332 /* section_ID */
4333 {
4334 int *new_section_ID = (int *)HECMW_malloc(sizeof(int) * mesh->n_elem_gross);
4335 if (new_section_ID == NULL) {
4336 HECMW_set_error(errno, "");
4337 return HECMW_ERROR;
4338 }
4339 for (i = 0; i < mesh->n_elem_gross; i++) {
4340 int old = new2old[i];
4341 new_section_ID[i] = mesh->section_ID[old];
4342 }
4344 mesh->section_ID = new_section_ID;
4345 }
4346 /* elem_mat_ID_index, elem_mat_ID_item */
4347 {
4348 int *new_emID_index, *new_emID_item;
4349 new_emID_index =
4350 (int *)HECMW_malloc(sizeof(int) * (mesh->n_elem_gross + 1));
4351 if (new_emID_index == NULL) {
4352 HECMW_set_error(errno, "");
4353 return HECMW_ERROR;
4354 }
4355 new_emID_item = (int *)HECMW_malloc(
4356 sizeof(int) * mesh->elem_mat_ID_index[mesh->n_elem_gross]);
4357 if (new_emID_item == NULL) {
4358 HECMW_set_error(errno, "");
4359 return HECMW_ERROR;
4360 }
4361 new_emID_index[0] = 0;
4362 for (i = 0; i < mesh->n_elem_gross; i++) {
4363 int old, js_old, je_old, len, js_new, j;
4364 old = new2old[i];
4365 js_old = mesh->elem_mat_ID_index[old];
4366 je_old = mesh->elem_mat_ID_index[old + 1];
4367 len = je_old - js_old;
4368 js_new = new_emID_index[i];
4369 new_emID_index[i + 1] = js_new + len;
4370 for (j = 0; j < len; j++) {
4371 new_emID_item[js_new + j] = mesh->elem_mat_ID_item[js_old + j];
4372 }
4373 }
4376 mesh->elem_mat_ID_index = new_emID_index;
4377 mesh->elem_mat_ID_item = new_emID_item;
4378 }
4379
4380 /*
4381 * Update using old2new
4382 */
4383 /* elem_internal_list */
4384 for (i = 0; i < mesh->ne_internal; i++) {
4385 int old = mesh->elem_internal_list[i];
4386 int new = old2new[old - 1] + 1;
4387 mesh->elem_internal_list[i] = new;
4388 }
4389 /* shared_item */
4390 for (i = 0; i < mesh->shared_index[mesh->n_neighbor_pe]; i++) {
4391 int old = mesh->shared_item[i];
4392 int new = old2new[old - 1] + 1;
4393 mesh->shared_item[i] = new;
4394 }
4395 /* elem_groups->grp_item */
4396 for (i = 0; i < mesh->elem_group->grp_index[mesh->elem_group->n_grp]; i++) {
4397 int old = mesh->elem_group->grp_item[i];
4398 int new = old2new[old - 1] + 1;
4399 mesh->elem_group->grp_item[i] = new;
4400 }
4401 /* surf_groups->grp_item */
4402 for (i = 0; i < mesh->surf_group->grp_index[mesh->surf_group->n_grp]; i++) {
4403 int old = mesh->surf_group->grp_item[2 * i];
4404 int new = old2new[old - 1] + 1;
4405 mesh->surf_group->grp_item[2 * i] = new;
4406 }
4407 return HECMW_SUCCESS;
4408}
4409
4410static int renumber_elements(struct hecmwST_local_mesh *mesh) {
4411 if (mesh->n_elem_gross == mesh->n_elem) return HECMW_SUCCESS;
4412
4413 HECMW_log(HECMW_LOG_DEBUG, "Started renumbering elements...\n");
4414
4415 if (reorder_elems(mesh, mesh->elem_old2new, mesh->elem_new2old) !=
4416 HECMW_SUCCESS) {
4417 return HECMW_ERROR;
4418 }
4419
4420 /* if (delete_external_items( mesh->n_elem, mesh->n_neighbor_pe, */
4421 /* mesh->shared_index, mesh->shared_item, 1 ) != */
4422 /* HECMW_SUCCESS) { */
4423 /* return HECMW_ERROR; */
4424 /* } */
4425 if (delete_external_items(mesh->n_elem, mesh->elem_group->n_grp,
4427 1) != HECMW_SUCCESS) {
4428 return HECMW_ERROR;
4429 }
4430 if (delete_external_items(mesh->n_elem, mesh->surf_group->n_grp,
4432 2) != HECMW_SUCCESS) {
4433 return HECMW_ERROR;
4434 }
4435
4436 HECMW_log(HECMW_LOG_DEBUG, "Finished renumbering elements.\n");
4437 return HECMW_SUCCESS;
4438}
4439
4440static int renumber_back_elements(struct hecmwST_local_mesh *mesh) {
4441 if (mesh->n_elem_gross == mesh->n_elem) return HECMW_SUCCESS;
4442
4443 HECMW_log(HECMW_LOG_DEBUG, "Started renumbering back elements...\n");
4444
4445 if (reorder_elems(mesh, mesh->elem_new2old, mesh->elem_old2new) !=
4446 HECMW_SUCCESS) {
4447 return HECMW_ERROR;
4448 }
4449
4450 HECMW_log(HECMW_LOG_DEBUG, "Finished renumbering back elements.\n");
4451 return HECMW_SUCCESS;
4452}
4453
4454/*============================================================================*/
4455/* */
4456/* refine distributed mesh ONCE */
4457/* */
4458/*============================================================================*/
4459static int refine_dist_mesh(struct hecmwST_local_mesh *mesh,
4460 struct hecmwST_local_mesh *ref_mesh) {
4461 HECMW_log(HECMW_LOG_DEBUG, "Started refining mesh...\n");
4462
4463 if (copy_unchanging_info(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
4464 if (call_refiner(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
4465 if (rebuild_info(mesh, ref_mesh) != HECMW_SUCCESS) return HECMW_ERROR;
4466
4467 HECMW_log(HECMW_LOG_DEBUG, "Finished refining mesh.\n");
4468 return HECMW_SUCCESS;
4469}
4470
4471/*============================================================================*/
4472/* */
4473/* refine HEC-MW distributed mesh data */
4474/* */
4475/*============================================================================*/
4476int HECMW_dist_refine(struct hecmwST_local_mesh **mesh, int refine,
4477 const char *cad_filename, const char *part_filename) {
4478 int error_flag = HECMW_SUCCESS;
4479 int i;
4480
4481 if (refine <= 0) return HECMW_SUCCESS;
4482
4483 HECMW_log(HECMW_LOG_DEBUG, "Refinement requested; starting...\n");
4484
4485 if ((*mesh)->n_refine > 0) {
4486 if (renumber_back_elements(*mesh) != HECMW_SUCCESS) return HECMW_ERROR;
4487 if (renumber_back_nodes(*mesh) != HECMW_SUCCESS) return HECMW_ERROR;
4488 }
4489
4490 if (prepare_refiner(*mesh, cad_filename, part_filename) != HECMW_SUCCESS)
4491 return HECMW_ERROR;
4492
4493 for (i = 0; i < refine; i++) {
4494 struct hecmwST_local_mesh *ref_mesh = HECMW_dist_alloc();
4495 if (ref_mesh == NULL) {
4496 error_flag = HECMW_ERROR;
4497 break;
4498 }
4499
4500 if (i > 0) {
4501 clear_refiner();
4502 }
4503
4504 HECMW_log(HECMW_LOG_DEBUG, "Refining(%d)...\n", i + 1);
4505
4506 if (refine_dist_mesh(*mesh, ref_mesh) != HECMW_SUCCESS) {
4507 HECMW_log(HECMW_LOG_ERROR, "Refinement failed\n");
4508 error_flag = HECMW_ERROR;
4509 break;
4510 }
4511
4512 HECMW_dist_free(*mesh);
4513 *mesh = ref_mesh;
4514 }
4515
4516 if (terminate_refiner(*mesh) != HECMW_SUCCESS) return HECMW_ERROR;
4517
4518 if (error_flag != HECMW_SUCCESS) return HECMW_ERROR;
4519
4520 if (renumber_nodes(*mesh) != HECMW_SUCCESS) return HECMW_ERROR;
4521 if (renumber_elements(*mesh) != HECMW_SUCCESS) return HECMW_ERROR;
4522
4523 HECMW_log(HECMW_LOG_DEBUG, "Refinement finished.\n");
4524
4525 return HECMW_SUCCESS;
4526}
4527
4528#else /* HECMW_WITH_REFINER */
4529
4530int HECMW_dist_refine(struct hecmwST_local_mesh **mesh, int refine,
4531 const char *cad_filename, const char *part_filename) {
4532 if (refine > 0) {
4533 HECMW_log(HECMW_LOG_WARN, "Refiner not enabled; ignoring...\n");
4534 }
4535 return HECMW_SUCCESS;
4536}
4537#endif /* HECMW_WITH_REFINER */
int HECMW_bit_array_init(struct hecmw_bit_array *ba, size_t len)
int HECMW_bit_array_get(struct hecmw_bit_array *ba, size_t index)
void HECMW_bit_array_set(struct hecmw_bit_array *ba, size_t index)
void HECMW_bit_array_finalize(struct hecmw_bit_array *ba)
int HECMW_Isend(void *buffer, int count, HECMW_Datatype datatype, int dest, int tag, HECMW_Comm comm, HECMW_Request *request)
Definition: hecmw_comm.c:278
int HECMW_Waitall(int count, HECMW_Request *array_of_requests, HECMW_Status *array_of_statuses)
Definition: hecmw_comm.c:132
HECMW_Comm HECMW_comm_get_comm(void)
Definition: hecmw_comm.c:699
int HECMW_Allreduce(void *sendbuf, void *recvbuf, int count, HECMW_Datatype datatype, HECMW_Op op, HECMW_Comm comm)
Definition: hecmw_comm.c:364
int HECMW_comm_get_rank(void)
Definition: hecmw_comm.c:707
int HECMW_Recv(void *buffer, int count, HECMW_Datatype datatype, int source, int tag, HECMW_Comm comm, HECMW_Status *status)
Definition: hecmw_comm.c:235
#define HECMW_ETYPE_PRI1
#define HECMW_ETYPE_PRI2
#define HECMW_ETYPE_LN24
#define HECMW_ETYPE_MST1
#define HECMW_ETYPE_TET1_4
#define HECMW_ETYPE_BEM3
#define HECMW_ETYPE_LN54
#define HECMW_ETYPE_LN51
#define HECMW_ETYPE_LN15
#define HECMW_ETYPE_PYR1
#define HECMW_ETYPE_LN66
#define HECMW_ETYPE_LN65
#define HECMW_ETYPE_LN11
#define HECMW_ETYPE_LN43
#define HECMW_ETYPE_QUA2
#define HECMW_ETYPE_LN42
#define HECMW_ETYPE_LN32
#define HECMW_ETYPE_LN55
#define HECMW_ETYPE_MSQ2
#define HECMW_ETYPE_LN61
#define HECMW_ETYPE_JTQ2
#define HECMW_ETYPE_LN22
#define HECMW_ETYPE_SHQ1
#define HECMW_ETYPE_SHT1
#define HECMW_ETYPE_TET2
#define HECMW_ETYPE_PTQ1
#define HECMW_ETYPE_JTT2
#define HECMW_ETYPE_PTT1
#define HECMW_ETYPE_LN26
#define HECMW_ETYPE_LN33
#define HECMW_ETYPE_LN53
#define HECMW_ETYPE_LN64
#define HECMW_ETYPE_LN13
#define HECMW_ETYPE_LN46
#define HECMW_ETYPE_SHT2
#define HECMW_ETYPE_LN34
#define HECMW_ETYPE_TRI1
#define HECMW_ETYPE_SHT6
#define HECMW_ETYPE_LN41
#define HECMW_ETYPE_ROD2
#define HECMW_ETYPE_ROD31
#define HECMW_MAX_NODE_MAX
#define HECMW_ETYPE_HEX1
#define HECMW_ETYPE_JTT1
#define HECMW_ETYPE_ROD1
#define HECMW_ETYPE_LN62
#define HECMW_ETYPE_LN45
#define HECMW_ETYPE_JTQ1
#define HECMW_ETYPE_BEM1
#define HECMW_ETYPE_QUA1
#define HECMW_ETYPE_BEM2
#define HECMW_ETYPE_SHQ8
#define HECMW_ETYPE_LN35
#define HECMW_ETYPE_LN63
#define HECMW_ETYPE_LN25
#define HECMW_ETYPE_PYR2
#define HECMW_ETYPE_LN36
#define HECMW_ETYPE_TET1
#define HECMW_ETYPE_LN52
#define HECMW_ETYPE_LN16
#define HECMW_ETYPE_HEX1_4
#define HECMW_ETYPE_LN12
#define HECMW_ETYPE_MSQ1
#define HECMW_ETYPE_LN23
#define HECMW_ETYPE_HEX2
#define HECMW_ETYPE_TRI2
#define HECMW_ETYPE_PTT2
#define HECMW_ETYPE_LN31
#define HECMW_ETYPE_LN21
#define HECMW_ETYPE_PTQ2
#define HECMW_ETYPE_LN56
#define HECMW_ETYPE_LN44
#define HECMW_ETYPE_MST2
#define HECMW_ETYPE_SHQ2
#define HECMW_ETYPE_LN14
#define HECMW_INT
Definition: hecmw_config.h:48
MPI_Status HECMW_Status
Definition: hecmw_config.h:36
MPI_Comm HECMW_Comm
Definition: hecmw_config.h:30
MPI_Request HECMW_Request
Definition: hecmw_config.h:34
#define HECMW_SUM
Definition: hecmw_config.h:58
#define HECMW_ERROR
Definition: hecmw_config.h:66
#define HECMW_SUCCESS
Definition: hecmw_config.h:64
struct hecmwST_local_mesh * HECMW_dist_alloc()
void HECMW_dist_free(struct hecmwST_local_mesh *mesh)
int HECMW_dist_refine(struct hecmwST_local_mesh **mesh, int refine, const char *cad_filename, const char *part_filename)
struct hecmwST_local_mesh * mesh
Definition: hecmw_repart.h:71
int HECMW_set_error(int errorno, const char *fmt,...)
Definition: hecmw_error.c:37
int HECMW_get_max_node(int etype)
Definition: hecmw_etype.c:409
int HECMW_is_etype_33struct(int etype)
Definition: hecmw_etype.c:2007
int HECMW_is_etype_link(int etype)
Definition: hecmw_etype.c:1964
int HECMW_is_etype_surface(int etype)
Definition: hecmw_etype.c:1891
int HECMW_is_etype_truss(int etype)
Definition: hecmw_etype.c:2017
int HECMW_is_etype_patch(int etype)
Definition: hecmw_etype.c:2025
char * HECMW_get_ucd_label(int etype)
Definition: hecmw_etype.c:1216
int HECMW_is_etype_beam(int etype)
Definition: hecmw_etype.c:1940
int HECMW_is_etype_solid(int etype)
Definition: hecmw_etype.c:1903
int HECMW_is_etype_shell(int etype)
Definition: hecmw_etype.c:1950
int HECMW_is_etype_rod(int etype)
Definition: hecmw_etype.c:1882
#define NULL
int HECMW_log(int loglv, const char *fmt,...)
Definition: hecmw_log.c:260
#define HECMW_LOG_ERROR
Definition: hecmw_log.h:15
#define HECMW_LOG_WARN
Definition: hecmw_log.h:17
#define HECMW_LOG_DEBUG
Definition: hecmw_log.h:21
#define HECMW_calloc(nmemb, size)
Definition: hecmw_malloc.h:21
#define HECMW_free(ptr)
Definition: hecmw_malloc.h:24
#define HECMW_strdup(s)
Definition: hecmw_malloc.h:23
#define HECMW_malloc(size)
Definition: hecmw_malloc.h:20
struct result_list * elem_list
int nelem
int HECMW_set_int_iter_next(struct hecmw_set_int *set, int *value)
void HECMW_set_int_finalize(struct hecmw_set_int *set)
Definition: hecmw_set_int.c:36
int HECMW_set_int_add(struct hecmw_set_int *set, int value)
Definition: hecmw_set_int.c:60
int HECMW_set_int_init(struct hecmw_set_int *set)
Definition: hecmw_set_int.c:16
void HECMW_set_int_iter_init(struct hecmw_set_int *set)
size_t HECMW_set_int_check_dup(struct hecmw_set_int *set)
Definition: hecmw_set_int.c:78
#define HECMW_FLAG_PARTTYPE_UNKNOWN
Definition: hecmw_struct.h:142
#define HECMW_FLAG_PARTTYPE_NODEBASED
Definition: hecmw_struct.h:143
void HECMW_abort(HECMW_Comm comm)
Definition: hecmw_util.c:88
#define HECMW_assert(cond)
Definition: hecmw_util.h:40
int HECMW_varray_int_rmdup(struct hecmw_varray_int *varray)
size_t HECMW_varray_int_nval(const struct hecmw_varray_int *varray)
const int * HECMW_varray_int_get_cv(const struct hecmw_varray_int *varray)
int * HECMW_varray_int_get_v(struct hecmw_varray_int *varray)
int HECMW_varray_int_get(const struct hecmw_varray_int *varray, size_t index)
int HECMW_varray_int_init(struct hecmw_varray_int *varray)
void HECMW_varray_int_finalize(struct hecmw_varray_int *varray)
int HECMW_varray_int_append(struct hecmw_varray_int *varray, int value)
int HECMW_varray_int_cat(struct hecmw_varray_int *varray, const struct hecmw_varray_int *varray2)
int HECMW_varray_int_resize(struct hecmw_varray_int *varray, size_t len)
int * amp_type_definition
Definition: hecmw_struct.h:61
double * amp_table
Definition: hecmw_struct.h:72
struct hecmwST_section * section
Definition: hecmw_struct.h:244
struct hecmwST_amplitude * amp
Definition: hecmw_struct.h:247
struct hecmwST_material * material
Definition: hecmw_struct.h:245
struct hecmwST_refine_origin * refine_origin
Definition: hecmw_struct.h:252
struct hecmwST_mpc * mpc
Definition: hecmw_struct.h:246
struct hecmwST_node_grp * node_group
Definition: hecmw_struct.h:248
double * node_init_val_item
Definition: hecmw_struct.h:180
struct hecmwST_contact_pair * contact_pair
Definition: hecmw_struct.h:251
struct hecmwST_surf_grp * surf_group
Definition: hecmw_struct.h:250
char gridfile[HECMW_FILENAME_LEN+1]
Definition: hecmw_struct.h:153
char header[HECMW_HEADER_LEN+1]
Definition: hecmw_struct.h:156
HECMW_Comm HECMW_COMM
Definition: hecmw_struct.h:208
struct hecmwST_elem_grp * elem_group
Definition: hecmw_struct.h:249
int * when_i_was_refined_node
Definition: hecmw_struct.h:226
int * when_i_was_refined_elem
Definition: hecmw_struct.h:227
int * mat_subitem_index
Definition: hecmw_struct.h:42
double * mat_val
Definition: hecmw_struct.h:44
double * mat_temp
Definition: hecmw_struct.h:45
int * mpc_dof
Definition: hecmw_struct.h:52
double * mpc_val
Definition: hecmw_struct.h:53
double * mpc_const
Definition: hecmw_struct.h:54
int * mpc_index
Definition: hecmw_struct.h:50
int * mpc_item
Definition: hecmw_struct.h:51
double * sect_R_item
Definition: hecmw_struct.h:32
int * sect_mat_ID_index
Definition: hecmw_struct.h:27
int * sect_mat_ID_item
Definition: hecmw_struct.h:28