00001 #include "siftop.h"
00002 #include "galias.h"
00003
00004
00005 #include <stdlib.h>
00006 #include <string.h>
00007
00008 const gtypel operator ++(gtypel &a, int)
00009 {
00010 gtypel old;
00011 switch (a)
00012 {
00013 case isolinear1d:
00014 old = a;
00015 a = isoquadratic1d;
00016 return old;
00017 case isoquadratic1d:
00018 old = a;
00019 a = trianglelinear;
00020 return old;
00021 case trianglelinear:
00022 old = a;
00023 a = trianglequadratic;
00024 return old;
00025 case trianglequadratic:
00026 old = a;
00027 a = isolinear2d;
00028 return old;
00029 case isolinear2d:
00030 old = a;
00031 a = isoquadratic2d;
00032 return old;
00033 case isoquadratic2d:
00034 old = a;
00035 a = tetrahedronlinear;
00036 return old;
00037 case tetrahedronlinear:
00038 old = a;
00039 a = tetrahedronquadratic;
00040 return old;
00041 case tetrahedronquadratic:
00042 old = a;
00043 a = pyramidelinear;
00044 return old;
00045 case pyramidelinear:
00046 old = a;
00047 a = pyramidequadratic;
00048 return old;
00049 case pyramidequadratic:
00050 old = a;
00051 a = vedgelinear;
00052 return old;
00053 case vedgelinear:
00054 old = a;
00055 a = vedgequadratic;
00056 return old;
00057 case vedgequadratic:
00058 old = a;
00059 a = isolinear3d;
00060 return old;
00061 case isolinear3d:
00062 old = a;
00063 a = isoquadratic3d;
00064 return old;
00065 case isoquadratic3d:
00066 old = a;
00067 a = all_element_types;
00068 return old;
00069 default:
00070 print_err("increasing of invalid gtypel value is required", __FILE__, __LINE__, __func__);
00071 }
00072 return all_element_types;
00073 }
00074
00075
00076
00077
00078
00079
00080
00081
00082 snode::snode()
00083 {
00084 x = y = z = 0.0;
00085 nprop = 0;
00086 entid = NULL;
00087 prop = NULL;
00088 }
00089
00090
00091
00092
00093
00094
00095
00096
00097 snode::~snode()
00098 {
00099 delete [] entid;
00100 delete [] prop;
00101 }
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111 void snode::alloc(long num_prop)
00112 {
00113 nprop = num_prop;
00114 entid = new entityp[nprop];
00115 prop = new long[nprop];
00116 }
00117
00118
00119
00120
00121
00122
00123
00124
00125
00126
00127
00128
00129 void snode::add_prop(entityp ent, long propid)
00130 {
00131 entityp *auxent = entid;
00132 long *auxprop = prop;
00133 long n = nprop;
00134 long i;
00135
00136 for(i=0; i<nprop; i++)
00137 {
00138 if ((entid[i] == ent) && (prop[i] == propid))
00139 return;
00140 }
00141
00142 nprop += 1;
00143 entid = new entityp[nprop];
00144 prop = new long[nprop];
00145 for(i=0; i<n; i++)
00146 {
00147 entid[i] = auxent[i];
00148 prop[i] = auxprop[i];
00149 }
00150 entid[n] = ent;
00151 prop[n] = propid;
00152 delete [] auxent;
00153 delete [] auxprop;
00154 }
00155
00156
00157
00158
00159
00160
00161
00162
00163
00164
00165
00166
00167
00168
00169
00170
00171
00172
00173
00174 long snode::searchprop(long aprop, gentity ent, long nid, sedges *edg, sfaces *sfc)
00175 {
00176 long i, j;
00177 for(i=0; i<nprop; i++)
00178 {
00179 switch(ent)
00180 {
00181 case gvertex:
00182 case gcurve:
00183 case gregion:
00184 if ((int(ent) == int(entid[i])) && (prop[i] == aprop)) return 1;
00185 break;
00186 case gsurface:
00187 if (((entid[i] == esurface) || (entid[i] == epatch) ||
00188 (entid[i] == eshell)) && (prop[i] == aprop)) return 1;
00189 break;
00190 default:
00191 print_err("Unknown type of entity is required", __FILE__, __LINE__, __func__);
00192 abort();
00193 }
00194 }
00195
00196 if ((ent == gcurve) && (edg))
00197 {
00198
00199 for(i=0; i<edg->nprop; i++)
00200 {
00201 if (edg->proplist[i] == aprop)
00202 {
00203 for(j=0; j<edg->prop_list_nod_num[i]; j++)
00204 {
00205 if (edg->prop_list_nod[i][j] == nid) return 1;
00206 if (edg->prop_list_nod[i][j] > nid) break;
00207 }
00208 }
00209 }
00210 }
00211
00212 if ((ent == gsurface) && (sfc))
00213 {
00214
00215 for(i=0; i<sfc->nprop; i++)
00216 {
00217 if (sfc->proplist[i] == aprop)
00218 {
00219 for(j=0; j<sfc->prop_list_nod_num[i]; j++)
00220 {
00221 if (sfc->prop_list_nod[i][j] == nid) return 1;
00222 if (sfc->prop_list_nod[i][j] > nid) break;
00223 }
00224 }
00225 }
00226 }
00227 return 0;
00228 }
00229
00230
00231
00232
00233
00234
00235
00236
00237 selement::selement()
00238 {
00239 type = gtypel(0);
00240 nne = 0L;
00241
00242 ned = tnned = 0L;
00243 nned = NULL;
00244 propedg = NULL;
00245 edgenod = NULL;
00246
00247 nsurf = tnnsurf = 0;
00248 nnsurf = NULL;
00249 surfnod = NULL;
00250 propsurf = NULL;
00251
00252 nodes = NULL;
00253 prop = 0L;
00254 }
00255
00256
00257
00258
00259
00260
00261
00262
00263 selement::~selement()
00264 {
00265 long i;
00266 delete [] nned;
00267 delete [] nnsurf;
00268 delete [] nodes;
00269 delete [] propedg;
00270 delete [] propsurf;
00271 for(i=0; i<ned; i++)
00272 delete [] edgenod[i];
00273 for(i=0; i<nsurf; i++)
00274 delete [] surfnod[i];
00275 delete [] edgenod;
00276 delete [] surfnod;
00277 }
00278
00279
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294 long selement::alloc(long alloc_bp)
00295 {
00296 long i;
00297 switch (type)
00298 {
00299 case isolinear1d:
00300 nne = 2;
00301 ned = 1;
00302 nned = new long[ned];
00303 edgenod = new long*[ned];
00304 for(i=0; i<ned; i++)
00305 {
00306 nned[i] = 2;
00307 edgenod[i] = new long[nned[i]];
00308 }
00309 edgenod[0][0] = 0;
00310 edgenod[0][1] = 1;
00311 nsurf = 0;
00312 tnned = 2;
00313 tnnsurf = 0;
00314 break;
00315 case isoquadratic1d:
00316 nne = 3;
00317 ned = 1;
00318 nned = new long[ned];
00319 edgenod = new long*[ned];
00320 for(i=0; i<ned; i++)
00321 {
00322 nned[i] = 3;
00323 edgenod[i] = new long[nned[i]];
00324 }
00325 edgenod[0][0] = 0;
00326 edgenod[0][1] = 1;
00327 edgenod[0][2] = 2;
00328 nsurf = 0;
00329 tnned = 3;
00330 tnnsurf = 0;
00331 break;
00332 case trianglelinear:
00333 nne = 3;
00334 ned = 3;
00335 nsurf = 1;
00336 nned = new long[ned];
00337 edgenod = new long*[ned];
00338 surfnod = new long*[nsurf];
00339 for(i=0; i<ned; i++)
00340 {
00341 nned[i] = 2;
00342 edgenod[i] = new long[nned[i]];
00343 }
00344 edgenod[0][0] = 0; edgenod[0][1] = 1;
00345 edgenod[1][0] = 1; edgenod[1][1] = 2;
00346 edgenod[2][0] = 2; edgenod[2][1] = 0;
00347 nnsurf = new long[1];
00348 for(i=0; i<nsurf; i++)
00349 {
00350 nnsurf[i] = 3;
00351 surfnod[i] = new long[nnsurf[i]];
00352 }
00353 surfnod[0][0] = 0; surfnod[0][1] = 1; surfnod[0][2] = 2;
00354 tnned = 6;
00355 tnnsurf = 3;
00356 break;
00357 case trianglequadratic:
00358 nne = 6;
00359 ned = 3;
00360 nsurf = 1;
00361 nned = new long[ned];
00362 edgenod = new long*[ned];
00363 surfnod = new long*[nsurf];
00364 for(i=0; i<ned; i++)
00365 {
00366 nned[i] = 3;
00367 edgenod[i] = new long[nned[i]];
00368 }
00369 edgenod[0][0] = 0; edgenod[0][1] = 1; edgenod[0][2] = 3;
00370 edgenod[1][0] = 1; edgenod[1][1] = 2; edgenod[1][2] = 4;
00371 edgenod[2][0] = 2; edgenod[2][1] = 0; edgenod[2][2] = 5;
00372 nnsurf = new long[nsurf];
00373 for(i=0; i<nsurf; i++)
00374 {
00375 nnsurf[i] = 6;
00376 surfnod[i] = new long[nnsurf[i]];
00377 }
00378 surfnod[0][0] = 0; surfnod[0][1] = 1; surfnod[0][2] = 2;
00379 surfnod[0][3] = 3; surfnod[0][4] = 4; surfnod[0][5] = 5;
00380 tnned = 9;
00381 tnnsurf = 6;
00382 break;
00383 case isolinear2d:
00384 nne = 4;
00385 ned = 4;
00386 nsurf = 1;
00387 nned = new long[ned];
00388 edgenod = new long*[ned];
00389 surfnod = new long*[nsurf];
00390 for(i=0; i<ned; i++)
00391 {
00392 nned[i] = 2;
00393 edgenod[i] = new long[nned[i]];
00394 }
00395 edgenod[0][0] = 0; edgenod[0][1] = 1;
00396 edgenod[1][0] = 1; edgenod[1][1] = 2;
00397 edgenod[2][0] = 2; edgenod[2][1] = 3;
00398 edgenod[3][0] = 3; edgenod[3][1] = 0;
00399 nnsurf = new long[nsurf];
00400 for(i=0; i<nsurf; i++)
00401 {
00402 nnsurf[i] = 4;
00403 surfnod[i] = new long[nnsurf[i]];
00404 }
00405 surfnod[0][0] = 0; surfnod[0][1] = 1; surfnod[0][2] = 2; surfnod[0][3] = 3;
00406 tnned = 8;
00407 tnnsurf = 4;
00408 break;
00409 case isoquadratic2d:
00410 nne = 8;
00411 ned = 4;
00412 nsurf = 1;
00413 nned = new long[ned];
00414 edgenod = new long*[ned];
00415 surfnod = new long*[nsurf];
00416 for(i=0; i<ned; i++)
00417 {
00418 nned[i] = 3;
00419 edgenod[i] = new long[nned[i]];
00420 }
00421 edgenod[0][0] = 0; edgenod[0][1] = 1; edgenod[0][2] = 4;
00422 edgenod[1][0] = 1; edgenod[1][1] = 2; edgenod[1][2] = 5;
00423 edgenod[2][0] = 2; edgenod[2][1] = 3; edgenod[2][2] = 6;
00424 edgenod[3][0] = 3; edgenod[3][1] = 0; edgenod[3][2] = 7;
00425 nnsurf = new long[nsurf];
00426 for(i=0; i<nsurf; i++)
00427 {
00428 nnsurf[i] = 8;
00429 surfnod[i] = new long[nnsurf[i]];
00430 }
00431 surfnod[0][0] = 0; surfnod[0][1] = 1; surfnod[0][2] = 2; surfnod[0][3] = 3;
00432 surfnod[0][4] = 4; surfnod[0][5] = 5; surfnod[0][6] = 6; surfnod[0][7] = 7;
00433 tnned = 12;
00434 tnnsurf = 8;
00435 break;
00436 case tetrahedronlinear:
00437 nne = 4;
00438 ned = 6;
00439 nsurf = 4;
00440 nned = new long[ned];
00441 edgenod = new long*[ned];
00442 surfnod = new long*[nsurf];
00443 for(i=0; i<ned; i++)
00444 {
00445 nned[i] = 2;
00446 edgenod[i] = new long[nned[i]];
00447 }
00448 edgenod[0][0]=0; edgenod[0][1]=1;
00449 edgenod[1][0]=1; edgenod[1][1]=2;
00450 edgenod[2][0]=2; edgenod[2][1]=0;
00451 edgenod[3][0]=1; edgenod[3][1]=3;
00452 edgenod[4][0]=3; edgenod[4][1]=2;
00453 edgenod[5][0]=3; edgenod[5][1]=0;
00454 nnsurf = new long[nsurf];
00455 for(i=0; i<nsurf; i++)
00456 {
00457 nnsurf[i] = 3;
00458 surfnod[i] = new long[nnsurf[i]];
00459 }
00460 surfnod[0][0]=2; surfnod[0][1]=1; surfnod[0][2]=3;
00461 surfnod[1][0]=0; surfnod[1][1]=2; surfnod[1][2]=3;
00462 surfnod[2][0]=1; surfnod[2][1]=0; surfnod[2][2]=3;
00463 surfnod[3][0]=0; surfnod[3][1]=1; surfnod[3][2]=2;
00464 tnned = 12;
00465 tnnsurf = 12;
00466 break;
00467 case tetrahedronquadratic:
00468 nne = 10;
00469 ned = 6;
00470 nsurf = 4;
00471 nned = new long[ned];
00472 edgenod = new long*[ned];
00473 surfnod = new long*[nsurf];
00474 for(i=0; i<ned; i++)
00475 {
00476 nned[i] = 3;
00477 edgenod[i] = new long[nned[i]];
00478 }
00479 edgenod[0][0]=0; edgenod[0][1]=1; edgenod[0][2]=4;
00480 edgenod[1][0]=1; edgenod[1][1]=2; edgenod[1][2]=5;
00481 edgenod[2][0]=2; edgenod[2][1]=0; edgenod[2][2]=6;
00482 edgenod[3][0]=1; edgenod[3][1]=3; edgenod[3][2]=8;
00483 edgenod[4][0]=3; edgenod[4][1]=2; edgenod[4][2]=9;
00484 edgenod[5][0]=3; edgenod[5][1]=0; edgenod[5][2]=7;
00485 nnsurf = new long[nsurf];
00486 for(i=0; i<nsurf; i++)
00487 {
00488 nnsurf[i] = 6;
00489 surfnod[i] = new long[nnsurf[i]];
00490 }
00491 surfnod[0][0]=2; surfnod[0][1]=1; surfnod[0][2]=3;
00492 surfnod[0][3]=5; surfnod[0][4]=8; surfnod[0][5]=9;
00493 surfnod[1][0]=0; surfnod[1][1]=2; surfnod[1][2]=3;
00494 surfnod[1][3]=6; surfnod[1][4]=9; surfnod[1][5]=7;
00495 surfnod[2][0]=1; surfnod[2][1]=0; surfnod[2][2]=3;
00496 surfnod[2][3]=4; surfnod[2][4]=7; surfnod[2][5]=8;
00497 surfnod[3][0]=0; surfnod[3][1]=1; surfnod[3][2]=2;
00498 surfnod[3][3]=4; surfnod[3][4]=5; surfnod[3][5]=6;
00499 tnned = 18;
00500 tnnsurf = 24;
00501 break;
00502 case pyramidelinear:
00503 nne = 5;
00504 ned = 8;
00505 nsurf = 5;
00506 nned = new long[ned];
00507 edgenod = new long*[ned];
00508 surfnod = new long*[nsurf];
00509 for(i=0; i<ned; i++)
00510 {
00511 nned[i] = 2;
00512 edgenod[i] = new long[nned[i]];
00513 }
00514 nnsurf = new long[nsurf];
00515 for(i=0; i<nsurf-1; i++)
00516 {
00517 nnsurf[i] = 3;
00518 surfnod[i] = new long[nnsurf[i]];
00519 }
00520 nnsurf[nsurf-1] = 4;
00521 surfnod[nsurf-1] = new long[nnsurf[nsurf-1]];
00522 tnned = 16;
00523 tnnsurf = 16;
00524 break;
00525 case pyramidequadratic:
00526 nne = 13;
00527 ned = 9;
00528 nsurf = 5;
00529 nned = new long[ned];
00530 edgenod = new long*[ned];
00531 surfnod = new long*[nsurf];
00532 for(i=0; i<ned; i++)
00533 {
00534 nned[i] = 3;
00535 edgenod[i] = new long[nned[i]];
00536 }
00537 nnsurf = new long[nsurf];
00538 for(i=0; i<nsurf-1; i++)
00539 {
00540 nnsurf[i] = 6;
00541 surfnod[i] = new long[nnsurf[i]];
00542 }
00543 nnsurf[nsurf-1] = 8;
00544 surfnod[nsurf-1] = new long[nnsurf[nsurf-1]];
00545 tnned = 27;
00546 tnnsurf = 32;
00547 break;
00548 case vedgelinear:
00549 nne = 6;
00550 ned = 9;
00551 nsurf = 5;
00552 nned = new long[ned];
00553 edgenod = new long*[ned];
00554 surfnod = new long*[nsurf];
00555 for(i=0; i<ned; i++)
00556 {
00557 nned[i] = 2;
00558 edgenod[i] = new long[nned[i]];
00559 }
00560 nnsurf = new long[nsurf];
00561 for(i=0; i<nsurf-2; i++)
00562 {
00563 nnsurf[i] = 4;
00564 surfnod[i] = new long[nnsurf[i]];
00565 }
00566 nnsurf[nsurf-2] = 3;
00567 nnsurf[nsurf-1] = 3;
00568 surfnod[nsurf-2] = new long[nnsurf[nsurf-2]];
00569 surfnod[nsurf-1] = new long[nnsurf[nsurf-1]];
00570 tnned = 18;
00571 tnnsurf = 18;
00572 break;
00573 case vedgequadratic:
00574 nne = 15;
00575 ned = 9;
00576 nsurf = 5;
00577 nned = new long[ned];
00578 edgenod = new long*[ned];
00579 surfnod = new long*[nsurf];
00580 for(i=0; i<ned; i++)
00581 {
00582 nned[i] = 3;
00583 edgenod[i] = new long[nned[i]];
00584 }
00585 nnsurf = new long[nsurf];
00586 for(i=0; i<nsurf-2; i++)
00587 {
00588 nnsurf[i] = 8;
00589 surfnod[i] = new long[nnsurf[i]];
00590 }
00591 nnsurf[nsurf-2] = 6;
00592 nnsurf[nsurf-1] = 6;
00593 surfnod[nsurf-2] = new long[nnsurf[nsurf-2]];
00594 surfnod[nsurf-1] = new long[nnsurf[nsurf-1]];
00595 tnned = 27;
00596 tnnsurf = 36;
00597 break;
00598 case isolinear3d:
00599 nne = 8;
00600 ned = 12;
00601 nsurf = 6;
00602 nned = new long[ned];
00603 edgenod = new long*[ned];
00604 surfnod = new long*[nsurf];
00605 for(i=0; i<ned; i++)
00606 {
00607 nned[i] = 2;
00608 edgenod[i] = new long[nned[i]];
00609 }
00610 edgenod[0][0] =0; edgenod[0][1] =1;
00611 edgenod[1][0] =1; edgenod[1][1] =2;
00612 edgenod[2][0] =2; edgenod[2][1] =3;
00613 edgenod[3][0] =3; edgenod[3][1] =0;
00614 edgenod[4][0] =0; edgenod[4][1] =4;
00615 edgenod[5][0] =1; edgenod[5][1] =5;
00616 edgenod[6][0] =2; edgenod[6][1] =6;
00617 edgenod[7][0] =3; edgenod[7][1] =7;
00618 edgenod[8][0] =4; edgenod[8][1] =5;
00619 edgenod[9][0] =5; edgenod[9][1] =6;
00620 edgenod[10][0]=6; edgenod[10][1]=7;
00621 edgenod[11][0]=7; edgenod[11][1]=4;
00622 nnsurf = new long[nsurf];
00623 for(i=0; i<nsurf; i++)
00624 {
00625 nnsurf[i] = 4;
00626 surfnod[i] = new long[nnsurf[i]];
00627 linhexahedral_surfnod (surfnod[i], i);
00628 }
00629 tnned = 24;
00630 tnnsurf = 24;
00631 break;
00632 case isoquadratic3d:
00633 nne = 20;
00634 ned = 12;
00635 nsurf = 6;
00636 nned = new long[ned];
00637 edgenod = new long*[ned];
00638 surfnod = new long*[nsurf];
00639 for(i=0; i<ned; i++)
00640 {
00641 nned[i] = 3;
00642 edgenod[i] = new long[nned[i]];
00643 }
00644 edgenod[0][0] =0; edgenod[0][1] =1; edgenod[0][2] = 8;
00645 edgenod[1][0] =1; edgenod[1][1] =2; edgenod[1][2] = 9;
00646 edgenod[2][0] =2; edgenod[2][1] =3; edgenod[2][2] = 10;
00647 edgenod[3][0] =3; edgenod[3][1] =0; edgenod[3][2] = 11;
00648 edgenod[4][0] =0; edgenod[4][1] =4; edgenod[4][2] = 12;
00649 edgenod[5][0] =1; edgenod[5][1] =5; edgenod[5][2] = 13;
00650 edgenod[6][0] =2; edgenod[6][1] =6; edgenod[6][2] = 14;
00651 edgenod[7][0] =3; edgenod[7][1] =7; edgenod[7][2] = 15;
00652 edgenod[8][0] =4; edgenod[8][1] =5; edgenod[8][2] = 16;
00653 edgenod[9][0] =5; edgenod[9][1] =6; edgenod[9][2] = 17;
00654 edgenod[10][0]=6; edgenod[10][1]=7; edgenod[10][2]= 18;
00655 edgenod[11][0]=7; edgenod[11][1]=4; edgenod[11][2]= 19;
00656 nnsurf = new long[nsurf];
00657 for(i=0; i<nsurf; i++)
00658 {
00659 nnsurf[i] = 8;
00660 surfnod[i] = new long[nnsurf[i]];
00661 quadhexahedral_surfnod (surfnod[i],i);
00662 }
00663
00664 tnned = 36;
00665 tnnsurf = 48;
00666 break;
00667 default:
00668 return 1;
00669 }
00670 nodes = NULL;
00671 propedg = NULL;
00672 propsurf = NULL;
00673 if (nne)
00674 {
00675 nodes = new long[nne];
00676 memset(nodes, 0, sizeof(*nodes)*nne);
00677 }
00678 if (ned && alloc_bp)
00679 {
00680 propedg = new long[ned];
00681 memset(propedg, 0, sizeof(*propedg)*ned);
00682 }
00683 if (nsurf && alloc_bp)
00684 {
00685 propsurf = new long[nsurf];
00686 memset(propsurf, 0, sizeof(*propsurf)*nsurf);
00687 }
00688 return 0;
00689 }
00690
00691
00692
00693
00694
00695
00696
00697
00698
00699
00700
00701
00702
00703
00704
00705
00706
00707
00708
00709
00710
00711
00712
00713
00714 long selement::searchprop(long aprop, gentity ent, snode *tnodes, long eid, sedges *edg, sfaces *sfc, long *&entid, long &nentid)
00715 {
00716 long i, j, ret = 0;
00717 entid = NULL;
00718 nentid = 0;
00719
00720 switch(ent)
00721 {
00722 case gvertex:
00723 for (i=0; i<nne; i++)
00724 {
00725 ret = tnodes[nodes[i]].searchprop(aprop, ent, nodes[i], edg, sfc);
00726 if (ret)
00727 return 1;
00728 }
00729 break;
00730 case gcurve:
00731 if (propedg)
00732 {
00733 for (i=0; i<ned; i++)
00734 {
00735 if (propedg[i] == aprop)
00736 return 1;
00737 }
00738 }
00739 else
00740 ret=-1;
00741 break;
00742 case gsurface:
00743 if (propsurf)
00744 {
00745 for (i=0; i<nsurf; i++)
00746 {
00747 if (propsurf[i] == aprop)
00748 return 1;
00749 }
00750 }
00751 else
00752 ret=-1;
00753 break;
00754 case gregion:
00755 if (prop == aprop) return 1;
00756 break;
00757 default:
00758 print_err("Unknown type of entity is required", __FILE__, __LINE__, __func__);
00759 abort();
00760 }
00761
00762 if ((ent == gcurve) && (edg))
00763 {
00764 ret = 0;
00765
00766 for(i=0; i<edg->nprop; i++)
00767 {
00768 if (edg->proplist[i] == aprop)
00769 {
00770 for(j=0; j<edg->prop_list_elem_num[i]; j++)
00771 {
00772 if (edg->prop_list_elem[i][j] == eid)
00773 {
00774 entid = edg->prop_list_elem_edgid[i][j];
00775 nentid = edg->prop_list_elem_edgid_num[i][j];
00776 return 1;
00777 }
00778 if (edg->prop_list_elem[i][j] > eid) break;
00779 }
00780 }
00781 }
00782 }
00783
00784 if ((ent == gsurface) && (sfc))
00785 {
00786 ret = 0;
00787
00788 for(i=0; i<sfc->nprop; i++)
00789 {
00790 if (sfc->proplist[i] == aprop)
00791 {
00792 for(j=0; j<sfc->prop_list_elem_num[i]; j++)
00793 {
00794 if (sfc->prop_list_elem[i][j] == eid)
00795 {
00796 entid = sfc->prop_list_elem_sfid[i][j];
00797 nentid = sfc->prop_list_elem_sfid_num[i][j];
00798 return 1;
00799 }
00800 if (sfc->prop_list_elem[i][j] > eid) break;
00801 }
00802 }
00803 }
00804 }
00805 return ret;
00806 }
00807
00808
00809
00810
00811
00812
00813
00814
00815
00816
00817
00818
00819
00820
00821
00822
00823
00824
00825
00826
00827
00828 long selement::searchnodprop(long aprop, gentity ent, snode *tnodes, sedges *edg, sfaces *sfc, long *propent)
00829 {
00830 long i, aux, ret = 0;
00831 for (i=0; i<nne; i++)
00832 {
00833 propent[i] = -1;
00834
00835 aux = tnodes[nodes[i]].searchprop(aprop, ent, nodes[i], edg, sfc);
00836 if (aux)
00837 {
00838 propent[i] = aprop;
00839 ret = 1;
00840 }
00841 }
00842 return ret;
00843 }
00844
00845
00846
00847
00848
00849
00850
00851
00852
00853
00854
00855
00856
00857
00858
00859
00860
00861 long selement::compare_edg(long *nod, long n)
00862 {
00863 long i, j, k;
00864 long nf;
00865 long stop;
00866
00867 for (i=0; i<ned; i++)
00868 {
00869 nf = 0;
00870 for(j=0; j<n; j++)
00871 {
00872 stop = 1;
00873 for(k=0; k<nned[i]; k++)
00874 {
00875 if (nodes[edgenod[i][k]] == nod[j])
00876 {
00877 nf++;
00878 stop = 0;
00879 break;
00880 }
00881 }
00882 if (stop)
00883 break;
00884 }
00885 if (nf == n)
00886 return i;
00887 }
00888 return -1;
00889 }
00890
00891
00892
00893
00894
00895
00896
00897
00898
00899
00900
00901
00902
00903
00904
00905
00906
00907 long selement::compare_surf(long *nod, long n)
00908 {
00909 long i, j, k;
00910 long nf;
00911 long stop;
00912
00913 for (i=0; i<nsurf; i++)
00914 {
00915 nf = 0;
00916 for(j=0; j<n; j++)
00917 {
00918 stop = 1;
00919 for(k=0; k<nnsurf[i]; k++)
00920 {
00921 if (nodes[surfnod[i][k]] == nod[j])
00922 {
00923 nf++;
00924 stop = 0;
00925 break;
00926 }
00927 }
00928 if (stop)
00929 break;
00930 }
00931 if (nf == n)
00932 return i;
00933 }
00934 return -1;
00935 }
00936
00937
00938
00939
00940
00941
00942
00943
00944 sedge::sedge(void)
00945 {
00946 n1 = n2 = -1;
00947 prop = -1;
00948 }
00949
00950
00951
00952
00953
00954
00955
00956
00957 sedge::~sedge(void)
00958 {
00959
00960 }
00961
00962
00963
00964
00965
00966
00967
00968
00969
00970
00971
00972
00973
00974 long sedge::read(XFILE *in)
00975 {
00976 if (xfscanf(in, "%ld %ld %ld", &n1, &n2, &prop) < 3)
00977 return 1;
00978
00979 n1--;
00980 n2--;
00981 return 0;
00982 }
00983
00984
00985
00986
00987
00988
00989
00990
00991
00992
00993
00994
00995
00996
00997
00998
00999
01000
01001
01002 void sedge::print(FILE *out, long *lnn)
01003 {
01004 if (lnn)
01005 fprintf(out, "%ld %ld %ld\n", lnn[n1], lnn[n2], prop);
01006 else
01007 fprintf(out, "%ld %ld %ld\n", n1+1, n2+1, prop);
01008 }
01009
01010
01011
01012
01013
01014
01015
01016
01017
01018
01019
01020
01021
01022
01023
01024
01025
01026
01027 long sedge::cordomtest(long *lnn)
01028 {
01029 if (lnn[n1] && lnn[n2])
01030 return 1;
01031
01032 return 0;
01033 }
01034
01035
01036
01037
01038
01039
01040
01041
01042 sedges::sedges()
01043 {
01044 nedg = 0L;
01045 nprop = 0L;
01046 edges = NULL;
01047 proplist = NULL;
01048 prop_list_nod = NULL;
01049 prop_list_nod_num = NULL;
01050 prop_list_elem = NULL;
01051 prop_list_elem_num = NULL;
01052 prop_list_elem_edgid = NULL;
01053 prop_list_elem_edgid_num = NULL;
01054 }
01055
01056
01057
01058
01059
01060
01061
01062
01063 sedges::~sedges()
01064 {
01065 long i, j;
01066
01067 delete [] edges;
01068 for (i=0; i<nprop; i++)
01069 {
01070 for (j=0; j<prop_list_elem_num[i]; j++)
01071 delete [] prop_list_elem_edgid[i][j];
01072
01073 delete [] prop_list_nod[i];
01074 delete [] prop_list_elem[i];
01075 delete [] prop_list_elem_edgid[i];
01076 delete [] prop_list_elem_edgid_num[i];
01077 }
01078 if (nprop)
01079 {
01080 delete [] proplist;
01081 delete [] prop_list_nod;
01082 delete [] prop_list_nod_num;
01083 delete [] prop_list_elem;
01084 delete [] prop_list_elem_num;
01085 delete [] prop_list_elem_edgid;
01086 delete [] prop_list_elem_edgid_num;
01087 }
01088 }
01089
01090
01091
01092
01093
01094
01095
01096
01097
01098
01099
01100 sedges::sedges(long n)
01101 {
01102 nedg = n;
01103 nprop = 0L;
01104 edges = new sedge[n];
01105 proplist = NULL;
01106 prop_list_nod = NULL;
01107 prop_list_nod_num = NULL;
01108 prop_list_elem = NULL;
01109 prop_list_elem_num = NULL;
01110 prop_list_elem_edgid = NULL;
01111 prop_list_elem_edgid_num = NULL;
01112 }
01113
01114
01115
01116
01117
01118
01119
01120
01121
01122
01123
01124
01125
01126
01127 void sedges::read(XFILE *in)
01128 {
01129 long i;
01130 for (i=0; i<nedg; i++)
01131 {
01132 if (edges[i].read(in))
01133 {
01134 print_err("cannot read edge number %ld", __FILE__, __LINE__, __func__, i+1);
01135 abort();
01136 }
01137 }
01138 }
01139
01140
01141
01142
01143
01144
01145
01146
01147
01148
01149
01150
01151
01152
01153
01154
01155
01156
01157
01158 void sedges::print(FILE *out, long *lnn)
01159 {
01160 long i;
01161 long *aux, nedg_l;
01162
01163 if (lnn == NULL)
01164 {
01165 fprintf(out, "edges %ld\n", nedg);
01166 for (i=0; i<nedg; i++)
01167 edges[i].print(out, lnn);
01168 }
01169 else
01170 {
01171 nedg_l=0;
01172 aux = new long[nedg];
01173 memset(aux, 0, sizeof(*aux)*nedg);
01174
01175 for(i=0; i<nedg; i++)
01176 {
01177 if (edges[i].cordomtest(lnn))
01178 {
01179 aux[i]++;
01180 nedg_l++;
01181 }
01182 }
01183 fprintf(out, "edges %ld\n", nedg_l);
01184 for(i=0; i<nedg; i++)
01185 {
01186 if (aux[i])
01187 edges[i].print(out, lnn);
01188 }
01189 delete [] aux;
01190 }
01191 }
01192
01193
01194
01195
01196
01197
01198
01199
01200
01201
01202
01203
01204
01205
01206
01207
01208 long sedges::search_list_newprop(long aprop, long *propid, long n)
01209 {
01210 long i;
01211
01212 for(i=0; i<n; i++)
01213 {
01214 if (propid[i] == aprop)
01215 return 0;
01216 }
01217 return 1;
01218 }
01219
01220
01221
01222
01223
01224
01225
01226
01227
01228
01229
01230
01231
01232
01233
01234
01235
01236
01237
01238 void sedges::gen_list_edg(long nn, long ne, selement *elems, long *nadjelnod, long **adjelnod)
01239 {
01240 long i;
01241 long aux;
01242 long *propid;
01243
01244
01245 if (nedg > 0)
01246 {
01247 propid = new long[nedg];
01248 nprop = 0;
01249 propid[nprop] = aux = edges[0].prop;
01250 nprop++;
01251 for(i=0; i<nedg; i++)
01252 {
01253 if (edges[i].prop != aux)
01254 {
01255 aux = edges[i].prop;
01256 if (search_list_newprop(aux, propid, nprop))
01257 {
01258 propid[nprop] = aux;
01259 nprop++;
01260 }
01261 }
01262 }
01263 proplist = new long[nprop];
01264 for(i=0; i<nprop; i++)
01265 proplist[i] = propid[i];
01266 delete [] propid;
01267 prop_list_nod = new long*[nprop];
01268 prop_list_nod_num = new long[nprop];
01269 prop_list_elem = new long*[nprop];
01270 prop_list_elem_num = new long[nprop];
01271 prop_list_elem_edgid = new long**[nprop];
01272 prop_list_elem_edgid_num = new long*[nprop];
01273 memset(prop_list_nod, 0, sizeof(*prop_list_nod)*nprop);
01274 memset(prop_list_nod_num, 0, sizeof(*prop_list_nod_num)*nprop);
01275 memset(prop_list_elem, 0, sizeof(*prop_list_elem)*nprop);
01276 memset(prop_list_elem_num, 0, sizeof(*prop_list_elem_num)*nprop);
01277 memset(prop_list_elem_edgid, 0, sizeof(*prop_list_elem_edgid)*nprop);
01278 memset(prop_list_elem_edgid_num, 0, sizeof(*prop_list_elem_edgid_num)*nprop);
01279
01280 gen_list_node(nn);
01281 gen_list_elem(ne, elems, nadjelnod, adjelnod);
01282 }
01283 }
01284
01285
01286
01287
01288
01289
01290
01291
01292
01293
01294
01295
01296
01297
01298
01299
01300 void sedges::gen_list_node(long nn)
01301 {
01302 long i, j, k, nnodp;
01303 long *nodid;
01304
01305 nodid = new long[nn];
01306
01307 for(i=0; i<nprop; i++)
01308 {
01309 memset(nodid, 0, sizeof(*nodid)*nn);
01310 nnodp = 0;
01311 for(j=0; j<nedg; j++)
01312 {
01313 if (edges[j].prop == proplist[i])
01314 {
01315 if (nodid[edges[j].n1] == 0)
01316 {
01317 nodid[edges[j].n1] = 1;
01318 nnodp++;
01319 }
01320 if (nodid[edges[j].n2] == 0)
01321 {
01322 nodid[edges[j].n2] = 1;
01323 nnodp++;
01324 }
01325 }
01326 }
01327 prop_list_nod_num[i] = nnodp;
01328 prop_list_nod[i] = new long[nnodp];
01329 k = 0;
01330 for(j=0; j<nn; j++)
01331 {
01332 if (nodid[j])
01333 {
01334 prop_list_nod[i][k] = j;
01335 k++;
01336 }
01337 }
01338 }
01339 delete [] nodid;
01340 }
01341
01342
01343
01344
01345
01346
01347
01348
01349
01350
01351
01352
01353
01354
01355
01356
01357
01358
01359
01360
01361
01362
01363
01364
01365
01366
01367
01368
01369
01370
01371
01372
01373
01374
01375
01376
01377
01378
01379
01380
01381
01382
01383
01384
01385
01386
01387
01388
01389
01390
01391
01392
01393
01394
01395
01396
01397
01398
01399
01400
01401
01402
01403
01404
01405
01406
01407
01408
01409
01410
01411
01412
01413
01414
01415
01416
01417
01418
01419
01420
01421
01422
01423
01424
01425
01426
01427
01428
01429
01430
01431
01432
01433
01434
01435
01436
01437
01438
01439
01440
01441
01442
01443
01444
01445
01446
01447 void sedges::gen_list_elem(long ne, selement *elems, long *nadjelnod, long **adjelnod)
01448 {
01449 long i, j, k, l, nelemp;
01450 long aeid1, aeid2;
01451 long *elemid;
01452 long **elemedgid;
01453 long edgid;
01454 long nodes[2];
01455 long aux;
01456
01457 elemid = new long[ne];
01458 elemedgid = new long*[ne];
01459 for(i=0; i<ne; i++)
01460 elemedgid[i] = new long [elems[i].ned];
01461
01462 for(i=0; i<nprop; i++)
01463 {
01464 memset(elemid, 0, sizeof(*elemid)*ne);
01465 for(j=0; j<ne; j++)
01466 memset(elemedgid[j], 0, sizeof(*elemedgid[j])*elems[j].ned);
01467 nelemp = 0;
01468
01469 for(j=0; j<nedg; j++)
01470 {
01471 if (edges[j].prop == proplist[i])
01472 {
01473
01474
01475 for(k=0; k<nadjelnod[edges[j].n1]; k++)
01476 {
01477 aeid1 = adjelnod[edges[j].n1][k];
01478 aux=1;
01479 for(l=0; l<nadjelnod[edges[j].n2]; l++)
01480 {
01481 aeid2 = adjelnod[edges[j].n2][l];
01482 if (aeid2 == aeid1)
01483 {
01484 aux++;
01485 break;
01486 }
01487 if (aeid1 < aeid2)
01488 break;
01489 }
01490 if (aux==2)
01491 {
01492 nodes[0] = edges[j].n1;
01493 nodes[1] = edges[j].n2;
01494 edgid = elems[aeid1].compare_edg(nodes, 2);
01495 if (edgid >= 0)
01496 {
01497 if (elemid[aeid1] == 0)
01498 {
01499 nelemp++;
01500 elemid[aeid1] = 1;
01501 }
01502 elemedgid[aeid1][edgid] = 1;
01503 }
01504 }
01505 }
01506 }
01507 }
01508 prop_list_elem_num[i] = nelemp;
01509 prop_list_elem[i] = new long[nelemp];
01510 prop_list_elem_edgid_num[i] = new long[nelemp];
01511 prop_list_elem_edgid[i] = new long*[nelemp];
01512 k = 0;
01513 for(j=0; j<ne; j++)
01514 {
01515 if (elemid[j] == 0)
01516 continue;
01517
01518 prop_list_elem[i][k] = j;
01519 aux = 0;
01520 for(l = 0; l<elems[j].ned; l++)
01521 {
01522 if (elemedgid[j][l])
01523 aux++;
01524 }
01525 prop_list_elem_edgid_num[i][k] = aux;
01526 if (aux)
01527 prop_list_elem_edgid[i][k] = new long[aux];
01528 else
01529 prop_list_elem_edgid[i][k] = NULL;
01530 aux = 0;
01531 for(l = 0; l<elems[j].ned; l++)
01532 {
01533 if (elemedgid[j][l])
01534 {
01535 prop_list_elem_edgid[i][k][aux] = l;
01536 aux++;
01537 }
01538 }
01539 k++;
01540 }
01541 }
01542 delete [] elemid;
01543 for(i=0; i<ne; i++)
01544 delete [] elemedgid[i];
01545 delete [] elemedgid;
01546 }
01547
01548
01549
01550
01551
01552
01553
01554 sface::sface(void)
01555 {
01556 nnod = 0;
01557 nodes = NULL;
01558 prop = -1;
01559 }
01560
01561
01562
01563
01564
01565
01566
01567
01568
01569
01570
01571 sface::sface(long n)
01572 {
01573 nnod = n;
01574 nodes = new long[n];
01575 memset(nodes, 0, sizeof(*nodes)*n);
01576 prop = -1;
01577 }
01578
01579
01580
01581
01582
01583
01584
01585
01586 sface::~sface(void)
01587 {
01588 delete [] nodes;
01589 }
01590
01591
01592
01593
01594
01595
01596
01597
01598
01599
01600
01601
01602
01603 long sface::read(XFILE *in)
01604 {
01605 long i;
01606
01607 if (xfscanf(in, "%ld", &nnod) == 0)
01608 return 1;
01609
01610 nodes = new long[nnod];
01611 memset(nodes, 0, sizeof(*nodes)*nnod);
01612
01613 for(i=0; i<nnod; i++)
01614 {
01615 if (xfscanf(in, "%ld", nodes+i) == 0)
01616 return 1;
01617
01618 nodes[i]--;
01619
01620 }
01621 if (xfscanf(in, "%ld", &prop) == 0)
01622 return 1;
01623
01624 return 0;
01625 }
01626
01627
01628
01629
01630
01631
01632
01633
01634
01635
01636
01637
01638
01639
01640
01641
01642
01643
01644
01645 void sface::print(FILE *out, long *lnn)
01646 {
01647 long i;
01648
01649 fprintf(out, "%ld", nnod);
01650 for(i=0; i<nnod; i++)
01651 {
01652 if (lnn)
01653 fprintf(out, " %ld", lnn[nodes[i]]);
01654 else
01655 fprintf(out, " %ld", nodes[i]+1);
01656 }
01657
01658 fprintf(out, " %ld\n", prop);
01659 }
01660
01661
01662
01663
01664
01665
01666
01667
01668
01669
01670
01671
01672
01673
01674
01675
01676
01677
01678 long sface::cordomtest(long *lnn)
01679 {
01680 long i, aux=0;
01681
01682 for(i=0; i<nnod; i++)
01683 {
01684 if (lnn[nodes[i]])
01685 aux++;
01686 }
01687 if (aux == nnod)
01688 return 1;
01689
01690 return 0;
01691 }
01692
01693
01694
01695
01696
01697
01698
01699
01700 sfaces::sfaces()
01701 {
01702 nfcs = 0L;
01703 nprop = 0L;
01704 faces = NULL;
01705 proplist = NULL;
01706 prop_list_nod = NULL;
01707 prop_list_nod_num = NULL;
01708 prop_list_elem = NULL;
01709 prop_list_elem_num = NULL;
01710 prop_list_elem_sfid = NULL;
01711 prop_list_elem_sfid_num = NULL;
01712 }
01713
01714
01715
01716
01717
01718
01719
01720
01721 sfaces::~sfaces()
01722 {
01723 long i, j;
01724
01725 delete [] faces;
01726 delete [] proplist;
01727 for (i=0; i<nprop; i++)
01728 {
01729 for (j=0; j<prop_list_elem_num[i]; j++)
01730 delete [] prop_list_elem_sfid[i][j];
01731
01732 delete [] prop_list_nod[i];
01733 delete [] prop_list_elem[i];
01734 delete [] prop_list_elem_sfid[i];
01735 delete [] prop_list_elem_sfid_num[i];
01736 }
01737 if (nprop)
01738 {
01739 delete [] prop_list_nod;
01740 delete [] prop_list_nod_num;
01741 delete [] prop_list_elem;
01742 delete [] prop_list_elem_num;
01743 delete [] prop_list_elem_sfid;
01744 delete [] prop_list_elem_sfid_num;
01745 }
01746 }
01747
01748
01749
01750
01751
01752
01753
01754
01755
01756
01757
01758 sfaces::sfaces(long n)
01759 {
01760 nfcs = n;
01761 nprop = 0L;
01762 faces = new sface[n];
01763 proplist = NULL;
01764 prop_list_nod = NULL;
01765 prop_list_nod_num = NULL;
01766 prop_list_elem = NULL;
01767 prop_list_elem_num = NULL;
01768 prop_list_elem_sfid = NULL;
01769 prop_list_elem_sfid_num = NULL;
01770 }
01771
01772
01773
01774
01775
01776
01777
01778
01779
01780
01781
01782
01783
01784
01785 void sfaces::read(XFILE *in)
01786 {
01787 long i;
01788
01789 for(i=0; i<nfcs; i++)
01790 {
01791 if (faces[i].read(in))
01792 {
01793 print_err("cannot read face number %ld\n", __FILE__, __LINE__, __func__, i+1);
01794 abort();
01795 }
01796 }
01797 }
01798
01799
01800
01801
01802
01803
01804
01805
01806
01807
01808
01809
01810
01811
01812
01813
01814
01815
01816
01817
01818 void sfaces::print(FILE *out, long *lnn)
01819 {
01820 long i;
01821 long *aux;
01822 long nfcs_l;
01823
01824 if (lnn == NULL)
01825 {
01826 fprintf(out, "faces %ld\n", nfcs);
01827 for(i=0; i<nfcs; i++)
01828 faces[i].print(out, lnn);
01829 }
01830 else
01831 {
01832 nfcs_l=0;
01833 aux = new long[nfcs];
01834 memset(aux, 0, sizeof(*aux)*nfcs);
01835
01836 for(i=0; i<nfcs; i++)
01837 {
01838 if (faces[i].cordomtest(lnn))
01839 {
01840 aux[i]++;
01841 nfcs_l++;
01842 }
01843 }
01844 fprintf(out, "faces %ld\n", nfcs_l);
01845 for(i=0; i<nfcs; i++)
01846 {
01847 if (aux[i])
01848 faces[i].print(out, lnn);
01849 }
01850 delete [] aux;
01851 }
01852 }
01853
01854
01855
01856
01857
01858
01859
01860
01861
01862
01863
01864
01865
01866
01867
01868
01869 long sfaces::search_list_newprop(long prop, long *propid, long n)
01870 {
01871 long i;
01872
01873 for(i=0; i<n; i++)
01874 {
01875 if (propid[i] == prop)
01876 return 0;
01877 }
01878 return 1;
01879 }
01880
01881
01882
01883
01884
01885
01886
01887
01888
01889
01890
01891
01892
01893
01894
01895
01896
01897
01898
01899 void sfaces::gen_list_surf(long nn, long ne, selement *elems, long *nadjelnod, long **adjelnod)
01900 {
01901 long i;
01902 long aux;
01903 long *propid;
01904
01905
01906 if (nfcs > 0)
01907 {
01908 propid = new long[nfcs];
01909 nprop = 0;
01910 propid[nprop] = aux = faces[0].prop;
01911 nprop++;
01912 for(i=0; i<nfcs; i++)
01913 {
01914 if (faces[i].prop != aux)
01915 {
01916 aux = faces[i].prop;
01917 if (search_list_newprop(aux, propid, nprop))
01918 {
01919 propid[nprop] = aux;
01920 nprop++;
01921 }
01922 }
01923 }
01924 proplist = new long[nprop];
01925 for(i=0; i<nprop; i++)
01926 proplist[i] = propid[i];
01927 delete [] propid;
01928 prop_list_nod = new long*[nprop];
01929 prop_list_nod_num = new long[nprop];
01930 prop_list_elem = new long*[nprop];
01931 prop_list_elem_num = new long[nprop];
01932 prop_list_elem_sfid = new long**[nprop];
01933 prop_list_elem_sfid_num = new long*[nprop];
01934 memset(prop_list_nod, 0, sizeof(*prop_list_nod)*nprop);
01935 memset(prop_list_nod_num, 0, sizeof(*prop_list_nod_num)*nprop);
01936 memset(prop_list_elem, 0, sizeof(*prop_list_elem)*nprop);
01937 memset(prop_list_elem_num, 0, sizeof(*prop_list_elem_num)*nprop);
01938 memset(prop_list_elem_sfid, 0, sizeof(*prop_list_elem_sfid)*nprop);
01939 memset(prop_list_elem_sfid_num, 0, sizeof(*prop_list_elem_sfid_num)*nprop);
01940
01941 gen_list_node(nn);
01942 gen_list_elem(ne, elems, nadjelnod, adjelnod);
01943 }
01944 }
01945
01946
01947
01948
01949
01950
01951
01952
01953
01954
01955
01956
01957
01958
01959
01960
01961 void sfaces::gen_list_node(long nn)
01962 {
01963 long i, j, k, nnodp;
01964 long *nodid;
01965
01966 nodid = new long[nn];
01967
01968 for(i=0; i<nprop; i++)
01969 {
01970 memset(nodid, 0, sizeof(*nodid)*nn);
01971 nnodp = 0;
01972 for(j=0; j<nfcs; j++)
01973 {
01974 if (faces[j].prop == proplist[i])
01975 {
01976 for(k=0; k<faces[j].nnod; k++)
01977 {
01978 if (nodid[faces[j].nodes[k]] == 0)
01979 {
01980 nodid[faces[j].nodes[k]] = 1;
01981 nnodp++;
01982 }
01983 }
01984 }
01985 }
01986 prop_list_nod_num[i] = nnodp;
01987 prop_list_nod[i] = new long[nnodp];
01988 k = 0;
01989 for(j=0; j<nn; j++)
01990 {
01991 if (nodid[j])
01992 {
01993 prop_list_nod[i][k] = j;
01994 k++;
01995 }
01996 }
01997 }
01998 delete [] nodid;
01999 }
02000
02001
02002
02003
02004
02005
02006
02007
02008
02009
02010
02011
02012
02013
02014
02015
02016
02017
02018
02019
02020
02021
02022
02023
02024
02025 void sfaces::gen_list_elem(long ne, selement *elems, long *nadjelnod, long **adjelnod)
02026 {
02027 long i, j, k, l, m, kmin, nelemp, nadjelmin;
02028 long aeid_min, aeid;
02029 long *elemid;
02030 long **elemsfcid;
02031 long sfid, aux;
02032
02033 elemid = new long[ne];
02034 elemsfcid = new long*[ne];
02035 for(i=0; i<ne; i++)
02036 elemsfcid[i] = new long [elems[i].nsurf];
02037
02038 for(i=0; i<nprop; i++)
02039 {
02040 memset(elemid, 0, sizeof(*elemid)*ne);
02041 for(j=0; j<ne; j++)
02042 memset(elemsfcid[j], 0, sizeof(*elemsfcid[j])*elems[j].nsurf);
02043 nelemp = 0;
02044
02045 for(j=0; j<nfcs; j++)
02046 {
02047 if (faces[j].prop == proplist[i])
02048 {
02049
02050 nadjelmin = ne;
02051 kmin = 0;
02052 for (k=0; k<faces[j].nnod; k++)
02053 {
02054 if (nadjelnod[faces[j].nodes[k]] < nadjelmin)
02055 {
02056 kmin = k;
02057 nadjelmin = nadjelnod[faces[j].nodes[k]];
02058 }
02059 }
02060
02061
02062
02063 for(k=0; k<nadjelmin; k++)
02064 {
02065 aeid_min = adjelnod[faces[j].nodes[kmin]][k];
02066 aux=1;
02067 for(l=0; l<faces[j].nnod; l++)
02068 {
02069 if(kmin == l)
02070 continue;
02071
02072 for(m=0; m<nadjelnod[faces[j].nodes[l]]; m++)
02073 {
02074 aeid = adjelnod[faces[j].nodes[l]][m];
02075 if (aeid == aeid_min)
02076 {
02077 aux++;
02078 break;
02079 }
02080 if (aeid_min < aeid)
02081 {
02082 aux = -1;
02083 break;
02084 }
02085 }
02086 if (aux < 0)
02087 break;
02088 }
02089 if (aux==faces[j].nnod)
02090 {
02091 sfid = elems[aeid_min].compare_surf(faces[j].nodes, faces[j].nnod);
02092 if (sfid >= 0)
02093 {
02094 if (elemid[aeid_min] == 0)
02095 {
02096 nelemp++;
02097 elemid[aeid_min] = 1;
02098 }
02099 elemsfcid[aeid_min][sfid] = 1;
02100 }
02101 }
02102 }
02103 }
02104 }
02105 prop_list_elem_num[i] = nelemp;
02106 prop_list_elem[i] = new long[nelemp];
02107 prop_list_elem_sfid_num[i] = new long[nelemp];
02108 prop_list_elem_sfid[i] = new long*[nelemp];
02109 k = 0;
02110 for(j=0; j<ne; j++)
02111 {
02112 if (elemid[j] == 0)
02113 continue;
02114
02115 prop_list_elem[i][k] = j;
02116 aux = 0;
02117 for(l = 0; l<elems[j].nsurf; l++)
02118 {
02119 if (elemsfcid[j][l])
02120 aux++;
02121 }
02122 prop_list_elem_sfid_num[i][k] = aux;
02123 if (aux)
02124 prop_list_elem_sfid[i][k] = new long[aux];
02125 else
02126 prop_list_elem_sfid[i][k] = NULL;
02127 aux = 0;
02128 for(l = 0; l<elems[j].nsurf; l++)
02129 {
02130 if (elemsfcid[j][l])
02131 {
02132 prop_list_elem_sfid[i][k][aux] = l;
02133 aux++;
02134 }
02135 }
02136 k++;
02137 }
02138 }
02139 delete [] elemid;
02140 for(i=0; i<ne; i++)
02141 delete [] elemsfcid[i];
02142 delete [] elemsfcid;
02143 }
02144
02145
02146
02147
02148
02149
02150
02151 siftop :: siftop (void)
02152 {
02153 elements = NULL;
02154 nodes = NULL;
02155 gnn = NULL;
02156 edges = NULL;
02157 surfaces = NULL;
02158 nadjelnod = NULL;
02159 adjelnod = NULL;
02160 nn = 0;
02161 ne = 0;
02162 memset(npet, 0, sizeof(npet));
02163 }
02164
02165
02166
02167
02168
02169
02170
02171
02172 siftop :: ~siftop(void)
02173 {
02174 long i;
02175
02176 if (elements)
02177 delete [] elements;
02178
02179 if (nodes)
02180 delete [] nodes;
02181
02182 if (gnn)
02183 delete [] gnn;
02184
02185 if (edges)
02186 delete edges;
02187
02188 if (surfaces)
02189 delete surfaces;
02190
02191 if (nadjelnod)
02192 {
02193 for (i=0; i<nn; i++)
02194 delete [] adjelnod[i];
02195 delete [] nadjelnod;
02196 delete [] adjelnod;
02197 }
02198 }
02199
02200
02201
02202
02203
02204
02205
02206
02207
02208
02209
02210
02211
02212
02213
02214
02215
02216
02217
02218 long siftop::read(XFILE *in, long rgnn, long rebp)
02219 {
02220 long i, j, eid, nid;
02221 long nfaces, nedges;
02222 char emsg[1001];
02223
02224
02225 xfscanf(in, "%k%ld", "num_nodes", &nn);
02226 nodes = new snode[nn];
02227 for (i = 0; i < nn; i++)
02228 {
02229 xfscanf(in, "%k%ld", "node_id", &nid);
02230 if ((nid < 1) && (nid > nn))
02231 {
02232 sprintf(emsg, "node number %ld is out of range <1, %ld> (file %s, line=%ld, col=%ld)\n", nid, nn, in->fname, in->line, in->col);
02233 print_err(emsg, __FILE__, __LINE__, __func__);
02234 return (1);
02235 }
02236 nid--;
02237 xfscanf(in, "%k%le %k%le %k%le %k%ld", "x", &nodes[nid].x, "y", &nodes[nid].y, "z", &nodes[nid].z, "numprop", &nodes[nid].nprop);
02238 if (nodes[nid].nprop < 1)
02239 {
02240 sprintf(emsg, "node %ld must have at least one property (file %s, line=%ld, col=%ld)\n", nid, in->fname, in->line, in->col);
02241 print_err(emsg, __FILE__, __LINE__, __func__);
02242 return (1);
02243 }
02244 nodes[nid].alloc(nodes[nid].nprop);
02245 for (j=0; j < nodes[nid].nprop; j++)
02246 xfscanf(in, "%k%m%ld", "prop", &entityp_kwdset, nodes[nid].entid+j, nodes[nid].prop+j);
02247 }
02248
02249
02250
02251 xfscanf(in, "%k%ld", "num_elements", &ne);
02252 elements = new selement[ne];
02253 for (i = 0; i < ne; i++)
02254 {
02255 xfscanf(in, "%k%ld", "elem_id", &eid);
02256 eid--;
02257 if ((eid < 0) && (eid >= ne))
02258 {
02259 sprintf(emsg, "element number %ld is out of range <1, %ld> (file %s, line=%ld, col=%ld)", eid+1, nn, in->fname, in->line, in->col);
02260 print_err(emsg, __FILE__, __LINE__, __func__);
02261 return (2);
02262 }
02263 xfscanf(in, "%k%m", "eltype", >ypel_kwdset, &elements[eid].type);
02264 if (elements[eid].alloc(rebp))
02265 {
02266 sprintf(emsg, "unknown type on element %ld is required (file %s, line=%ld, col=%ld)", eid+1, in->fname, in->line, in->col);
02267 print_err(emsg, __FILE__, __LINE__, __func__);
02268 }
02269 xfscanf(in, "%k", "enodes");
02270 for (j = 0; j < elements[eid].nne; j++)
02271 {
02272 xfscanf(in, "%ld", &elements[eid].nodes[j]);
02273 elements[eid].nodes[j]--;
02274 }
02275 xfscanf(in, "%k%ld", "eprop", &elements[eid].prop);
02276 if (rebp)
02277 {
02278 if (elements[eid].ned)
02279 xfscanf(in, "%k", "propedg");
02280 for(j = 0; j < elements[eid].ned; j++)
02281 xfscanf(in, "%ld", elements[eid].propedg+j);
02282 if (elements[eid].nsurf)
02283 xfscanf(in, "%k", "propsurf");
02284 for(j = 0; j < elements[eid].nsurf; j++)
02285 xfscanf(in, "%ld", elements[eid].propsurf+j);
02286 }
02287 else{
02288 xfscanf(in, "%*[^\n]");
02289 }
02290 }
02291
02292
02293 if (rgnn==1){
02294
02295 gnn = new long[nn];
02296 memset(gnn, 0, sizeof(*gnn)*nn);
02297 for (i = 0; i < nn; i++){
02298 xfscanf(in, "%k%ld", "node_id", &j);
02299 if ((j < 1) && (j > nn)){
02300 sprintf(emsg, "node number %ld is out of range <1, %ld> (file %s, line=%ld, col=%ld)\n", j, nn, in->fname, in->line, in->col);
02301 print_err(emsg, __FILE__, __LINE__, __func__);
02302 return (3);
02303 }
02304 j--;
02305 xfscanf(in, "%k%ld", "glob_id", gnn+j);
02306 }
02307 }
02308
02309 if (rgnn==2){
02310
02311
02312
02313
02314
02315
02316 xfscanf (in,"%ld %ld",&meshtype,&nsd);
02317
02318
02319 nnsd = new long [nsd];
02320
02321 for (i=0;i<nsd;i++){
02322 xfscanf (in,"%ld",nnsd+i);
02323 }
02324
02325
02326 gnn = new long[nn];
02327 memset(gnn, 0, sizeof(*gnn)*nn);
02328 for (i = 0; i < nn; i++){
02329 xfscanf(in, "%k%ld", "node_id", &j);
02330 if ((j < 1) && (j > nn)){
02331 sprintf(emsg, "node number %ld is out of range <1, %ld> (file %s, line=%ld, col=%ld)\n", j, nn, in->fname, in->line, in->col);
02332 print_err(emsg, __FILE__, __LINE__, __func__);
02333 return (3);
02334 }
02335 j--;
02336 xfscanf(in, "%k%ld", "glob_id", gnn+j);
02337 }
02338 }
02339
02340 kwd_handling bmode = in->kwdmode;
02341 in->kwdmode = sequent_mode;
02342
02343 if (xfscanf(in, " %+k", "faces"))
02344 {
02345 in->kwdmode = bmode;
02346 xfscanf(in, "%ld", &nfaces);
02347 surfaces = new sfaces(nfaces);
02348 surfaces->read(in);
02349 if (nadjelnod == NULL)
02350 gen_adjelnod();
02351 surfaces->gen_list_surf(nn, ne, elements, nadjelnod, adjelnod);
02352 }
02353
02354 in->kwdmode = sequent_mode;
02355
02356 if (xfscanf(in, " %+k", "edges"))
02357 {
02358 in->kwdmode = bmode;
02359 xfscanf(in, "%ld", &nedges);
02360 edges = new sedges(nedges);
02361 edges->read(in);
02362 if (nadjelnod == NULL)
02363 gen_adjelnod();
02364 edges->gen_list_edg(nn, ne, elements, nadjelnod, adjelnod);
02365 }
02366 in->kwdmode = bmode;
02367
02368 return(0);
02369 }
02370
02371
02372
02373
02374
02375
02376
02377
02378
02379
02380
02381 void siftop::exporttop (gtopology *top)
02382 {
02383 long i, j;
02384 top->alloc_nodes(nn);
02385 top->alloc_elements(ne);
02386 for (i = 0; i < nn; i++)
02387 {
02388 top->gnodes[i].x = nodes[i].x;
02389 top->gnodes[i].y = nodes[i].y;
02390 top->gnodes[i].z = nodes[i].z;
02391 }
02392 for (i = 0; i < ne; i++)
02393 {
02394 top->gelements[i].nne = elements[i].nne;
02395 top->gelements[i].nodes = new long[elements[i].nne];
02396 for (j = 0; j < elements[i].nne; j++)
02397 top->gelements[i].nodes[j] = elements[i].nodes[j];
02398 switch(elements[i].type){
02399 case isolinear1d:{
02400 top->gelements[i].get=linbar;
02401 break;
02402 }
02403 case isoquadratic1d:{
02404 top->gelements[i].get=quadbar;
02405 break;
02406 }
02407 case trianglelinear:{
02408 top->gelements[i].get=lintriag;
02409 break;
02410 }
02411 case trianglequadratic:{
02412 top->gelements[i].get=quadtriag;
02413 break;
02414 }
02415 case isolinear2d:{
02416 top->gelements[i].get=linquad;
02417 break;
02418 }
02419 case isoquadratic2d:{
02420 top->gelements[i].get=quadquad;
02421 break;
02422 }
02423 case tetrahedronlinear:{
02424 top->gelements[i].get=lintetra;
02425 break;
02426 }
02427 case tetrahedronquadratic:{
02428 top->gelements[i].get=quadtetra;
02429 break;
02430 }
02431 case pyramidelinear:{
02432 top->gelements[i].get=noelem;
02433 break;
02434 }
02435 case pyramidequadratic:{
02436 top->gelements[i].get=noelem;
02437 break;
02438 }
02439 case vedgelinear:{
02440 top->gelements[i].get=noelem;
02441 break;
02442 }
02443 case vedgequadratic:{
02444 top->gelements[i].get=noelem;
02445 break;
02446 }
02447 case isolinear3d:{
02448 top->gelements[i].get=linhexa;
02449 break;
02450 }
02451 case isoquadratic3d:{
02452 top->gelements[i].get=quadhexa;
02453 break;
02454 }
02455 case all_element_types:{
02456 top->gelements[i].get=noelem;
02457 break;
02458 }
02459 default:{
02460 print_err("unknown type of element is required", __FILE__, __LINE__, __func__);
02461 }
02462 }
02463 }
02464
02465 }
02466
02467
02468
02469
02470
02471
02472
02473
02474
02475
02476
02477
02478
02479
02480 void siftop::import_t3d (XFILE *in, long paral)
02481 {
02482 long i;
02483 long type, renum, mode, temp;
02484 long outmode[9];
02485
02486
02487 xfscanf(in, "%ld", &type);
02488 xfscanf(in, "%ld", °ree);
02489 xfscanf(in, "%ld", &renum);
02490 xfscanf(in, "%ld", &mode);
02491
02492 if (paral==1)
02493 xfscanf(in,"%ld", &temp);
02494
02495
02496 for (i = 0; i < 9; i++)
02497 outmode[i]=(mode>>i) & 1;
02498
02499 xfscanf(in, "%ld", &nn);
02500
02501 if (nodes)
02502 delete [] nodes;
02503 nodes = new snode[nn];
02504 if (paral)
02505 {
02506 if (gnn)
02507 delete [] gnn;
02508 gnn = new long[nn];
02509 memset(gnn, 0, sizeof(*gnn)*nn);
02510 }
02511 switch (type)
02512
02513 {
02514 case 3:
02515
02516 xfscanf(in, "%ld%ld%ld", npet+isolinear1d+degree-2, npet+trianglelinear+degree-2, npet+tetrahedronlinear+degree-2);
02517 break ;
02518 case 4:
02519
02520 xfscanf(in, "%ld%ld%ld", npet+isolinear1d+degree-2, npet+isolinear2d+degree-2, npet+isolinear3d+degree-2);
02521 break ;
02522 case 7:
02523
02524 xfscanf(in, "%ld%ld%ld", npet+isolinear1d+degree-2, npet+trianglelinear+degree-2, npet+isolinear2d+degree-2);
02525
02526 xfscanf(in, "%ld%ld%ld%ld", npet+tetrahedronlinear+degree-2, npet+pyramidelinear+degree-2, npet+vedgelinear+degree-2, npet+isolinear3d+degree-2);
02527 break ;
02528 default:
02529 print_err("unknown type of generator is required", __FILE__, __LINE__, __func__);
02530 }
02531
02532 ne=0;
02533 for (i = 0; i < all_element_types-1; i++)
02534 ne+=npet[i];
02535
02536
02537 import_t3d_nodes(in, paral, outmode) ;
02538
02539
02540 if (elements)
02541 delete [] elements;
02542 elements = new selement[ne];
02543
02544
02545 import_t3d_elements(in, outmode);
02546
02547 return;
02548 }
02549
02550
02551
02552
02553
02554
02555
02556
02557
02558
02559
02560
02561
02562
02563
02564
02565 long siftop::import_t3d_nodes (XFILE *in, long paral, long *outmode)
02566 {
02567 long i, j, nid;
02568 long nent;
02569 entityp ent;
02570 long eprop;
02571 char emsg[1001];
02572
02573 for (i = 0; i < nn; i++)
02574 {
02575 xfscanf(in, "%ld", &nid);
02576 if ((nid < 1) || (nid > nn))
02577 {
02578 sprintf(emsg, "node number %ld is out of range <1, %ld> (file %s, line=%ld, col=%ld)\n", nid, nn, in->fname, in->line, in->col);
02579 print_err(emsg, __FILE__, __LINE__, __func__);
02580 return (1);
02581 }
02582 nid--;
02583 if (paral)
02584 xfscanf(in, "%ld", gnn+nid);
02585 xfscanf(in, "%le %le %le %m %*ld %ld", &nodes[nid].x, &nodes[nid].y, &nodes[nid].z, &entityp_kwdset, &ent, &eprop);
02586 if (outmode[omt_sec_clasif_nod])
02587 {
02588 xfscanf(in, "%ld", &nent);
02589 nent++;
02590 nodes[nid].alloc(nent);
02591 nodes[nid].entid[0] = ent;
02592 nodes[nid].prop[0] = eprop;
02593 for (j=1; j < nent; j++)
02594 xfscanf(in, "%m %*ld %ld", &entityp_kwdset, nodes[nid].entid+j, nodes[nid].prop+j);
02595 }
02596 else
02597 {
02598 nent = 1;
02599 nodes[nid].alloc(nent);
02600 nodes[nid].entid[0] = ent;
02601 nodes[nid].prop[0] = eprop;
02602 }
02603 switch (nodes[nid].entid[0])
02604 {
02605 case evertex:
02606 case eregion:
02607 break ;
02608 case ecurve:
02609 if (outmode[omt_parameters] || outmode[omt_tangent])
02610 skipline(in);
02611 break ;
02612 case esurface:
02613 case epatch:
02614 case eshell:
02615 if (outmode[omt_parameters] || outmode[omt_normal])
02616 skipline(in);
02617 break ;
02618 }
02619 }
02620 return(0);
02621 }
02622
02623
02624
02625
02626
02627
02628
02629
02630
02631
02632
02633
02634
02635
02636 long siftop::import_t3d_elements (XFILE *in, long *outmode)
02637 {
02638 long i, j, eid;
02639 gtypel etype;
02640 char emsg[1001];
02641
02642 if (degree == 1)
02643 etype=isolinear1d;
02644 else
02645 etype=isoquadratic1d;
02646
02647 for ( ; etype<all_element_types; etype++, etype++)
02648 {
02649 for (i=0; i<npet[etype-1]; i++)
02650 {
02651 xfscanf(in, "%ld", &eid);
02652 eid--;
02653 if ((eid < 0) || (eid >= ne))
02654 {
02655 sprintf(emsg, "element number %ld is out of range <1, %ld> (file %s, line=%ld, col=%ld)", eid+1, nn, in->fname, in->line, in->col);
02656 print_err(emsg, __FILE__, __LINE__, __func__);
02657 return(1);
02658 }
02659
02660 elements[eid].type = etype;
02661 if (outmode[omt_neighbourhood])
02662 elements[eid].alloc(1);
02663 else
02664 elements[eid].alloc(1);
02665
02666 for(j=0; j<elements[eid].nne; j++)
02667 {
02668 xfscanf(in, "%ld", &elements[eid].nodes[j]);
02669 elements[eid].nodes[j]--;
02670 }
02671
02672 xfscanf(in, "%*ld %*ld %ld", &elements[eid].prop);
02673 switch (etype)
02674 {
02675 case isolinear1d:
02676 case isoquadratic1d:
02677 break;
02678 case trianglelinear:
02679 case trianglequadratic:
02680 case isolinear2d:
02681 case isoquadratic2d:
02682 {
02683 if (outmode[omt_neighbourhood])
02684 skipline(in);
02685 if (outmode[omt_boundary_entities])
02686 {
02687 for (j=0; j<elements[eid].ned; j++)
02688 xfscanf(in, "%*ld");
02689 for (j=0; j<elements[eid].ned; j++)
02690 xfscanf(in, "%ld", elements[eid].propedg+j);
02691 skipline(in);
02692 }
02693 else
02694 skipline(in);
02695 break;
02696 }
02697 case tetrahedronlinear:
02698 case tetrahedronquadratic:
02699 case pyramidelinear:
02700 case pyramidequadratic:
02701 case vedgelinear:
02702 case vedgequadratic:
02703 case isolinear3d:
02704 case isoquadratic3d:
02705 {
02706 if (outmode[omt_neighbourhood] || outmode[omt_associated_elements])
02707 skipline(in);
02708 if (outmode[omt_boundary_entities])
02709 {
02710 for (j=0; j<elements[eid].nsurf; j++)
02711 xfscanf(in, "%*ld");
02712 for (j=0; j<elements[eid].nsurf; j++)
02713 xfscanf(in, "%*ld");
02714 for (j=0; j<elements[eid].nsurf; j++)
02715 xfscanf(in, "%ld", elements[eid].propsurf+j);
02716 skipline(in);
02717 }
02718 else
02719 skipline(in);
02720 break;
02721 }
02722 default:
02723 sprintf(emsg, "unknown type on element %ld is required (file %s, line=%ld, col=%ld)", eid+1, in->fname, in->line, in->col);
02724 print_err(emsg, __FILE__, __LINE__, __func__);
02725 break;
02726 }
02727 }
02728 }
02729 return(0);
02730 }
02731
02732
02733
02734 enum bsecgid {begsecgid_coord=0, begsecgid_elem=1};
02735 const enumstr bsecgid_str[] = {{"Coordinates",0}, {"Elements",1}};
02736 const kwdset bsecgid_kwdset(sizeof(bsecgid_str)/sizeof(*bsecgid_str), bsecgid_str);
02737
02738
02739 enum esecgid {endsecgid_coord=0, endsecgid_elem=1};
02740 const enumstr esecgid_str[] = {{"end Coordinates",0}, {"end Elements",1}};
02741 const kwdset esecgid_kwdset(sizeof(esecgid_str)/sizeof(*esecgid_str), esecgid_str);
02742
02743
02744
02745
02746
02747
02748
02749
02750
02751
02752
02753
02754 void siftop::import_gid (XFILE *in)
02755 {
02756
02757
02758
02759
02760
02761
02762
02763
02764
02765
02766
02767
02768
02769
02770
02771
02772
02773
02774
02775
02776
02777
02778
02779
02780
02781
02782
02783
02784
02785
02786
02787
02788
02789 return;
02790 }
02791
02792
02793
02794
02795
02796
02797
02798
02799 void siftop::print (FILE *out)
02800 {
02801 long i, j;
02802
02803 fprintf (out,"%ld\n",nn);
02804
02805 for (i = 0; i < nn; i++)
02806 {
02807 fprintf (out,"%ld %15.10le %15.10le %15.10le %ld",i+1, nodes[i].x, nodes[i].y, nodes[i].z, nodes[i].nprop);
02808 for(j=0; j<nodes[i].nprop; j++)
02809 {
02810 fprintf(out, " %d %ld", nodes[i].entid[j], nodes[i].prop[j]);
02811 }
02812 fprintf(out, "\n");
02813 }
02814
02815 fprintf (out,"\n%ld\n", ne);
02816
02817 for (i = 0; i < ne; i++)
02818 {
02819 fprintf (out,"%ld %d",i+1, elements[i].type);
02820 for (j = 0; j < elements[i].nne; j++)
02821 fprintf (out," %ld",elements[i].nodes[j]+1);
02822 fprintf (out," %ld",elements[i].prop);
02823 if (elements[i].propedg)
02824 {
02825 fprintf(out, " ");
02826 for (j=0; j < elements[i].ned; j++)
02827 fprintf (out," %ld",elements[i].propedg[j]);
02828 }
02829 if (elements[i].propsurf)
02830 {
02831 fprintf(out, " ");
02832 for (j=0; j < elements[i].nsurf; j++)
02833 fprintf (out," %ld",elements[i].propsurf[j]);
02834 }
02835 fprintf(out, "\n");
02836 }
02837 if (edges)
02838 edges->print(out, NULL);
02839 if (surfaces)
02840 surfaces->print(out, NULL);
02841 }
02842
02843
02844
02845
02846
02847
02848
02849
02850
02851
02852
02853
02854
02855 void siftop::shift_print_nodes (FILE *out,long shift)
02856 {
02857 long i,j;
02858 for (i = 0; i < nn; i++){
02859 fprintf (out,"%ld %15.10le %15.10le %15.10le %ld",i+1+shift, nodes[i].x, nodes[i].y, nodes[i].z, nodes[i].nprop);
02860 for(j=0; j<nodes[i].nprop; j++)
02861 {
02862 fprintf(out, " %d %ld", nodes[i].entid[j], nodes[i].prop[j]);
02863 }
02864 fprintf(out, "\n");
02865 }
02866 }
02867
02868
02869
02870
02871
02872
02873
02874
02875
02876
02877
02878
02879 void siftop::shift_print_elements (FILE *out,long shiftnn,long shiftne){
02880
02881 long i,j;
02882 for (i = 0; i < ne; i++){
02883 fprintf (out,"%ld %d",i+1+shiftne, elements[i].type);
02884 for (j = 0; j < elements[i].nne; j++)
02885 fprintf (out," %ld",elements[i].nodes[j]+1+shiftnn);
02886 fprintf (out," %ld",elements[i].prop);
02887 if (elements[i].propedg)
02888 {
02889 fprintf(out, " ");
02890 for (j=0; j < elements[i].ned; j++)
02891 fprintf (out," %ld",elements[i].propedg[j]);
02892 }
02893 if (elements[i].propsurf)
02894 {
02895 fprintf(out, " ");
02896 for (j=0; j < elements[i].nsurf; j++)
02897 fprintf (out," %ld",elements[i].propsurf[j]);
02898 }
02899 fprintf(out, "\n");
02900 }
02901 }
02902
02903
02904
02905
02906
02907
02908
02909
02910
02911
02912 void siftop::gen_adjelnod(void)
02913 {
02914 long i, j, nne, *enodes;
02915
02916
02917 if (nadjelnod==NULL)
02918 nadjelnod = new long [nn];
02919
02920 memset (nadjelnod,0,nn*sizeof(*nadjelnod));
02921
02922
02923 for (i=0; i<ne; i++)
02924 {
02925 nne = elements[i].nne;
02926 enodes = elements[i].nodes;
02927
02928 for (j=0;j<nne;j++)
02929 nadjelnod[enodes[j]]++;
02930 }
02931
02932
02933 if (adjelnod==NULL)
02934 {
02935 adjelnod = new long*[nn];
02936
02937 for (i=0;i<nn;i++)
02938 adjelnod[i] = new long[nadjelnod[i]];
02939 }
02940
02941 memset (nadjelnod,0,nn*sizeof(*nadjelnod));
02942
02943
02944
02945
02946 for (i=0;i<ne;i++)
02947 {
02948 nne = elements[i].nne;
02949 enodes = elements[i].nodes;
02950
02951 for (j=0;j<nne;j++)
02952 adjelnod[enodes[j]][nadjelnod[enodes[j]]++] = i;
02953 }
02954
02955 return;
02956 }
02957
02958
02959
02960
02961
02962
02963
02964
02965
02966
02967
02968
02969
02970
02971
02972
02973
02974 long siftop::get_propent_nodes(long prop, gentity ent, long *setnodes)
02975 {
02976 long i, j, ret=0;
02977 for (i=0; i<nn; i++)
02978 {
02979 setnodes[i] = -1;
02980 if (nodes[i].searchprop(prop, ent, i, NULL, NULL))
02981 {
02982 setnodes[i] = i;
02983 ret++;
02984 }
02985 }
02986 if ((ent == gcurve) && (edges))
02987 {
02988
02989 for(i=0; i<edges->nprop; i++)
02990 {
02991 if (edges->proplist[i] == prop)
02992 {
02993 for(j=0; j<edges->prop_list_nod_num[i]; j++)
02994 {
02995 setnodes[edges->prop_list_nod[i][j]] = edges->prop_list_nod[i][j];
02996 ret++;
02997 }
02998 }
02999 }
03000 }
03001
03002 if ((ent == gsurface) && (surfaces))
03003 {
03004
03005 for(i=0; i<surfaces->nprop; i++)
03006 {
03007 if (surfaces->proplist[i] == prop)
03008 {
03009 for(j=0; j<surfaces->prop_list_nod_num[i]; j++)
03010 {
03011 setnodes[surfaces->prop_list_nod[i][j]] = surfaces->prop_list_nod[i][j];
03012 ret++;
03013 }
03014 }
03015 }
03016 }
03017 return ret;
03018 }
03019
03020
03021
03022
03023
03024
03025
03026
03027
03028
03029
03030
03031
03032
03033
03034
03035
03036
03037
03038
03039 long siftop::get_ndofn(long prop, gentity ent, long *ndofn, long *setnodes)
03040 {
03041 long i, j, id;
03042 long ndof = -1;
03043 for (i=0; i<nn; i++)
03044 {
03045 setnodes[i] = -1;
03046 if (nodes[i].searchprop(prop, ent, i, NULL, NULL))
03047 {
03048 setnodes[i] = i;
03049 if (ndof == -1)
03050 ndof = ndofn[i];
03051 else
03052 {
03053 if (ndof != ndofn[i])
03054 return 0;
03055 }
03056 }
03057 }
03058
03059
03060 if ((ent == gcurve) && (edges))
03061 {
03062
03063 for(i=0; i<edges->nprop; i++)
03064 {
03065 if (edges->proplist[i] == prop)
03066 {
03067 for(j=0; j<edges->prop_list_nod_num[i]; j++)
03068 {
03069 id = edges->prop_list_nod[i][j];
03070 setnodes[id] = id;
03071 if (ndof == -1)
03072 ndof = ndofn[id];
03073 else
03074 {
03075 if (ndof != ndofn[id])
03076 return 0;
03077 }
03078 }
03079 }
03080 }
03081 }
03082
03083 if ((ent == gsurface) && (surfaces))
03084 {
03085
03086 for(i=0; i<surfaces->nprop; i++)
03087 {
03088 if (surfaces->proplist[i] == prop)
03089 {
03090 for(j=0; j<surfaces->prop_list_nod_num[i]; j++)
03091 {
03092 id = surfaces->prop_list_nod[i][j];
03093 setnodes[id] = id;
03094 if (ndof == -1)
03095 ndof = ndofn[id];
03096 else
03097 {
03098 if (ndof != ndofn[id])
03099 return 0;
03100 }
03101 }
03102 }
03103 }
03104 }
03105 return ndof;
03106 }
03107
03108
03109
03110
03111
03112
03113
03114
03115
03116
03117
03118
03119
03120
03121
03122 void siftop::export_gid_mesh(FILE *out, long idn1, long ide1)
03123 {
03124 long i, print_header, print_coord = 1;
03125
03126 print_header = 1;
03127 for (i=0; i < ne; i++)
03128 {
03129 if (elements[i].type == isolinear1d)
03130 {
03131 if (print_header)
03132 {
03133 fprintf(out, "MESH \"%ld-beams2\" dimension 3 Elemtype Linear Nnode 2\n", ide1);
03134 print_header = 0;
03135 if (print_coord)
03136 {
03137 write_gid_nodes(out, idn1);
03138 print_coord = 0;
03139 }
03140 fprintf(out, "Elements\n");
03141 }
03142 write_gid_element(out, i, idn1, ide1);
03143 }
03144 }
03145 if (print_header == 0)
03146 fprintf(out, "end Elements\n");
03147
03148 print_header = 1;
03149 for (i=0; i < ne; i++)
03150 {
03151 if (elements[i].type == isoquadratic1d)
03152 {
03153 if (print_header)
03154 {
03155 fprintf(out, "MESH \"%ld-beams3\" dimension 3 Elemtype Linear Nnode 3\n", ide1);
03156 print_header = 0;
03157 if (print_coord)
03158 {
03159 write_gid_nodes(out, idn1);
03160 print_coord = 0;
03161 }
03162 fprintf(out, "Elements\n");
03163 }
03164 write_gid_element(out, i, idn1, ide1);
03165 }
03166 }
03167 if (print_header == 0)
03168 fprintf(out, "end Elements\n");
03169
03170 print_header = 1;
03171 for (i=0; i < ne; i++)
03172 {
03173 if (elements[i].type == trianglelinear)
03174 {
03175 if (print_header)
03176 {
03177 fprintf(out, "MESH \"%ld-trias3\" dimension 3 Elemtype Triangle Nnode 3\n", ide1);
03178 print_header = 0;
03179 if (print_coord)
03180 {
03181 write_gid_nodes(out, idn1);
03182 print_coord = 0;
03183 }
03184 fprintf(out, "Elements\n");
03185 }
03186 write_gid_element(out, i, idn1, ide1);
03187 }
03188 }
03189 if (print_header == 0)
03190 fprintf(out, "end Elements\n");
03191
03192 print_header = 1;
03193 for (i=0; i < ne; i++)
03194 {
03195 if (elements[i].type == trianglequadratic)
03196 {
03197 if (print_header)
03198 {
03199 fprintf(out, "MESH \"%ld-trias6\" dimension 3 Elemtype Triangle Nnode 6\n", ide1);
03200 print_header = 0;
03201 if (print_coord)
03202 {
03203 write_gid_nodes(out, idn1);
03204 print_coord = 0;
03205 }
03206 fprintf(out, "Elements\n");
03207 }
03208 write_gid_element(out, i, idn1, ide1);
03209 }
03210 }
03211 if (print_header == 0)
03212 fprintf(out, "end Elements\n");
03213
03214 print_header = 1;
03215 for (i=0; i < ne; i++)
03216 {
03217 if (elements[i].type == isolinear2d)
03218 {
03219 if (print_header)
03220 {
03221 fprintf(out, "MESH \"%ld-Quads4\" dimension 3 Elemtype Quadrilateral Nnode 4\n", ide1);
03222 print_header = 0;
03223 if (print_coord)
03224 {
03225 write_gid_nodes(out, idn1);
03226 print_coord = 0;
03227 }
03228 fprintf(out, "Elements\n");
03229 }
03230 write_gid_element(out, i, idn1, ide1);
03231 }
03232 }
03233 if (print_header == 0)
03234 fprintf(out, "end Elements\n");
03235
03236 print_header = 1;
03237 for (i=0; i < ne; i++)
03238 {
03239 if (elements[i].type == isoquadratic2d)
03240 {
03241 if (print_header)
03242 {
03243 fprintf(out, "MESH \"%ld-Quads8\" dimension 3 Elemtype Quadrilateral Nnode 8\n", ide1);
03244 print_header = 0;
03245 if (print_coord)
03246 {
03247 write_gid_nodes(out, idn1);
03248 print_coord = 0;
03249 }
03250 fprintf(out, "Elements\n");
03251 }
03252 write_gid_element(out, i, idn1, ide1);
03253 }
03254 }
03255 if (print_header == 0)
03256 fprintf(out, "end Elements\n");
03257
03258 print_header = 1;
03259 for (i=0; i < ne; i++)
03260 {
03261 if (elements[i].type == tetrahedronlinear)
03262 {
03263 if (print_header)
03264 {
03265 fprintf(out, "MESH \"%ld-Tetras4\" dimension 3 Elemtype Tetrahedra Nnode 4\n", ide1);
03266 print_header = 0;
03267 if (print_coord)
03268 {
03269 write_gid_nodes(out, idn1);
03270 print_coord = 0;
03271 }
03272 fprintf(out, "Elements\n");
03273 }
03274 write_gid_element(out, i, idn1, ide1);
03275 }
03276 }
03277 if (print_header == 0)
03278 fprintf(out, "end Elements\n");
03279
03280 print_header = 1;
03281 for (i=0; i < ne; i++)
03282 {
03283 if (elements[i].type == tetrahedronquadratic)
03284 {
03285 if (print_header)
03286 {
03287 fprintf(out, "MESH \"%ld-Tetras10\" dimension 3 Elemtype Tetrahedra Nnode 10\n", ide1);
03288 print_header = 0;
03289 if (print_coord)
03290 {
03291 write_gid_nodes(out, idn1);
03292 print_coord = 0;
03293 }
03294 fprintf(out, "Elements\n");
03295 }
03296 write_gid_element(out, i, idn1, ide1);
03297 }
03298 }
03299 if (print_header == 0)
03300 fprintf(out, "end Elements\n");
03301
03302 print_header = 1;
03303 for (i=0; i < ne; i++)
03304 {
03305 if (elements[i].type == isolinear3d)
03306 {
03307 if (print_header)
03308 {
03309 fprintf(out, "MESH \"%ld-Brick8\" dimension 3 Elemtype Hexahedra Nnode 8\n", ide1);
03310 print_header = 0;
03311 if (print_coord)
03312 {
03313 write_gid_nodes(out, idn1);
03314 print_coord = 0;
03315 }
03316 fprintf(out, "Elements\n");
03317 }
03318 write_gid_element(out, i, idn1, ide1);
03319 }
03320 }
03321 if (print_header == 0)
03322 fprintf(out, "end Elements\n");
03323
03324 print_header = 1;
03325 for (i=0; i < ne; i++)
03326 {
03327 if (elements[i].type == isoquadratic3d)
03328 {
03329 if (print_header)
03330 {
03331 fprintf(out, "MESH \"%ld-Brick20\" dimension 3 Elemtype Hexahedra Nnode 20\n", ide1);
03332 print_header = 0;
03333 if (print_coord)
03334 {
03335 write_gid_nodes(out, idn1);
03336 print_coord = 0;
03337 }
03338 fprintf(out, "Elements\n");
03339 }
03340 write_gid_element(out, i, idn1, ide1);
03341 }
03342 }
03343 if (print_header == 0)
03344 fprintf(out, "end Elements\n");
03345 }
03346
03347
03348
03349
03350
03351
03352
03353
03354
03355
03356
03357
03358
03359 void siftop::write_gid_nodes(FILE *out, long idn1)
03360 {
03361 long i;
03362 fprintf(out, "Coordinates\n");
03363 for (i=0; i<nn; i++)
03364 fprintf(out, "%ld %e %e %e\n", i+idn1, nodes[i].x, nodes[i].y, nodes[i].z);
03365 fprintf(out, "end Coordinates\n");
03366 }
03367
03368
03369
03370
03371
03372
03373
03374
03375
03376
03377
03378
03379
03380
03381
03382 void siftop::write_gid_element(FILE *out, long i, long idn1, long ide1)
03383 {
03384 long j;
03385 fprintf(out, "%ld ", i+ide1);
03386 for (j=0; j<elements[i].nne; j++)
03387 fprintf(out, "%ld ", elements[i].nodes[j]+idn1);
03388 fprintf(out, "%ld\n", elements[i].prop);
03389 }
03390
03391
03392
03393
03394
03395
03396
03397
03398
03399
03400 void siftop::t3d2sifelformat()
03401 {
03402 long i,j;
03403 long *help;
03404 help = new long[20];
03405
03406 for(i = 0; i < ne; i++){
03407 switch(elements[i].type){
03408 case tetrahedronlinear:{
03409 for(j = 0; j < 4; j++){
03410 help[j] = elements[i].nodes[j];
03411 }
03412 elements[i].nodes[2] = help[3];
03413 elements[i].nodes[3] = help[2];
03414
03415 for(j = 0; j < 4; j++){
03416 help[j] = elements[i].propsurf[j];
03417 }
03418 elements[i].propsurf[0]= help[2];
03419 elements[i].propsurf[1]= help[3];
03420 elements[i].propsurf[2]= help[0];
03421 elements[i].propsurf[3]= help[1];
03422
03423 break;
03424 }
03425 case tetrahedronquadratic:{
03426
03427 for(j = 0; j < 10; j++){
03428 help[j] = elements[i].nodes[j];
03429 }
03430
03431 elements[i].nodes[2] = help[3];
03432 elements[i].nodes[3] = help[2];
03433 elements[i].nodes[5] = help[8];
03434 elements[i].nodes[6] = help[7];
03435 elements[i].nodes[7] = help[6];
03436 elements[i].nodes[8] = help[5];
03437
03438 for(j = 0; j < 4; j++){
03439 help[j] = elements[i].propsurf[j];
03440 }
03441 elements[i].propsurf[0]= help[2];
03442 elements[i].propsurf[1]= help[3];
03443 elements[i].propsurf[2]= help[0];
03444 elements[i].propsurf[3]= help[1];
03445
03446 break;
03447 }
03448 case isolinear3d:{
03449 for(j = 0; j < 8; j++){
03450 help[j] = elements[i].nodes[j];
03451 }
03452
03453 elements[i].nodes[0] = help[3];
03454 elements[i].nodes[1] = help[0];
03455 elements[i].nodes[2] = help[4];
03456 elements[i].nodes[3] = help[7];
03457 elements[i].nodes[4] = help[2];
03458 elements[i].nodes[5] = help[1];
03459 elements[i].nodes[6] = help[5];
03460 elements[i].nodes[7] = help[6];
03461
03462
03463 for(j = 0; j < 6; j++){
03464 help[j] = elements[i].propsurf[j];
03465 }
03466
03467 elements[i].propsurf[0] = help[0];
03468 elements[i].propsurf[1] = help[3];
03469 elements[i].propsurf[2] = help[1];
03470 elements[i].propsurf[3] = help[5];
03471 elements[i].propsurf[4] = help[2];
03472 elements[i].propsurf[5] = help[4];
03473
03474
03475 break;
03476 }
03477 case isoquadratic3d:{
03478 for(j = 0; j < 20; j++){
03479 help[j] = elements[i].nodes[j];
03480 }
03481
03482 elements[i].nodes[0] = help[3];
03483 elements[i].nodes[1] = help[0];
03484 elements[i].nodes[2] = help[4];
03485 elements[i].nodes[3] = help[7];
03486 elements[i].nodes[4] = help[2];
03487 elements[i].nodes[5] = help[1];
03488 elements[i].nodes[6] = help[5];
03489 elements[i].nodes[7] = help[6];
03490 elements[i].nodes[8] = help[11];
03491 elements[i].nodes[9] = help[16];
03492 elements[i].nodes[10] = help[15];
03493 elements[i].nodes[11] = help[19];
03494 elements[i].nodes[12] = help[10];
03495 elements[i].nodes[13] = help[8];
03496 elements[i].nodes[14] = help[12];
03497 elements[i].nodes[15] = help[14];
03498 elements[i].nodes[16] = help[9];
03499 elements[i].nodes[17] = help[17];
03500 elements[i].nodes[18] = help[13];
03501 elements[i].nodes[19] = help[18];
03502
03503 for(j = 0; j < 6; j++){
03504 help[j] = elements[i].propsurf[j];
03505 }
03506 elements[i].propsurf[0] = help[0];
03507 elements[i].propsurf[1] = help[3];
03508 elements[i].propsurf[2] = help[1];
03509 elements[i].propsurf[3] = help[5];
03510 elements[i].propsurf[4] = help[2];
03511 elements[i].propsurf[5] = help[4];
03512
03513 break;
03514 }
03515 default:{
03516 break;
03517 }
03518 }
03519 }
03520 delete []help;
03521
03522 }