00001 #include "mechprint.h"
00002 #include "global.h"
00003 #include "globmat.h"
00004 #include "alias.h"
00005 #include "gtopology.h"
00006 #include "meshtransfer.h"
00007 #include "loadcase.h"
00008 #include "element.h"
00009 #include "node.h"
00010 #include "intpoints.h"
00011 #include "plelemqq.h"
00012 #include "mathem.h"
00013 #include "genfile.h"
00014 #include "siftop_element_types.h"
00015 #include "probdesc.h"
00016 #include "siftop.h"
00017 #include "elemhead.h"
00018 #include "intp.h"
00019 #include "vector.h"
00020 #include "vecttens.h"
00021 #include "sequent.h"
00022
00023
00024
00025 #include <stdlib.h>
00026 #include <math.h>
00027 #include <string.h>
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038
00039
00040
00041
00042
00043
00044
00045
00046
00047
00048 void print_init(long istep, const char *mode, long idn1, long ide1)
00049 {
00050 char fname[FNAMELEN+20];
00051 char *path;
00052 char *name;
00053 char *suffix;
00054 long i;
00055 long sl;
00056
00057
00058 if ((Outdm->outf == NULL) && (Outdm->textout==1))
00059 {
00060 filename_decomposition (Outdm->outfn,path,name,suffix);
00061
00062 if (St == NULL)
00063 {
00064 if (istep < 0)
00065 sprintf(fname, "%s%s%s", path, name, suffix);
00066 else
00067 sprintf(fname, "%s%s.%ld%s", path, name, istep+1, suffix);
00068 }
00069 else
00070 {
00071 if (istep < 0)
00072 sprintf(fname, "%s%s.%ld%s", path, name, Mp->ns+1, suffix);
00073 else
00074 sprintf(fname, "%s%s.%ld.%ld%s", path, name, istep+1, Mp->ns+1, suffix);
00075 }
00076 Outdm->outf = fopen(fname, mode);
00077 if (Outdm->outf == NULL)
00078 {
00079 print_err("unable to open output text file '%s'", __FILE__, __LINE__, __func__, fname);
00080 abort();
00081 }
00082 fseek(Outdm->outf, 0, SEEK_END);
00083 if (ftell(Outdm->outf) == 0)
00084 Outdm->print_header(Outdm->outf);
00085
00086 delete [] path;
00087 delete [] name;
00088 delete [] suffix;
00089 }
00090
00091 if ((Outdm->gf != grfmt_no) && (Outdm->outgr == NULL))
00092 {
00093 filename_decomposition (Outdm->outgrfn,path,name,suffix);
00094 if (St == NULL)
00095 {
00096 if (istep < 0)
00097 sprintf(fname, "%s%s%s", path, name, suffix);
00098 else
00099 sprintf(fname, "%s%s.%ld%s", path, name, istep+1, suffix);
00100 }
00101 else
00102 {
00103 if (istep < 0)
00104 sprintf(fname, "%s%s.%ld%s", path, name, Mp->ns+1, suffix);
00105 else
00106 sprintf(fname, "%s%s.%ld.%ld%s", path, name, istep+1, Mp->ns+1, suffix);
00107 }
00108 if ((Outdm->gf == grfmt_gid) || (Outdm->gf == grfmt_gid_sep))
00109 {
00110 sl = strlen(fname);
00111 sprintf(fname+sl, ".msh");
00112 }
00113 if(Outdm->gf != grfmt_vtk)
00114 Outdm->outgr = fopen(fname, mode);
00115 if (Outdm->outgr == NULL && Outdm->gf != grfmt_vtk)
00116 {
00117 print_err("unable to open graphics file '%s'", __FILE__, __LINE__, __func__, fname);
00118 abort();
00119 }
00120 switch(Outdm->gf)
00121 {
00122 case grfmt_no:
00123 break;
00124 case grfmt_open_dx:
00125 break;
00126 case grfmt_femcad:
00127 export_femcad(Outdm->outgr);
00128 break;
00129 case grfmt_gid:
00130 fseek(Outdm->outgr, 0, SEEK_END);
00131 if (ftell(Outdm->outgr) == 0)
00132 export_gid_mesh(Outdm->outgr, idn1, ide1);
00133 fclose(Outdm->outgr);
00134 if (Outdm->ncut > 0)
00135 {
00136 sprintf(fname+sl, "2d.msh");
00137 Outdm->outgr = fopen(fname, mode);
00138 fseek(Outdm->outgr, 0, SEEK_END);
00139 if (ftell(Outdm->outgr) == 0)
00140 export_gid_2dmesh(Outdm->outgr, Outdm->ncut, idn1, ide1);
00141 fclose(Outdm->outgr);
00142 }
00143 sprintf(fname+sl, ".res");
00144 Outdm->outgr = fopen(fname, mode);
00145 if (Outdm->outgr == NULL)
00146 {
00147 print_err("unable to open graphics file '%s'", __FILE__, __LINE__, __func__, fname);
00148 abort();
00149 }
00150 if (strchr(mode, 'w') != NULL)
00151 {
00152 fseek(Outdm->outgr, 0, SEEK_END);
00153 if (ftell(Outdm->outgr) == 0)
00154 {
00155 fprintf(Outdm->outgr, "GiD Post Results File 1.0\n");
00156 export_gid_gauss_pt(Outdm->outgr, ide1);
00157 }
00158 }
00159 break;
00160 case grfmt_gid_sep:
00161 fseek(Outdm->outgr, 0, SEEK_END);
00162 if (ftell(Outdm->outgr) == 0)
00163 export_gid_mesh(Outdm->outgr, idn1, ide1);
00164 fclose(Outdm->outgr);
00165 if (Outdm->ncut > 0)
00166 {
00167 sprintf(fname+sl, "2d.msh");
00168 Outdm->outgr = fopen(fname, mode);
00169 fseek(Outdm->outgr, 0, SEEK_END);
00170 if (ftell(Outdm->outgr) == 0)
00171 export_gid_2dmesh(Outdm->outgr, Outdm->ncut, idn1, ide1);
00172 fclose(Outdm->outgr);
00173 }
00174 strncpy(Outdm->outgrfngs, fname, sl);
00175 Outdm->outgrfngs[sl] = 0;
00176 Outdm->outgr = NULL;
00177 Outdm->create_files_gidsp(mode);
00178 break;
00179 case grfmt_vtk:
00180 break;
00181 default:
00182 print_err("unknown type of graphics format is required", __FILE__, __LINE__, __func__);
00183 }
00184 delete [] path;
00185 delete [] name;
00186 delete [] suffix;
00187 }
00188 if ((Outdm->ndiag > 0) && (Outdm->outdiagf[0] == NULL))
00189 {
00190 filename_decomposition (Outdm->outdiagfn,path,name,suffix);
00191 for (i=0; i<Outdm->ndiag; i++)
00192 {
00193 if (St == NULL)
00194 {
00195 if (Outdm->ndiag > 1)
00196 sprintf(fname, "%s%s.%ld%s", path, name, i+1, suffix);
00197 else
00198 sprintf(fname, "%s%s%s", path, name, suffix);
00199 }
00200 else
00201 {
00202 if (Outdm->ndiag > 1)
00203 sprintf(fname, "%s%s.%ld.%ld%s", path, name, i+1, Mp->ns+1, suffix);
00204 else
00205 sprintf(fname, "%s%s.%ld%s", path, name, Mp->ns+1, suffix);
00206 }
00207 if (Outdm->outdiagf[i] == NULL)
00208 {
00209 Outdm->outdiagf[i] = fopen(fname, mode);
00210 if (Outdm->outdiagf[i] == NULL)
00211 {
00212 print_err("unable to open diagram file '%s'", __FILE__, __LINE__, __func__, fname);
00213 abort();
00214 }
00215 }
00216 }
00217 delete [] path;
00218 delete [] name;
00219 delete [] suffix;
00220 }
00221 Outdm->idn1 = idn1;
00222 Outdm->ide1 = ide1;
00223
00224
00225 }
00226
00227
00228
00229
00230
00231
00232
00233
00234
00235
00236
00237
00238
00239
00240
00241
00242
00243
00244 void print_step(long lcid, long istep, double lambda, double *fi)
00245 {
00246
00247
00248 switch(Mp->tprob)
00249 {
00250 case linear_statics:
00251 case load_balancing:
00252 case layered_linear_statics:
00253 Outdm->print_out(Outdm->outf, lcid, istep, lambda);
00254 Outdm->print_graphics(Outdm->outgr, lcid, lambda, istep, fi);
00255 break;
00256 case lin_floating_subdomain:
00257 case mat_nonlinear_statics:
00258 case geom_nonlinear_statics:
00259 case earth_pressure:
00260 case forced_dynamics:
00261 case mech_timedependent_prob:
00262 case nonlin_floating_subdomain:
00263 case growing_mech_structure:
00264 Outdm->print_newstep(Outdm->outf, lcid, istep, lambda);
00265 Outdm->print_out(Outdm->outf, lcid, istep, lambda);
00266 Outdm->print_diags(lcid, lambda, istep, fi);
00267 Outdm->print_graphics(Outdm->outgr, lcid, lambda, istep, fi);
00268 break;
00269 case eigen_dynamics:
00270 Outdm->print_newstep(Outdm->outf, lcid, istep, lambda);
00271 Outdm->print_out(Outdm->outf, lcid, istep, lambda);
00272 Outdm->print_graphics(Outdm->outgr, lcid, lambda, istep, fi);
00273 break;
00274 default:
00275 print_err("unsupported problem type is required", __FILE__, __LINE__, __func__);
00276 break;
00277 }
00278 }
00279
00280
00281
00282
00283
00284
00285
00286
00287
00288
00289
00290
00291
00292
00293
00294
00295
00296
00297
00298 void print_step_forced(long lcid, long istep, double lambda, double *fi)
00299 {
00300
00301
00302 switch(Mp->tprob)
00303 {
00304 case linear_statics:
00305 case load_balancing:
00306 case layered_linear_statics:
00307 Outdm->print_out_forced(Outdm->outf, lcid, istep, lambda);
00308 Outdm->print_graphics_forced(Outdm->outgr, lcid, lambda, istep, fi);
00309 break;
00310 case lin_floating_subdomain:
00311 case mat_nonlinear_statics:
00312 case geom_nonlinear_statics:
00313 case earth_pressure:
00314 case forced_dynamics:
00315 case mech_timedependent_prob:
00316 case nonlin_floating_subdomain:
00317 case growing_mech_structure:
00318 Outdm->print_newstep(Outdm->outf, lcid, istep, lambda);
00319 Outdm->print_out_forced(Outdm->outf, lcid, istep, lambda);
00320 Outdm->print_diags_forced(lcid, lambda, istep, fi);
00321 Outdm->print_graphics_forced(Outdm->outgr, lcid, lambda, istep, fi);
00322 break;
00323 case eigen_dynamics:
00324 Outdm->print_newstep(Outdm->outf, lcid, istep, lambda);
00325 Outdm->print_out_forced(Outdm->outf, lcid, istep, lambda);
00326 Outdm->print_graphics_forced(Outdm->outgr, lcid, lambda, istep, fi);
00327 break;
00328 default:
00329 print_err("unsupported problem type is required", __FILE__, __LINE__, __func__);
00330 break;
00331 }
00332 }
00333
00334
00335
00336
00337
00338
00339
00340
00341
00342
00343
00344
00345
00346 void print_flush()
00347 {
00348 long i;
00349 if (Outdm->outf)
00350 fflush(Outdm->outf);
00351 if (Outdm->outgr)
00352 fflush(Outdm->outgr);
00353 for (i=0; i<Outdm->ndiag; i++)
00354 {
00355 if (Outdm->outdiagf[i])
00356 fflush(Outdm->outdiagf[i]);
00357 }
00358 }
00359
00360
00361
00362
00363
00364
00365
00366
00367
00368
00369 void print_close()
00370 {
00371 long i;
00372 if (Outdm->outf)
00373 fclose(Outdm->outf);
00374 Outdm->outf = NULL;
00375 if (Outdm->outgr)
00376 fclose(Outdm->outgr);
00377 Outdm->outgr = NULL;
00378 for (i=0; i<Outdm->ndiag; i++)
00379 {
00380 if (Outdm->outdiagf[i])
00381 fclose(Outdm->outdiagf[i]);
00382 Outdm->outdiagf[i] = NULL;
00383 }
00384 }
00385
00386
00387
00388
00389
00390
00391
00392
00393
00394
00395
00396
00397
00398
00399 void export_gid_mesh(FILE *out, long idn1, long ide1)
00400 {
00401 long i, print_header, print_coord = 1;
00402
00403 print_header = 1;
00404 for (i=0; i < Mt->ne; i++)
00405 {
00406 if (Gtm->leso[i]==0)
00407 continue;
00408
00409 switch(Mt->elements[i].te)
00410 {
00411 case bar2d:
00412 case bar3d:
00413 case beam2d:
00414 case beam3d:
00415 case subsoilbeam:
00416 if (print_header)
00417 {
00418 fprintf(out, "MESH \"%ld-beams2\" dimension 3 Elemtype Linear Nnode 2\n", ide1);
00419 print_header = 0;
00420 if (print_coord)
00421 {
00422 write_gid_nodes(out, idn1);
00423 print_coord = 0;
00424 }
00425 fprintf(out, "Elements\n");
00426 }
00427 write_gid_element(out, i, idn1, ide1);
00428 break;
00429 default:
00430 break;
00431 }
00432 }
00433 if (print_header == 0)
00434 fprintf(out, "end Elements\n");
00435 print_header = 1;
00436 for (i=0; i < Mt->ne; i++)
00437 {
00438 if (Gtm->leso[i]==0)
00439 continue;
00440
00441 switch(Mt->elements[i].te)
00442 {
00443 case planequadcontact:
00444 if (print_header)
00445 {
00446 fprintf(out, "MESH \"%ld-plcontact2l\" dimension 3 Elemtype Linear Nnode 2\n", ide1);
00447 print_header = 0;
00448 if (print_coord)
00449 {
00450 write_gid_nodes(out, idn1);
00451 print_coord = 0;
00452 }
00453 fprintf(out, "Elements\n");
00454 }
00455 write_gid_contact_element(out, i, idn1, ide1);
00456 break;
00457 default:
00458 break;
00459 }
00460 }
00461 if (print_header == 0)
00462 fprintf(out, "end Elements\n");
00463 print_header = 1;
00464 for (i=0; i < Mt->ne; i++)
00465 {
00466 if (Gtm->leso[i]==0)
00467 continue;
00468
00469 switch(Mt->elements[i].te)
00470 {
00471 case barq2d:
00472 case barq3d:
00473 if (print_header)
00474 {
00475 fprintf(out, "MESH \"%ld-beams3\" dimension 3 Elemtype Linear Nnode 3\n", ide1);
00476 print_header = 0;
00477 if (print_coord)
00478 {
00479 write_gid_nodes(out, idn1);
00480 print_coord = 0;
00481 }
00482 fprintf(out, "Elements\n");
00483 }
00484 write_gid_element(out, i, idn1, ide1);
00485 break;
00486 default:
00487 break;
00488 }
00489 }
00490 if (print_header == 0)
00491 fprintf(out, "end Elements\n");
00492 print_header = 1;
00493 for (i=0; i < Mt->ne; i++)
00494 {
00495 if (Gtm->leso[i]==0)
00496 continue;
00497
00498 switch(Mt->elements[i].te)
00499 {
00500 case cctel:
00501 case dktel:
00502 case dstel:
00503 case subsoilplatetr:
00504 case axisymmlt:
00505 case planeelementlt:
00506 case planeelementrotlt:
00507 if (print_header)
00508 {
00509 fprintf(out, "MESH \"%ld-trias3\" dimension 3 Elemtype Triangle Nnode 3\n", ide1);
00510 print_header = 0;
00511 if (print_coord)
00512 {
00513 write_gid_nodes(out, idn1);
00514 print_coord = 0;
00515 }
00516 fprintf(out, "Elements\n");
00517 }
00518 write_gid_element(out, i, idn1, ide1);
00519 break;
00520 default:
00521 break;
00522 }
00523 }
00524 if (print_header == 0)
00525 fprintf(out, "end Elements\n");
00526 print_header = 1;
00527 for (i=0; i < Mt->ne; i++)
00528 {
00529 if (Gtm->leso[i]==0)
00530 continue;
00531
00532 switch(Mt->elements[i].te)
00533 {
00534 case planeelementqt:
00535 if (print_header)
00536 {
00537 fprintf(out, "MESH \"%ld-trias6\" dimension 3 Elemtype Triangle Nnode 6\n", ide1);
00538 print_header = 0;
00539 if (print_coord)
00540 {
00541 write_gid_nodes(out, idn1);
00542 print_coord = 0;
00543 }
00544 fprintf(out, "Elements\n");
00545 }
00546 write_gid_element(out, i, idn1, ide1);
00547 break;
00548 default:
00549 break;
00550 }
00551 }
00552 if (print_header == 0)
00553 fprintf(out, "end Elements\n");
00554 print_header = 1;
00555 for (i=0; i < Mt->ne; i++)
00556 {
00557 if (Gtm->leso[i]==0)
00558 continue;
00559
00560 switch(Mt->elements[i].te)
00561 {
00562 case axisymmlq:
00563 case planeelementlq:
00564 case planeelementrotlq:
00565 case q4plateel:
00566 case subsoilplateq:
00567 if (print_header)
00568 {
00569 fprintf(out, "MESH \"%ld-Quads4\" dimension 3 Elemtype Quadrilateral Nnode 4\n", ide1);
00570 print_header = 0;
00571 if (print_coord)
00572 {
00573 write_gid_nodes(out, idn1);
00574 print_coord = 0;
00575 }
00576 fprintf(out, "Elements\n");
00577 }
00578 write_gid_element(out, i, idn1, ide1);
00579 break;
00580 default:
00581 break;
00582 }
00583 }
00584 if (print_header == 0)
00585 fprintf(out, "end Elements\n");
00586 print_header = 1;
00587 for (i=0; i < Mt->ne; i++)
00588 {
00589 if (Gtm->leso[i]==0)
00590 continue;
00591
00592 switch(Mt->elements[i].te)
00593 {
00594 case planeelementqq:
00595 case axisymmqq:
00596 if (print_header)
00597 {
00598 fprintf(out, "MESH \"%ld-Quads8\" dimension 3 Elemtype Quadrilateral Nnode 8\n", ide1);
00599 print_header = 0;
00600 if (print_coord)
00601 {
00602 write_gid_nodes(out, idn1);
00603 print_coord = 0;
00604 }
00605 fprintf(out, "Elements\n");
00606 }
00607 write_gid_element(out, i, idn1, ide1);
00608 break;
00609 default:
00610 break;
00611 }
00612 }
00613 if (print_header == 0)
00614 fprintf(out, "end Elements\n");
00615 print_header = 1;
00616 for (i=0; i < Mt->ne; i++)
00617 {
00618 if (Gtm->leso[i]==0)
00619 continue;
00620
00621 switch(Mt->elements[i].te)
00622 {
00623 case lineartet:
00624 if (print_header)
00625 {
00626 fprintf(out, "MESH \"%ld-Tetras4\" dimension 3 Elemtype Tetrahedra Nnode 4\n", ide1);
00627 print_header = 0;
00628 if (print_coord)
00629 {
00630 write_gid_nodes(out, idn1);
00631 print_coord = 0;
00632 }
00633 fprintf(out, "Elements\n");
00634 }
00635 write_gid_element(out, i, idn1, ide1);
00636 break;
00637 default:
00638 break;
00639 }
00640 }
00641 if (print_header == 0)
00642 fprintf(out, "end Elements\n");
00643 print_header = 1;
00644 for (i=0; i < Mt->ne; i++)
00645 {
00646 if (Gtm->leso[i]==0)
00647 continue;
00648
00649 switch(Mt->elements[i].te)
00650 {
00651 case quadrtet:
00652 if (print_header)
00653 {
00654 fprintf(out, "MESH \"%ld-Tetras10\" dimension 3 Elemtype Tetrahedra Nnode 10\n", ide1);
00655 print_header = 0;
00656 if (print_coord)
00657 {
00658 write_gid_nodes(out, idn1);
00659 print_coord = 0;
00660 }
00661 fprintf(out, "Elements\n");
00662 }
00663 write_gid_element(out, i, idn1, ide1);
00664 break;
00665 default:
00666 break;
00667 }
00668 }
00669 if (print_header == 0)
00670 fprintf(out, "end Elements\n");
00671 print_header = 1;
00672 for (i=0; i < Mt->ne; i++)
00673 {
00674 if (Gtm->leso[i]==0)
00675 continue;
00676
00677 switch(Mt->elements[i].te)
00678 {
00679 case linearhex:
00680 if (print_header)
00681 {
00682 fprintf(out, "MESH \"%ld-Brick8\" dimension 3 Elemtype Hexahedra Nnode 8\n", ide1);
00683 print_header = 0;
00684 if (print_coord)
00685 {
00686 write_gid_nodes(out, idn1);
00687 print_coord = 0;
00688 }
00689 fprintf(out, "Elements\n");
00690 }
00691 write_gid_element(out, i, idn1, ide1);
00692 break;
00693 default:
00694 break;
00695 }
00696 }
00697 if (print_header == 0)
00698 fprintf(out, "end Elements\n");
00699 print_header = 1;
00700 for (i=0; i < Mt->ne; i++)
00701 {
00702 if (Gtm->leso[i]==0)
00703 continue;
00704
00705 switch(Mt->elements[i].te)
00706 {
00707 case quadrhex:
00708 if (print_header)
00709 {
00710 fprintf(out, "MESH \"%ld-Brick20\" dimension 3 Elemtype Hexahedra Nnode 20\n", ide1);
00711 print_header = 0;
00712 if (print_coord)
00713 {
00714 write_gid_nodes(out, idn1);
00715 print_coord = 0;
00716 }
00717 fprintf(out, "Elements\n");
00718 }
00719 write_gid_element(out, i, idn1, ide1);
00720 break;
00721 default:
00722 break;
00723 }
00724 }
00725 if (print_header == 0)
00726 fprintf(out, "end Elements\n");
00727 }
00728
00729
00730
00731
00732
00733
00734
00735
00736
00737
00738
00739
00740
00741
00742 void export_gid_gauss_pt(FILE *out, long ide1)
00743 {
00744 long i, j, k, ii, jj, brk;
00745 vector gp1, gp2, gp3, gp, w;
00746
00747
00748 brk = 0;
00749 for (i=0; i < Mt->ne; i++)
00750 {
00751 if (Gtm->leso[i]==0)
00752 continue;
00753
00754 switch(Mt->elements[i].te)
00755 {
00756 case bar2d:
00757 case bar3d:
00758 fprintf(out, "GaussPoints \"Lin_1D\" Elemtype Linear \"%ld-beams2\"\n", ide1);
00759 fprintf(out, "Number Of Gauss Points: 1\n");
00760 fprintf(out, "Nodes not included\n");
00761 fprintf(out, "Natural coordinates: internal\n");
00762 fprintf(out, "end GaussPoints\n\n");
00763 brk = 1;
00764 break;
00765 default:
00766 break;
00767 }
00768 if (brk)
00769 break;
00770 }
00771
00772 for (i=0; i < Mt->ne; i++)
00773 {
00774 if (Gtm->leso[i]==0)
00775 continue;
00776
00777 if (Mt->elements[i].te == planequadcontact)
00778 {
00779 fprintf(out, "GaussPoints \"%d\" Elemtype Linear \"%ld-plcontact2l\"\n", Mt->elements[i].te, ide1);
00780 fprintf(out, "Number Of Gauss Points: 3\n");
00781 fprintf(out, "Nodes not included\n");
00782 fprintf(out, "Natural coordinates: internal\n");
00783
00784
00785
00786
00787
00788
00789
00790
00791
00792
00793
00794
00795
00796
00797
00798
00799
00800 fprintf(out, "end GaussPoints\n\n");
00801 break;
00802 }
00803 }
00804
00805 brk = 0;
00806 for (i=0; i < Mt->ne; i++)
00807 {
00808 if (Gtm->leso[i]==0)
00809 continue;
00810
00811 switch(Mt->elements[i].te)
00812 {
00813 case beam2d:
00814 case beam3d:
00815 case subsoilbeam:
00816 fprintf(out, "GaussPoints \"Beam_1D\" Elemtype Linear \"%ld-beams2\"\n", ide1);
00817 fprintf(out, "Number Of Gauss Points: 3\n");
00818 fprintf(out, "Nodes included\n");
00819 fprintf(out, "Natural coordinates: internal\n");
00820 fprintf(out, "end GaussPoints\n\n");
00821 brk = 1;
00822 break;
00823 default:
00824 break;
00825 }
00826 if (brk)
00827 break;
00828 }
00829
00830 brk = 0;
00831 for (i=0; i < Mt->ne; i++)
00832 {
00833 if (Gtm->leso[i]==0)
00834 continue;
00835
00836 switch(Mt->elements[i].te)
00837 {
00838 case barq2d:
00839 case barq3d:
00840 fprintf(out, "GaussPoints \"Quad_1D\" Elemtype Linear \"%ld-beams3\"\n", ide1);
00841 fprintf(out, "Number Of Gauss Points: 3\n");
00842 fprintf(out, "Nodes not included\n");
00843 fprintf(out, "Natural coordinates: internal\n");
00844 fprintf(out, "end GaussPoints\n\n");
00845 brk = 1;
00846 break;
00847 default:
00848 break;
00849 }
00850 if (brk)
00851 break;
00852 }
00853
00854 for (i=0; i < Mt->ne; i++)
00855 {
00856 if (Gtm->leso[i]==0)
00857 continue;
00858
00859 if (Mt->elements[i].te == axisymmlt)
00860 {
00861 fprintf(out, "GaussPoints \"%d\" Elemtype Triangle \"%ld-trias3\"\n", axisymmlt, ide1);
00862 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
00863 fprintf(out, "Nodes not included\n");
00864 fprintf(out, "Natural coordinates: given\n");
00865 for (ii=0;ii<Asymlt->nb;ii++)
00866 {
00867 for (jj=0;jj<Asymlt->nb;jj++)
00868 {
00869 if (Asymlt->intordsm[ii][jj]==0)
00870 continue;
00871 allocv (Asymlt->intordsm[ii][jj],w);
00872 allocv (Asymlt->intordsm[ii][jj],gp1);
00873 allocv (Asymlt->intordsm[ii][jj],gp2);
00874 gauss_points_tr (gp1.a,gp2.a,w.a,Asymlt->intordsm[ii][jj]);
00875 for (i=0;i<Asymlt->intordsm[ii][jj];i++)
00876 fprintf(out, "%le %le\n", gp1[i], gp2[i]);
00877 destrv(gp1); destrv(gp2); destrv(w);
00878 }
00879 }
00880 fprintf(out, "end GaussPoints\n\n");
00881 break;
00882 }
00883 }
00884
00885 for (i=0; i < Mt->ne; i++)
00886 {
00887 if (Gtm->leso[i]==0)
00888 continue;
00889
00890 if (Mt->elements[i].te == cctel)
00891 {
00892 fprintf(out, "GaussPoints \"%d\" Elemtype Triangle \"%ld-trias3\"\n", cctel, ide1);
00893 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
00894 fprintf(out, "Nodes not included\n");
00895 fprintf(out, "Natural coordinates: given\n");
00896 for (ii=0;ii<Cct->nb;ii++)
00897 {
00898 for (jj=0;jj<Cct->nb;jj++)
00899 {
00900 if (Cct->intordsm[ii][jj]==0)
00901 continue;
00902 allocv (Cct->intordsm[ii][jj],w);
00903 allocv (Cct->intordsm[ii][jj],gp1);
00904 allocv (Cct->intordsm[ii][jj],gp2);
00905 gauss_points_tr (gp1.a,gp2.a,w.a,Cct->intordsm[ii][jj]);
00906 for (i=0;i<Cct->intordsm[ii][jj];i++)
00907 fprintf(out, "%le %le\n", gp1[i], gp2[i]);
00908 destrv(gp1); destrv(gp2); destrv(w);
00909 }
00910 }
00911 fprintf(out, "end GaussPoints\n\n");
00912 break;
00913 }
00914 }
00915
00916 for (i=0; i < Mt->ne; i++)
00917 {
00918 if (Gtm->leso[i]==0)
00919 continue;
00920
00921 if (Mt->elements[i].te == dktel)
00922 {
00923 fprintf(out, "GaussPoints \"%d\" Elemtype Triangle \"%ld-trias3\"\n", dktel, ide1);
00924 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
00925 fprintf(out, "Nodes not included\n");
00926 fprintf(out, "Natural coordinates: given\n");
00927 for (ii=0;ii<Dkt->nb;ii++)
00928 {
00929 for (jj=0;jj<Dkt->nb;jj++)
00930 {
00931 if (Dkt->intordsm[ii][jj]==0)
00932 continue;
00933 allocv (Dkt->intordsm[ii][jj],w);
00934 allocv (Dkt->intordsm[ii][jj],gp1);
00935 allocv (Dkt->intordsm[ii][jj],gp2);
00936 gauss_points_tr (gp1.a,gp2.a,w.a,Dkt->intordsm[ii][jj]);
00937 for (i=0;i<Dkt->intordsm[ii][jj];i++)
00938 fprintf(out, "%le %le\n", gp1[i], gp2[i]);
00939 destrv(gp1); destrv(gp2); destrv(w);
00940 }
00941 }
00942 fprintf(out, "end GaussPoints\n\n");
00943 break;
00944 }
00945 }
00946
00947 for (i=0; i < Mt->ne; i++)
00948 {
00949 if (Gtm->leso[i]==0)
00950 continue;
00951
00952 if (Mt->elements[i].te == subsoilplatetr)
00953 {
00954 fprintf(out, "GaussPoints \"%d\" Elemtype Triangle \"%ld-trias3\"\n", subsoilplatetr, ide1);
00955
00956 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i)-1);
00957
00958 fprintf(out, "Nodes not included\n");
00959 fprintf(out, "Natural coordinates: given\n");
00960 for (ii=0;ii<Spltr->nb;ii++)
00961 {
00962 for (jj=0;jj<Spltr->nb;jj++)
00963 {
00964 if (Spltr->intordsm[ii][jj]==0)
00965 continue;
00966 allocv (Spltr->intordsm[ii][jj],w);
00967 allocv (Spltr->intordsm[ii][jj],gp1);
00968 allocv (Spltr->intordsm[ii][jj],gp2);
00969 gauss_points_tr (gp1.a,gp2.a,w.a,Spltr->intordsm[ii][jj]);
00970
00971 for (i=1;i<Spltr->intordsm[ii][jj];i++)
00972 fprintf(out, "%le %le\n", gp1[i], gp2[i]);
00973 destrv(gp1); destrv(gp2); destrv(w);
00974 }
00975 }
00976 fprintf(out, "end GaussPoints\n\n");
00977 break;
00978 }
00979 }
00980
00981 for (i=0; i < Mt->ne; i++)
00982 {
00983 if (Gtm->leso[i]==0)
00984 continue;
00985
00986 if (Mt->elements[i].te == planeelementlt)
00987 {
00988 fprintf(out, "GaussPoints \"%d\" Elemtype Triangle\n", planeelementlt);
00989 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
00990 fprintf(out, "Nodes not included\n");
00991 fprintf(out, "Natural coordinates: internal\n");
00992 fprintf(out, "end GaussPoints\n\n");
00993 break;
00994 }
00995 }
00996
00997 for (i=0; i < Mt->ne; i++)
00998 {
00999 if (Gtm->leso[i]==0)
01000 continue;
01001
01002 if (Mt->elements[i].te == planeelementrotlt)
01003 {
01004 fprintf(out, "GaussPoints \"%d\" Elemtype Triangle \"%ld-trias3\"\n", planeelementrotlt, ide1);
01005 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01006 fprintf(out, "Nodes not included\n");
01007 fprintf(out, "Natural coordinates: given\n");
01008 for (ii=0;ii<Perlt->nb;ii++)
01009 {
01010 for (jj=0;jj<Perlt->nb;jj++)
01011 {
01012 if (Perlt->intordsm[ii][jj]==0)
01013 continue;
01014 allocv (Perlt->intordsm[ii][jj],w);
01015 allocv (Perlt->intordsm[ii][jj],gp1);
01016 allocv (Perlt->intordsm[ii][jj],gp2);
01017 gauss_points_tr (gp1.a,gp2.a,w.a,Perlt->intordsm[ii][jj]);
01018 for (i=0;i<Perlt->intordsm[ii][jj];i++)
01019 fprintf(out, "%le %le\n", gp1[i], gp2[i]);
01020 destrv(gp1); destrv(gp2); destrv(w);
01021 }
01022 }
01023 fprintf(out, "end GaussPoints\n\n");
01024 break;
01025 }
01026 }
01027
01028 for (i=0; i < Mt->ne; i++)
01029 {
01030 if (Gtm->leso[i]==0)
01031 continue;
01032
01033 if (Mt->elements[i].te == planeelementqt)
01034 {
01035 fprintf(out, "GaussPoints \"%d\" Elemtype Triangle \"%ld-trias6\"\n", planeelementqt, ide1);
01036 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01037 fprintf(out, "Nodes not included\n");
01038 fprintf(out, "Natural coordinates: given\n");
01039 for (ii=0;ii<Peqt->nb;ii++)
01040 {
01041 for (jj=0;jj<Peqt->nb;jj++)
01042 {
01043 if (Peqt->intordsm[ii][jj]==0)
01044 continue;
01045 allocv (Peqt->intordsm[ii][jj],w);
01046 allocv (Peqt->intordsm[ii][jj],gp1);
01047 allocv (Peqt->intordsm[ii][jj],gp2);
01048 gauss_points_tr (gp1.a,gp2.a,w.a,Peqt->intordsm[ii][jj]);
01049 for (i=0;i<Peqt->intordsm[ii][jj];i++)
01050 fprintf(out, "%le %le\n", gp1[i], gp2[i]);
01051 destrv(gp1); destrv(gp2); destrv(w);
01052 }
01053 }
01054 fprintf(out, "end GaussPoints\n\n");
01055 break;
01056 }
01057 }
01058
01059 for (i=0; i < Mt->ne; i++)
01060 {
01061 if (Gtm->leso[i]==0)
01062 continue;
01063
01064 if (Mt->elements[i].te == axisymmlq)
01065 {
01066 fprintf(out, "GaussPoints \"%d\" Elemtype Quadrilateral \"%ld-Quads4\"\n", axisymmlq, ide1);
01067 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01068 fprintf(out, "Nodes not included\n");
01069 fprintf(out, "Natural coordinates: given\n");
01070 for (ii=0;ii<Asymlq->nb;ii++)
01071 {
01072 for (jj=0;jj<Asymlq->nb;jj++)
01073 {
01074 if (Asymlq->intordsm[ii][jj]==0)
01075 continue;
01076 allocv (Asymlq->intordsm[ii][jj],w);
01077 allocv (Asymlq->intordsm[ii][jj],gp);
01078 gauss_points (gp.a,w.a,Asymlq->intordsm[ii][jj]);
01079 for (i=0;i<Asymlq->intordsm[ii][jj];i++)
01080 for (j=0; j<Asymlq->intordsm[ii][jj]; j++)
01081 fprintf(out, "%le %le\n", gp[i], gp[j]);
01082 destrv(gp); destrv(w);
01083 }
01084 }
01085 fprintf(out, "end GaussPoints\n\n");
01086 break;
01087 }
01088 }
01089
01090 for (i=0; i < Mt->ne; i++)
01091 {
01092 if (Gtm->leso[i]==0)
01093 continue;
01094
01095
01096
01097
01098
01099
01100
01101
01102
01103
01104
01105
01106
01107
01108
01109
01110
01111
01112
01113
01114
01115
01116
01117
01118
01119
01120
01121
01122
01123
01124 if (Mt->elements[i].te == planeelementlq)
01125 {
01126 fprintf(out, "GaussPoints \"%d\" Elemtype Quadrilateral \"%ld-Quads4\"\n", planeelementlq, ide1);
01127 fprintf(out, "Number Of Gauss Points: 4\n");
01128 fprintf(out, "Nodes not included\n");
01129 fprintf(out, "Natural coordinates: given\n");
01130 for (ii=0;ii<1;ii++)
01131 {
01132 for (jj=0;jj<1;jj++)
01133 {
01134 if (Pelq->intordsm[ii][jj]==0)
01135 continue;
01136 allocv (Pelq->intordsm[ii][jj],w);
01137 allocv (Pelq->intordsm[ii][jj],gp);
01138 gauss_points (gp.a,w.a,Pelq->intordsm[ii][jj]);
01139 for (i=0;i<Pelq->intordsm[ii][jj];i++)
01140 for (j=0; j<Pelq->intordsm[ii][jj]; j++)
01141 fprintf(out, "%le %le\n", -gp[i], -gp[j]);
01142 destrv(gp); destrv(w);
01143 }
01144 }
01145 fprintf(out, "end GaussPoints\n\n");
01146 break;
01147 }
01148 }
01149
01150 for (i=0; i < Mt->ne; i++)
01151 {
01152 if (Gtm->leso[i]==0)
01153 continue;
01154
01155 if (Mt->elements[i].te == planeelementrotlq)
01156 {
01157 fprintf(out, "GaussPoints \"%d\" Elemtype Quadrilateral \"%ld-Quads4\"\n", planeelementrotlq, ide1);
01158 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01159 fprintf(out, "Nodes not included\n");
01160 fprintf(out, "Natural coordinates: given\n");
01161 for (ii=0;ii<Perlq->nb; ii++)
01162 {
01163 for (jj=0;jj<Perlq->nb;jj++)
01164 {
01165 if (Perlq->intordsm[ii][jj]==0)
01166 continue;
01167 allocv (Perlq->intordsm[ii][jj],w);
01168 allocv (Perlq->intordsm[ii][jj],gp);
01169 gauss_points (gp.a,w.a,Perlq->intordsm[ii][jj]);
01170 for (i=0;i<Perlq->intordsm[ii][jj];i++)
01171 for (j=0; j<Perlq->intordsm[ii][jj]; j++)
01172 fprintf(out, "%le %le\n", -gp[i], -gp[j]);
01173 destrv(gp); destrv(w);
01174 }
01175 }
01176 fprintf(out, "end GaussPoints\n\n");
01177 break;
01178 }
01179 }
01180
01181 for (i=0; i < Mt->ne; i++)
01182 {
01183 if (Gtm->leso[i]==0)
01184 continue;
01185
01186 if (Mt->elements[i].te == planeelementqq)
01187 {
01188 fprintf(out, "GaussPoints \"%d\" Elemtype Quadrilateral \"%ld-Quads8\"\n", planeelementqq, ide1);
01189 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01190 fprintf(out, "Nodes not included\n");
01191 fprintf(out, "Natural coordinates: given\n");
01192 for (ii=0;ii<Peqq->nb;ii++)
01193 {
01194 for (jj=0;jj<Peqq->nb;jj++)
01195 {
01196 if (Peqq->intordsm[ii][jj]==0)
01197 continue;
01198 allocv (Peqq->intordsm[ii][jj],w);
01199 allocv (Peqq->intordsm[ii][jj],gp);
01200 gauss_points (gp.a,w.a,Peqq->intordsm[ii][jj]);
01201 for (i=0;i<Peqq->intordsm[ii][jj];i++)
01202 for (j=0; j<Peqq->intordsm[ii][jj]; j++)
01203 fprintf(out, "%le %le\n", -gp[i], -gp[j]);
01204 destrv(gp); destrv(w);
01205 }
01206 }
01207 fprintf(out, "end GaussPoints\n\n");
01208 break;
01209 }
01210 }
01211
01212 for (i=0; i < Mt->ne; i++)
01213 {
01214 if (Gtm->leso[i]==0)
01215 continue;
01216
01217 if (Mt->elements[i].te == axisymmqq)
01218 {
01219 fprintf(out, "GaussPoints \"%d\" Elemtype Quadrilateral \"%ld-Quads8\"\n", axisymmqq, ide1);
01220 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01221 fprintf(out, "Nodes not included\n");
01222 fprintf(out, "Natural coordinates: given\n");
01223 for (ii=0;ii<Asymqq->nb;ii++)
01224 {
01225 for (jj=0;jj<Asymqq->nb;jj++)
01226 {
01227 if (Asymqq->intordsm[ii][jj]==0)
01228 continue;
01229 allocv (Asymqq->intordsm[ii][jj],w);
01230 allocv (Asymqq->intordsm[ii][jj],gp);
01231 gauss_points (gp.a,w.a,Asymqq->intordsm[ii][jj]);
01232 for (i=0;i<Asymqq->intordsm[ii][jj];i++)
01233 for (j=0; j<Asymqq->intordsm[ii][jj]; j++)
01234 fprintf(out, "%le %le\n", gp[i], gp[j]);
01235 destrv(gp); destrv(w);
01236 }
01237 }
01238 fprintf(out, "end GaussPoints\n\n");
01239 break;
01240 }
01241 }
01242
01243 for (i=0; i < Mt->ne; i++)
01244 {
01245 if (Gtm->leso[i]==0)
01246 continue;
01247
01248 if (Mt->elements[i].te == lineartet)
01249 {
01250 fprintf(out, "GaussPoints \"%d\" Elemtype TetraHedra \"%ld-Tetras4\"\n", lineartet, ide1);
01251 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01252 fprintf(out, "Nodes not included\n");
01253 fprintf(out, "Natural coordinates: given\n");
01254 for (ii=0;ii<Ltet->nb;ii++)
01255 {
01256 for (jj=0;jj<Ltet->nb;jj++)
01257 {
01258 if (Ltet->intordsm[ii][jj]==0)
01259 continue;
01260 allocv (Ltet->intordsm[ii][jj],w);
01261 allocv (Ltet->intordsm[ii][jj],gp1);
01262 allocv (Ltet->intordsm[ii][jj],gp2);
01263 allocv (Ltet->intordsm[ii][jj],gp3);
01264 gauss_points_tet (gp1.a,gp2.a,gp3.a,w.a,Ltet->intordsm[ii][jj]);
01265 for (i=0;i<Ltet->intordsm[ii][jj];i++)
01266 fprintf(out, "%le %le %le\n", gp1[i], gp2[i], gp3[i]);
01267 destrv(gp1); destrv(gp2); destrv(gp3); destrv(w);
01268 }
01269 }
01270 fprintf(out, "end GaussPoints\n\n");
01271 break;
01272 }
01273 }
01274
01275 for (i=0; i < Mt->ne; i++)
01276 {
01277 if (Gtm->leso[i]==0)
01278 continue;
01279
01280 if (Mt->elements[i].te == quadrtet)
01281 {
01282 fprintf(out, "GaussPoints \"%d\" Elemtype TetraHedra \"%ld-Tetras10\"\n", quadrtet, ide1);
01283 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01284 fprintf(out, "Nodes not included\n");
01285 fprintf(out, "Natural coordinates: given\n");
01286 for (ii=0;ii<Qtet->nb;ii++)
01287 {
01288 for (jj=0;jj<Qtet->nb;jj++)
01289 {
01290 if (Qtet->intordsm[ii][jj]==0)
01291 continue;
01292 allocv (Qtet->intordsm[ii][jj],w);
01293 allocv (Qtet->intordsm[ii][jj],gp1);
01294 allocv (Qtet->intordsm[ii][jj],gp2);
01295 allocv (Qtet->intordsm[ii][jj],gp3);
01296 gauss_points_tet (gp1.a,gp2.a,gp3.a,w.a,Qtet->intordsm[ii][jj]);
01297 for (i=0;i<Qtet->intordsm[ii][jj];i++)
01298 fprintf(out, "%le %le %le\n", gp1[i], gp2[i], gp3[i]);
01299 destrv(gp1); destrv(gp2); destrv(gp3); destrv(w);
01300 }
01301 }
01302 fprintf(out, "end GaussPoints\n\n");
01303 break;
01304 }
01305 }
01306
01307 for (i=0; i < Mt->ne; i++)
01308 {
01309 if (Gtm->leso[i]==0)
01310 continue;
01311
01312 if (Mt->elements[i].te == linearhex)
01313 {
01314 fprintf(out, "GaussPoints \"%d\" Elemtype Hexahedra \"%ld-Brick8\"\n", linearhex, ide1);
01315 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01316 fprintf(out, "Nodes not included\n");
01317 fprintf(out, "Natural coordinates: given\n");
01318 for (ii=0;ii<Lhex->nb;ii++)
01319 {
01320 for (jj=0;jj<Lhex->nb;jj++)
01321 {
01322 if (Lhex->intordsm[ii][jj]==0)
01323 continue;
01324 allocv (Lhex->intordsm[ii][jj],w);
01325 allocv (Lhex->intordsm[ii][jj],gp);
01326 gauss_points(gp.a,w.a,Lhex->intordsm[ii][jj]);
01327 for (i=0;i<Lhex->intordsm[ii][jj];i++)
01328 {
01329 for (j=0;j<Lhex->intordsm[ii][jj];j++)
01330 {
01331 for (k=0;k<Lhex->intordsm[ii][jj];k++)
01332 fprintf(out, "%le %le %le\n", -gp[i], -gp[j], -gp[k]);
01333 }
01334 }
01335 destrv(gp); destrv(w);
01336 }
01337 }
01338 fprintf(out, "end GaussPoints\n\n");
01339 break;
01340 }
01341 }
01342
01343 for (i=0; i < Mt->ne; i++)
01344 {
01345 if (Gtm->leso[i]==0)
01346 continue;
01347
01348 if (Mt->elements[i].te == quadrhex)
01349 {
01350 fprintf(out, "GaussPoints \"%d\" Elemtype Hexahedra \"%ld-Brick20\"\n", quadrhex, ide1);
01351 fprintf(out, "Number Of Gauss Points: %ld\n", Mt->give_totnip(i));
01352 fprintf(out, "Nodes not included\n");
01353 fprintf(out, "Natural coordinates: given\n");
01354 for (ii=0;ii<Qhex->nb;ii++)
01355 {
01356 for (jj=0;jj<Qhex->nb;jj++)
01357 {
01358 if (Qhex->intordsm[ii][jj]==0)
01359 continue;
01360 allocv (Qhex->intordsm[ii][jj],w);
01361 allocv (Qhex->intordsm[ii][jj],gp);
01362 gauss_points(gp.a,w.a,Qhex->intordsm[ii][jj]);
01363 for (i=0;i<Qhex->intordsm[ii][jj];i++)
01364 {
01365 for (j=0;j<Qhex->intordsm[ii][jj];j++)
01366 {
01367 for (k=0;k<Qhex->intordsm[ii][jj];k++)
01368 fprintf(out, "%le %le %le\n", -gp[i], -gp[j], -gp[k]);
01369 }
01370 }
01371 destrv(gp); destrv(w);
01372 }
01373 }
01374 fprintf(out, "end GaussPoints\n\n");
01375 break;
01376 }
01377 }
01378 }
01379
01380
01381
01382
01383
01384
01385
01386
01387
01388
01389
01390
01391
01392
01393
01394 void export_gid_2dmesh(FILE *out, long icut, long idn1, long ide1)
01395 {
01396 long i, print_header, print_coord = 1;
01397 long range = Mt->ne/icut;
01398
01399 print_header = 1;
01400 for (i=0; i < Mt->ne; i++)
01401 {
01402 if (Gtm->leso[i]==0)
01403 continue;
01404
01405 switch(Mt->elements[i].te)
01406 {
01407 case linearhex:
01408 if (print_header)
01409 {
01410 fprintf(out, "MESH Quads4 dimension 3 Elemtype Quadrilateral Nnode 4\n");
01411 print_header = 0;
01412 if (print_coord)
01413 {
01414 write_gid_nodes(out, idn1);
01415 print_coord = 0;
01416 }
01417 fprintf(out, "Elements\n");
01418 }
01419 write_gid_2delement(out, i, 0, 4, i/range, 0, idn1, ide1);
01420 break;
01421 default:
01422 break;
01423 }
01424 }
01425 for (i=Mt->ne-range; i < Mt->ne; i++)
01426 {
01427 write_gid_2delement(out, i, 4, 8, icut, range, idn1, ide1);
01428 }
01429 if (print_header == 0)
01430 fprintf(out, "end Elements\n");
01431 }
01432
01433
01434
01435
01436
01437
01438
01439
01440
01441
01442
01443
01444
01445 void write_gid_nodes(FILE *out, long idn1)
01446 {
01447 long i;
01448 fprintf(out, "Coordinates\n");
01449 for (i=0; i<Mt->nn; i++)
01450 fprintf(out, "%ld %e %e %e\n", i+idn1, Gtm->gnodes[i].x, Gtm->gnodes[i].y, Gtm->gnodes[i].z);
01451 fprintf(out, "end Coordinates\n");
01452 }
01453
01454
01455
01456
01457
01458
01459
01460
01461
01462
01463
01464
01465
01466
01467
01468 void write_gid_element(FILE *out, long i, long idn1, long ide1)
01469 {
01470 long j;
01471 fprintf(out, "%ld ", i+ide1);
01472 for (j=0; j<Gtm->gelements[i].give_nne(); j++)
01473 fprintf(out, "%ld ", Gtm->gelements[i].nodes[j]+idn1);
01474 fprintf(out, "%d%ld\n", Mt->elements[i].tm[0], Mt->elements[i].idm[0]+1);
01475
01476
01477
01478
01479 }
01480
01481
01482
01483
01484
01485
01486
01487
01488
01489
01490
01491
01492
01493
01494
01495
01496 void write_gid_contact_element(FILE *out, long i, long idn1, long ide1)
01497 {
01498 long j;
01499 fprintf(out, "%ld ", i+ide1);
01500 for (j=0; j<Gtm->gelements[i].give_nne()/2; j++)
01501 fprintf(out, "%ld ", Gtm->gelements[i].nodes[j]+idn1);
01502 fprintf(out, "%ld\n", Mt->elements[i].idm[0]+1);
01503 }
01504
01505
01506
01507
01508
01509
01510
01511
01512
01513
01514
01515
01516
01517
01518
01519
01520
01521
01522
01523
01524
01525 void write_gid_2delement(FILE *out, long i, long id1, long nne, long icut, long di, long idn1, long ide1)
01526 {
01527 long j;
01528 fprintf(out, "%ld ", i+di+ide1);
01529 for (j=id1; j<nne; j++)
01530 fprintf(out, "%ld ", Gtm->gelements[i].nodes[j]+idn1);
01531 fprintf(out, "%ld%ld\n", icut+1, Mt->elements[i].idm[0]+1);
01532
01533
01534
01535 }
01536
01537
01538
01539
01540
01541
01542
01543
01544
01545
01546
01547
01548
01549
01550
01551 void write_gid_displ(FILE *out, long lcid, const char *desclcid)
01552 {
01553 long i, j;
01554 vector l,g;
01555 matrix tm(3,3);
01556 long ndofn, transf;
01557
01558 fprintf(out, "\nResult \"Displacements\" \"%ld\" %s Vector OnNodes\n", lcid, desclcid);
01559
01560 fprintf(out, "ComponentNames \"r_1\" \"r_2\" ");
01561 if (Mt->give_maxndofn() > 2)
01562 fprintf(out, "\"r_3\"");
01563 fprintf(out, "\nValues\n");
01564 for (i=0; i < Mt->nn; i++)
01565 {
01566 if (Gtm->lnso[i]==0)
01567 continue;
01568 reallocv (Mt->give_ndofn (i),l);
01569 noddispl (lcid, l.a, i);
01570
01571 transf = Mt->nodes[i].transf;
01572 if (transf>0){
01573
01574
01575
01576 ndofn = Mt->give_ndofn (i);
01577 reallocv(ndofn,g);
01578 reallocm(ndofn,ndofn,tm);
01579 fillm(0.0, tm);
01580
01581 tm[0][0]=Mt->nodes[i].e1[0];
01582 tm[1][0]=Mt->nodes[i].e1[1];
01583 if (transf == 3)
01584 tm[2][0]=Mt->nodes[i].e1[2];
01585
01586 tm[0][1]=Mt->nodes[i].e2[0];
01587 tm[1][1]=Mt->nodes[i].e2[1];
01588 if (transf == 3)
01589 tm[2][1]=Mt->nodes[i].e2[2];
01590
01591 if (transf == 3)
01592 {
01593 tm[0][2]=Mt->nodes[i].e3[0];
01594 tm[1][2]=Mt->nodes[i].e3[1];
01595 tm[2][2]=Mt->nodes[i].e3[2];
01596 }
01597 else
01598 {
01599 if (ndofn == 3)
01600 tm[2][2]=1.0;
01601 }
01602
01603 mxv(tm,l,g);
01604 copyv (g,l);
01605 }
01606
01607 fprintf(out, "%ld", i+Outdm->idn1);
01608 for (j = 0; j < Mt->give_ndofn (i); j++)
01609
01610 fprintf(out, " % e", l[j]);
01611 fprintf(out, "\n");
01612 }
01613 fprintf(out, "End Values\n");
01614 }
01615
01616
01617
01618
01619
01620
01621
01622
01623
01624
01625
01626
01627
01628
01629
01630
01631 void write_gid_nforces(FILE *out, long lcid, const char *desclcid, double *ifor)
01632 {
01633 long i, j, ndof;
01634 vector f,g;
01635 matrix tm(3,3);
01636 long ndofn, transf;
01637
01638 fprintf(out, "\nResult \"Forces\" \"%ld\" %s Vector OnNodes\n", lcid, desclcid);
01639 fprintf(out, "ComponentNames \"f_1\" \"f_2\" \"f_3\"\nValues\n");
01640 for (i=0; i < Mt->nn; i++)
01641 {
01642 if (Gtm->lnso[i]==0)
01643 continue;
01644
01645 ndof = Mt->give_ndofn(i);
01646 reallocv (ndof, f);
01647 nodforce (ifor, i, f);
01648
01649 transf = Mt->nodes[i].transf;
01650 if (transf>0){
01651
01652
01653
01654 ndofn = Mt->give_ndofn (i);
01655 reallocv(ndofn, g);
01656 reallocm(ndofn,ndofn,tm);
01657 fillm(0.0, tm);
01658
01659 tm[0][0]=Mt->nodes[i].e1[0];
01660 tm[1][0]=Mt->nodes[i].e1[1];
01661 if (transf == 3)
01662 tm[2][0]=Mt->nodes[i].e1[2];
01663
01664 tm[0][1]=Mt->nodes[i].e2[0];
01665 tm[1][1]=Mt->nodes[i].e2[1];
01666 if (transf == 3)
01667 tm[2][1]=Mt->nodes[i].e2[2];
01668
01669 if (transf == 3)
01670 {
01671 tm[0][2]=Mt->nodes[i].e3[0];
01672 tm[1][2]=Mt->nodes[i].e3[1];
01673 tm[2][2]=Mt->nodes[i].e3[2];
01674 }
01675 else
01676 {
01677 if (ndofn == 3)
01678 tm[2][2]=1.0;
01679 }
01680
01681 mxv(tm,f,g);
01682 copyv (g,f);
01683 }
01684
01685 fprintf(out, "%ld", i+Outdm->idn1);
01686 if (ndof > 3) ndof = 3;
01687 for (j = 0; j < ndof; j++)
01688 fprintf(out, " % e", f[j]);
01689
01690 for (j = ndof; j < 3; j++)
01691 fprintf(out, " 0.0");
01692 fprintf(out, "\n");
01693 }
01694 fprintf(out, "End Values\n");
01695 }
01696
01697
01698
01699
01700
01701
01702
01703
01704
01705
01706
01707
01708
01709
01710
01711
01712
01713 void write_gid_nodscalar(FILE *out, strastre scal, long lcid, long dir, const char *desclcid)
01714 {
01715 const char *sig = "";
01716 long i, ncompstr, ncompother;
01717 double aux;
01718 long ir;
01719 char *namevar[] = {"x", "y", "z"};
01720 matrix t;
01721 vector str, p;
01722
01723 if(Outdm->nlcs > 0)
01724 {
01725 allocv(3, p);
01726 allocm(3, 3, t);
01727 }
01728
01729 if (scal == stress)
01730 {
01731 if (Mm->max_ncompstrn < 6)
01732 {
01733 switch (dir)
01734 {
01735 case 0:
01736 sig = "sig_n_1";
01737 break;
01738 case 1:
01739 sig = "sig_n_2";
01740 break;
01741 case 2:
01742 sig = "sig_n_3";
01743 break;
01744 case 3:
01745 sig = "sig_n_4";
01746 break;
01747 case 4:
01748 sig = "sig_n_5";
01749 break;
01750 default:
01751 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
01752 }
01753 }
01754 if (Mm->max_ncompstrn == 6)
01755 {
01756 switch (dir)
01757 {
01758 case 0:
01759 sig = "sig_n_x";
01760 break;
01761 case 1:
01762 sig = "sig_n_y";
01763 break;
01764 case 2:
01765 sig = "sig_n_z";
01766 break;
01767 case 3:
01768 sig = "tau_n_yz";
01769 break;
01770 case 4:
01771 sig = "tau_n_xz";
01772 break;
01773 case 5:
01774 sig = "tau_n_xy";
01775 break;
01776 default:
01777 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
01778 }
01779 }
01780 }
01781 if (scal == strain)
01782 {
01783 if (Mm->max_ncompstrn < 6)
01784 {
01785 switch (dir)
01786 {
01787 case 0:
01788 sig = "eps_n_1";
01789 break;
01790 case 1:
01791 sig = "eps_n_2";
01792 break;
01793 case 2:
01794 sig = "eps_n_3";
01795 break;
01796 case 3:
01797 sig = "eps_n_4";
01798 break;
01799 case 4:
01800 sig = "eps_n_5";
01801 break;
01802 default:
01803 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
01804 }
01805 }
01806 if (Mm->max_ncompstrn == 6)
01807 {
01808 switch (dir)
01809 {
01810 case 0:
01811 sig = "eps_n_x";
01812 break;
01813 case 1:
01814 sig = "eps_n_y";
01815 break;
01816 case 2:
01817 sig = "eps_n_z";
01818 break;
01819 case 3:
01820 sig = "eps_n_yz";
01821 break;
01822 case 4:
01823 sig = "eps_n_xz";
01824 break;
01825 case 5:
01826 sig = "eps_n_xy";
01827 break;
01828 default:
01829 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
01830 }
01831 }
01832 }
01833 if (scal == pstress)
01834 {
01835 switch (dir)
01836 {
01837 case 0:
01838 sig = "psig_n_1";
01839 break;
01840 case 1:
01841 sig = "psig_n_2";
01842 break;
01843 case 2:
01844 sig = "psig_n_3";
01845 break;
01846 case 3:
01847 sig = "tau_n_max";
01848 break;
01849 default:
01850 print_err("unknown direction of principal stresses is required", __FILE__, __LINE__, __func__);
01851 }
01852 }
01853 if (scal == pstrain)
01854 {
01855 switch (dir)
01856 {
01857 case 0:
01858 sig = "peps_n_1";
01859 break;
01860 case 1:
01861 sig = "peps_n_2";
01862 break;
01863 case 3:
01864 sig = "peps_n_3";
01865 break;
01866 default:
01867 print_err("unknown direction of principal strains is required", __FILE__, __LINE__, __func__);
01868 }
01869 }
01870 if (scal == other)
01871 fprintf(out, "\nResult \"other_n_%ld\" \"%ld\" %s Scalar OnNodes\n", dir+1, lcid, desclcid);
01872 else
01873 fprintf(out, "\nResult \"%s\" \"%ld\" %s Scalar OnNodes\n", sig, lcid, desclcid);
01874 fprintf(out, "Values\n");
01875 switch (scal)
01876 {
01877 case stress:
01878 for (i = 0; i < Mt->nn; i++)
01879 {
01880 if (Gtm->lnso[i]==0)
01881 continue;
01882
01883 if (Outdm->nog.selnstre.presence_id(Outdm->nog.selstre, i, dir, ir) == 0)
01884 continue;
01885 ncompstr = Mt->nodes[i].ncompstr;
01886 if (Outdm->nog.transtre[ir] > 0)
01887 {
01888 Mt->give_nodal_coord(i, p);
01889 Outdm->lcs[Outdm->nog.transtre[ir]-1].evaluate_t(p, namevar, t);
01890 str.n = ncompstr;
01891 str.a = Mt->nodes[i].stress+ncompstr*lcid;
01892 if (ncompstr == 4)
01893 gl_comp_engvectortransf(str, aux, dir, t, planestrain, stress);
01894 if (ncompstr == 6)
01895 gl_comp_engvectortransf(str, aux, dir, t, spacestress, stress);
01896 str.a = NULL;
01897 fprintf(out, "%ld % e\n", i+Outdm->idn1, aux);
01898 }
01899 else
01900 fprintf(out, "%ld % e\n", i+Outdm->idn1, Mt->nodes[i].stress[ncompstr*lcid+dir]);
01901 }
01902 break;
01903 case strain:
01904 for (i = 0; i < Mt->nn; i++)
01905 {
01906 if (Gtm->lnso[i]==0)
01907 continue;
01908
01909 if (Outdm->nog.selnstra.presence_id(Outdm->nog.selstra, i, dir, ir) == 0)
01910 continue;
01911 ncompstr = Mt->nodes[i].ncompstr;
01912 if (Outdm->nog.transtra[ir] > 0)
01913 {
01914 Mt->give_nodal_coord(i, p);
01915 Outdm->lcs[Outdm->nog.transtra[ir]-1].evaluate_t(p, namevar, t);
01916 str.n = ncompstr;
01917 str.a = Mt->nodes[i].strain+ncompstr*lcid;
01918 if (ncompstr == 4)
01919 gl_comp_engvectortransf(str, aux, dir, t, planestrain, strain);
01920 if (ncompstr == 6)
01921 gl_comp_engvectortransf(str, aux, dir, t, spacestress, strain);
01922 str.a = NULL;
01923 fprintf(out, "%ld % e\n", i+Outdm->idn1, aux);
01924 }
01925 fprintf(out, "%ld % e\n", i+Outdm->idn1, Mt->nodes[i].strain[ncompstr*lcid+dir]);
01926 }
01927 break;
01928 case pstress:
01929 for (i = 0; i < Mt->nn; i++)
01930 {
01931 if (Gtm->lnso[i]==0)
01932 continue;
01933 if (Mt->nodes[i].pstre == NULL)
01934 continue;
01935
01936 if (Outdm->nog.selnstre.presence_id(Outdm->nog.selstre, i, dir) == 0)
01937 continue;
01938 if (dir == 3)
01939 fprintf(out, "%ld % e\n", i+Outdm->idn1, (Mt->nodes[i].pstre[3*lcid]+Mt->nodes[i].pstre[3*lcid+2])/2);
01940 else
01941 fprintf(out, "%ld % e\n", i+Outdm->idn1, Mt->nodes[i].pstre[3*lcid+dir]);
01942 }
01943 break;
01944 case pstrain:
01945 for (i = 0; i < Mt->nn; i++)
01946 {
01947 if (Gtm->lnso[i]==0)
01948 continue;
01949 if (Mt->nodes[i].pstra == NULL)
01950 continue;
01951
01952 if (Outdm->nog.selnstra.presence_id(Outdm->nog.selstra, i, dir) == 0)
01953 continue;
01954 ncompstr = Mt->nodes[i].ncompstr;
01955 fprintf(out, "%ld % e\n", i+Outdm->idn1, Mt->nodes[i].pstra[3*lcid+dir]);
01956 }
01957 break;
01958 case other:
01959 for (i = 0; i < Mt->nn; i++)
01960 {
01961 if (Gtm->lnso[i]==0)
01962 continue;
01963
01964 if (Outdm->nog.selnoth.presence_id(Outdm->nog.seloth, i, dir) == 0)
01965 continue;
01966 ncompother = Mt->nodes[i].ncompother;
01967 fprintf(out, "%ld % e\n", i+Outdm->idn1, Mt->nodes[i].other[ncompother*lcid+dir]);
01968 }
01969 break;
01970 default:
01971 print_err("unknown value type is required", __FILE__, __LINE__, __func__);
01972 }
01973 fprintf(out, "End Values\n");
01974 }
01975
01976
01977
01978
01979
01980
01981
01982
01983
01984
01985
01986
01987
01988
01989
01990
01991
01992 void write_gid_nodvector(FILE *out, strastre q, long lcid, long sid, const char *desclcid)
01993 {
01994 char sig[70];
01995 long i, j, ncompstr, ncompother;
01996 long asid, q_id1, q_n;
01997
01998 switch (q)
01999 {
02000 case strain:
02001 q_id1 = Outdm->nog.selstra[sid].id1[0];
02002 q_n = Outdm->nog.selstra[sid].ncomp[0];
02003 sprintf(sig, "eps_n_v%ld-%ld_s%ld", q_id1+1, q_n, sid+1);
02004 break;
02005 case stress:
02006 q_id1 = Outdm->nog.selstre[sid].id1[0];
02007 q_n = Outdm->nog.selstre[sid].ncomp[0];
02008 sprintf(sig, "sig_n_v%ld-%ld_s%ld", q_id1+1, q_n, sid+1);
02009 break;
02010 case other:
02011 q_id1 = Outdm->nog.seloth[sid].id1[0];
02012 q_n = Outdm->nog.seloth[sid].ncomp[0];
02013 sprintf(sig, "other_n_v%ld-%ld_s%ld", q_id1+1, q_n, sid+1);
02014 break;
02015 default:
02016 print_err("unsupported type of quantity (%d) is required", __FILE__, __LINE__, __func__, q);
02017 abort();
02018 }
02019 fprintf(out, "\nResult \"%s\" \"%ld\" %s Vector OnNodes\n", sig, lcid, desclcid);
02020 fprintf(out, "ComponentNames \"%s_1\" \"%s_2\" \"%s_3\"\nValues\n", sig, sig, sig);
02021 switch (q)
02022 {
02023 case strain:
02024 for (i = 0; i < Mt->nn; i++)
02025 {
02026 if (Gtm->lnso[i]==0)
02027 continue;
02028
02029 if (Outdm->nog.selnstra.presence_id(i, sid, asid) == 0)
02030 continue;
02031 ncompstr = Mt->nodes[i].ncompstr;
02032 fprintf(out, "%ld", i+Outdm->idn1);
02033 if (q_id1+q_n > ncompstr)
02034 {
02035 print_err("index of selected vector component exceeds total number of strain components\n"
02036 " (node=%ld, ncompstr=%ld, id1=%ld, n=%ld)\n", __FILE__, __LINE__, __func__, i+1, ncompstr, q_id1+1, q_n);
02037 abort();
02038 }
02039 for (j = q_id1; j < q_id1+q_n; j++)
02040 fprintf(out, " % e", Mt->nodes[i].strain[ncompstr*lcid+j]);
02041 for (j = 0; j < 3-q_n; j++)
02042 fprintf(out, " 0.0");
02043 fprintf(out, "\n");
02044 }
02045 break;
02046 case stress:
02047 for (i = 0; i < Mt->nn; i++)
02048 {
02049 if (Gtm->lnso[i]==0)
02050 continue;
02051
02052 if (Outdm->nog.selnstre.presence_id(i, sid, asid) == 0)
02053 continue;
02054 ncompstr = Mt->nodes[i].ncompstr;
02055 fprintf(out, "%ld", i+Outdm->idn1);
02056 if (q_id1+q_n > ncompstr)
02057 {
02058 print_err("index of selected vector component exceeds total number of stress component\n"
02059 " (node=%ld, ncompstr=%ld, id1=%ld, n=%ld)\n", __FILE__, __LINE__, __func__, i+1, ncompstr, q_id1+1, q_n);
02060 abort();
02061 }
02062 for (j = q_id1; j < q_id1+q_n; j++)
02063 fprintf(out, " % e", Mt->nodes[i].stress[ncompstr*lcid+j]);
02064 for (j = 0; j < 3-q_n; j++)
02065 fprintf(out, " 0.0");
02066 fprintf(out, "\n");
02067 }
02068 break;
02069 case other:
02070 for (i = 0; i < Mt->nn; i++)
02071 {
02072 if (Gtm->lnso[i]==0)
02073 continue;
02074
02075 if (Outdm->nog.selnoth.presence_id(i, sid, asid) == 0)
02076 continue;
02077 ncompother = Mt->nodes[i].ncompother;
02078 fprintf(out, "%ld", i+Outdm->idn1);
02079 if (q_id1+q_n > ncompother)
02080 {
02081 print_err("index of selected vector component exceeds total number of other values component\n"
02082 " (node=%ld, ncompother=%ld, id1=%ld, n=%ld)\n", __FILE__, __LINE__, __func__, i+1, ncompother, q_id1+1, q_n);
02083 abort();
02084 }
02085 for (j = q_id1; j < q_id1+q_n; j++)
02086 fprintf(out, " % e", Mt->nodes[i].other[ncompother*lcid+j]);
02087 for (j = 0; j < 3-q_n; j++)
02088 fprintf(out, " 0.0");
02089 fprintf(out, "\n");
02090 }
02091 break;
02092 default:
02093 print_err("unknown quantity type is required", __FILE__, __LINE__, __func__);
02094 abort();
02095 }
02096 fprintf(out, "End Values\n");
02097 }
02098
02099
02100
02101
02102
02103
02104
02105
02106
02107
02108
02109
02110
02111
02112
02113
02114
02115
02116
02117 void write_gid_nodtensor(FILE *out, strastre q, long lcid, long sid, const char *desclcid)
02118 {
02119 char sig[70];
02120 long i, ncompstr, ncompother;
02121 long asid, q_id1, q_n;
02122 char *namevar[] = {"x", "y", "z"};
02123 matrix t(3,3);
02124 matrix tmat;
02125 vector str, p;
02126 vector v;
02127 strastrestate ssst;
02128
02129 if(Outdm->nlcs > 0)
02130 {
02131 allocv(3, p);
02132 allocm(3, 3, tmat);
02133 }
02134
02135 switch (q)
02136 {
02137 case strain:
02138 q_id1 = Outdm->nog.selstra[sid].id1[0];
02139 q_n = 0;
02140 sprintf(sig, "eps_n_m_s%ld", sid+1);
02141 break;
02142 case stress:
02143 q_id1 = Outdm->nog.selstre[sid].id1[0];
02144 q_n = 0;
02145 sprintf(sig, "sig_n_m_s%ld", sid+1);
02146 break;
02147 case other:
02148 q_id1 = Outdm->nog.seloth[sid].id1[0];
02149 q_n = Outdm->nog.seloth[sid].ncomp[0];
02150 sprintf(sig, "other_n_m_s%ld", sid+1);
02151 break;
02152 default:
02153 print_err("unsupported type of quantity (%d) is required", __FILE__, __LINE__, __func__, q);
02154 abort();
02155 }
02156
02157 if (q != other)
02158 {
02159 for (i = 0; i < Mt->nn; i++)
02160 {
02161 if (Gtm->lnso[i]==0)
02162 continue;
02163 switch (q)
02164 {
02165 case strain:
02166
02167 if (Outdm->nog.selnstra.presence_id(i, sid, asid) == 0)
02168 continue;
02169 if (q_n < Mt->nodes[i].ncompstr)
02170 q_n = Mt->nodes[i].ncompstr;
02171 break;
02172 case stress:
02173
02174 if (Outdm->nog.selnstre.presence_id(i, sid, asid) == 0)
02175 continue;
02176 if (q_n < Mt->nodes[i].ncompstr)
02177 q_n = Mt->nodes[i].ncompstr;
02178 break;
02179 default:
02180 break;
02181 }
02182 }
02183 if (q_n == 0)
02184 return;
02185 }
02186 if ((q_n != 6) && (q_n != 4) && (q_n != 1))
02187 {
02188 print_err("wrong number of required tensor components (%ld)\n"
02189 " only 1 or 4 or 6 components are allowed.", __FILE__, __LINE__, __func__, q_n);
02190 abort();
02191 }
02192 if (q_n == 6)
02193 fprintf(out, "\nResult \"%s\" \"%ld\" %s Matrix OnNodes \n", sig, lcid, desclcid);
02194 else
02195 fprintf(out, "\nResult \"%s\" \"%ld\" %s PlainDeformationMatrix OnNodes \n", sig, lcid, desclcid);
02196 fprintf(out, "Values\n");
02197
02198 switch (q)
02199 {
02200 case strain:
02201 for (i = 0; i < Mt->nn; i++)
02202 {
02203 if (Gtm->lnso[i]==0)
02204 continue;
02205
02206 if (Outdm->nog.selnstra.presence_id(i, sid, asid) == 0)
02207 continue;
02208 ncompstr = Mt->nodes[i].ncompstr;
02209 fprintf(out, "%ld", i+Outdm->idn1);
02210 allocv(ncompstr, v);
02211 copyv(Mt->nodes[i].strain+ncompstr*lcid+q_id1, v);
02212 ssst = guess_ssst(ncompstr);
02213 vector_tensor (v, t, ssst, strain);
02214 if (Outdm->nog.transtra[sid] > 0)
02215 {
02216 Mt->give_nodal_coord(i, p);
02217 Outdm->lcs[Outdm->nog.transtra[sid]-1].evaluate_t(p, namevar, tmat);
02218 glmatrixtransf(t, tmat);
02219 }
02220 if (q_n == 6)
02221 fprintf(out, " % e % e % e % e % e % e\n", t[0][0], t[1][1], t[2][2], t[0][1], t[1][2], t[0][2]);
02222 else
02223 fprintf(out, " % e % e % e % e\n", t[0][0], t[1][1], t[0][1], t[2][2]);
02224 destrv(v);
02225 }
02226 break;
02227 case stress:
02228 for (i = 0; i < Mt->nn; i++)
02229 {
02230 if (Gtm->lnso[i]==0)
02231 continue;
02232
02233 if (Outdm->nog.selnstre.presence_id(i, sid, asid) == 0)
02234 continue;
02235 ncompstr = Mt->nodes[i].ncompstr;
02236 fprintf(out, "%ld", i+Outdm->idn1);
02237 allocv(ncompstr, v);
02238 copyv(Mt->nodes[i].stress+ncompstr*lcid+q_id1, v);
02239 ssst = guess_ssst(ncompstr);
02240 vector_tensor (v, t, ssst, strain);
02241 if (Outdm->nog.transtre[sid] > 0)
02242 {
02243 Mt->give_nodal_coord(i, p);
02244 Outdm->lcs[Outdm->nog.transtre[sid]-1].evaluate_t(p, namevar, tmat);
02245 glmatrixtransf(t, tmat);
02246 }
02247 if (q_n == 6)
02248 fprintf(out, " % e % e % e % e % e % e\n", t[0][0], t[1][1], t[2][2], t[0][1], t[1][2], t[0][2]);
02249 else
02250 fprintf(out, " % e % e % e % e\n", t[0][0], t[1][1], t[0][1], t[2][2]);
02251 destrv(v);
02252 }
02253 break;
02254 case other:
02255 for (i = 0; i < Mt->nn; i++)
02256 {
02257 if (Gtm->lnso[i]==0)
02258 continue;
02259
02260 if (Outdm->nog.selnoth.presence_id(i, sid, asid) == 0)
02261 continue;
02262 ncompother = Mt->nodes[i].ncompother;
02263 fprintf(out, "%ld", i+Outdm->idn1);
02264 allocv(q_n, v);
02265 copyv(Mt->nodes[i].other+ncompother*lcid+q_id1, v);
02266 ssst = guess_ssst(q_n);
02267 vector_tensor (v, t, ssst, strain);
02268 if (q_n == 6)
02269 fprintf(out, " % e % e % e % e % e % e\n", t[0][0], t[1][1], t[2][2], t[0][1], t[1][2], t[0][2]);
02270 else
02271 fprintf(out, " % e % e % e % e\n", t[0][0], t[1][1], t[0][1], t[2][2]);
02272 destrv(v);
02273 }
02274 break;
02275 default:
02276 print_err("unknown quantity type is required", __FILE__, __LINE__, __func__);
02277 abort();
02278 }
02279 fprintf(out, "End Values\n");
02280 }
02281
02282
02283
02284
02285
02286
02287
02288
02289
02290
02291
02292
02293
02294
02295
02296
02297
02298 void write_gid_elemscalar(FILE *out, strastre scal, long lcid, long dir, const char *desclcid)
02299 {
02300 if (Bar2d)
02301 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, bar2d);
02302
02303 if (Bar3d)
02304 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, bar3d);
02305
02306 if (Barq2d)
02307 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, barq2d);
02308
02309 if (Barq3d)
02310 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, barq3d);
02311
02312 if (Beam2d)
02313 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, beam2d);
02314
02315 if (Beam3d)
02316 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, beam3d);
02317
02318 if (Sbeam)
02319 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, subsoilbeam);
02320
02321 if (Spring)
02322 {
02323 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, spring_1);
02324 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, spring_2);
02325 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, spring_3);
02326 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, spring_4);
02327 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, spring_5);
02328 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, spring_6);
02329 }
02330
02331 if (Pqcon)
02332 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, planequadcontact);
02333
02334 if (Pelt)
02335 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, planeelementlt);
02336
02337 if (Peqt)
02338 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, planeelementqt);
02339
02340 if (Perlt)
02341 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, planeelementrotlt);
02342
02343 if (Pelq)
02344 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, planeelementlq);
02345
02346 if (Peqq)
02347 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, planeelementqq);
02348
02349 if (Perlq)
02350 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, planeelementrotlq);
02351
02352 if (Pesqt)
02353 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, planeelementsubqt);
02354
02355 if (Cct)
02356 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, cctel);
02357
02358 if (Dkt)
02359 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, dktel);
02360
02361 if (Dst)
02362 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, dstel);
02363
02364 if (Q4pl)
02365 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, q4plateel);
02366
02367 if (Spltr)
02368 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, subsoilplatetr);
02369
02370 if (Splq)
02371 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, subsoilplateq);
02372
02373 if (Asymlq)
02374 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, axisymmlq);
02375
02376 if (Asymlt)
02377 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, axisymmlt);
02378
02379 if (Asymqq)
02380 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, axisymmqq);
02381
02382 if (Shtr)
02383 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, shelltrelem);
02384
02385 if (Shq)
02386 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, shellqelem);
02387
02388 if (Ltet)
02389 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, lineartet);
02390
02391 if (Qtet)
02392 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, quadrtet);
02393
02394 if (Lhex)
02395 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, linearhex);
02396
02397 if (Qhex)
02398 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, quadrhex);
02399
02400 if (Lwed)
02401 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, linearwed);
02402
02403 if (Qwed)
02404 write_gid_elem_type_scalar(out, scal, lcid, dir, desclcid, quadrwed);
02405 }
02406
02407
02408
02409
02410
02411
02412
02413
02414
02415
02416
02417
02418
02419
02420
02421
02422
02423
02424
02425 void write_gid_elem_type_scalar(FILE *out, strastre scal, long lcid, long dir, const char *desclcid, elemtype te)
02426 {
02427 const char *sig = "";
02428 char gpname[1000];
02429 long i, j, ipp,ncompstr,ncompother,tnip, ir;
02430 long print_header = 1;
02431 long num_copy_ip;
02432 char *namevar[] = {"x", "y", "z"};
02433 matrix t;
02434 vector str, p;
02435 double aux;
02436
02437 if(Outdm->nlcs > 0)
02438 {
02439 allocv(3, p);
02440 allocm(3, 3, t);
02441 }
02442
02443
02444
02445 for (i = 0; i < Mt->ne; i++)
02446 {
02447 if (Gtm->leso[i]==0)
02448 continue;
02449
02450 tnip = Mt->give_totnip(i);
02451 if (te == planeelementlq)
02452 tnip = 4;
02453 num_copy_ip = 0;
02454 switch (te)
02455 {
02456 case beam3d:
02457 case subsoilbeam:
02458 num_copy_ip = 1;
02459 break;
02460 default:
02461 break;
02462 }
02463
02464
02465 if (Mt->elements[i].te != te)
02466 continue;
02467
02468
02469 if ((scal == strain) && (Outdm->eog.selestra.presence_id(Outdm->eog.selstra, i, dir, ir) == 0))
02470 continue;
02471 if ((scal == stress) && (Outdm->eog.selestre.presence_id(Outdm->eog.selstre, i, dir, ir) == 0))
02472 continue;
02473 if ((scal == other) && (Outdm->eog.seleoth.presence_id(Outdm->eog.seloth, i, dir, ir) == 0))
02474 continue;
02475
02476
02477
02478 ipp = Mt->elements[i].ipp[0][0];
02479 if ((scal == strain) && (Outdm->eog.selstra[ir].st == sel_all) && (dir >= Mm->ip[ipp].ncompstr))
02480 continue;
02481 if ((scal == stress) && (Outdm->eog.selstre[ir].st == sel_all) && (dir >= Mm->ip[ipp].ncompstr))
02482 continue;
02483 if ((scal == other) && (Outdm->eog.seloth[ir].st == sel_all) && (dir >= Mm->ip[ipp].ncompeqother))
02484 continue;
02485
02486 if (print_header)
02487 {
02488 if (scal == stress)
02489 {
02490 switch (Mm->ip[ipp].ssst)
02491 {
02492 case bar:
02493 switch (dir)
02494 {
02495 case 0:
02496 sig = "sig_e_1";
02497 break;
02498 case 1:
02499 sig = "sig_e_2";
02500 break;
02501 case 2:
02502 sig = "sig_e_3";
02503 break;
02504 case 3:
02505 sig = "sig_e_4";
02506 break;
02507 case 4:
02508 sig = "sig_e_5";
02509 break;
02510 case 5:
02511 sig = "sig_e_6";
02512 break;
02513 default:
02514 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
02515 }
02516 break;
02517 case plbeam:
02518 switch (dir)
02519 {
02520 case 0:
02521 sig = "N";
02522 break;
02523 case 1:
02524 sig = "V";
02525 break;
02526 case 2:
02527 sig = "M";
02528 break;
02529 default:
02530 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
02531 }
02532 break;
02533 case planestress:
02534 case planestrain:
02535 switch (dir)
02536 {
02537 case 0:
02538 sig = "sig_e_x";
02539 break;
02540 case 1:
02541 sig = "sig_e_y";
02542 break;
02543 case 2:
02544 sig = "tau_e_xy";
02545 break;
02546 case 3:
02547 sig = "sig_e_z";
02548 break;
02549 default:
02550 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
02551 }
02552 break;
02553
02554 case platek:{
02555 switch (dir){
02556 case 0:
02557 sig = "m_x";
02558 break;
02559 case 1:
02560 sig = "m_y";
02561 break;
02562 case 2:
02563 sig = "m_xy";
02564 break;
02565 default:
02566 print_err("unknown direction of stress is required", __FILE__, __LINE__, __func__);
02567 }
02568 break;
02569 }
02570
02571 case plates:{
02572 switch (dir){
02573 case 0:
02574 sig = "m_x";
02575 break;
02576 case 1:
02577 sig = "m_y";
02578 break;
02579 case 2:
02580 sig = "m_xy";
02581 break;
02582 case 3:
02583 sig = "q_x";
02584 break;
02585 case 4:
02586 sig = "q_y";
02587 break;
02588 default:
02589 print_err("unknown direction of stress is required", __FILE__, __LINE__, __func__);
02590 }
02591 break;
02592 }
02593
02594 case axisymm:
02595 switch (dir)
02596 {
02597 case 0:
02598 sig = "sig_e_r";
02599 break;
02600 case 1:
02601 sig = "sig_e_z";
02602 break;
02603 case 2:
02604 sig = "sig_e_phi";
02605 break;
02606 case 3:
02607 sig = "sig_e_rz";
02608 break;
02609 default:
02610 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
02611 }
02612 break;
02613
02614 case shell:{
02615 switch (dir){
02616 case 0:
02617 sig = "sig_x";
02618 break;
02619 case 1:
02620 sig = "sig_y";
02621 break;
02622 case 2:
02623 sig = "tau_z";
02624 break;
02625 case 3:
02626 sig = "m_x";
02627 break;
02628 case 4:
02629 sig = "m_y";
02630 break;
02631 case 5:
02632 sig = "m_xy";
02633 break;
02634 default:
02635 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
02636 }
02637 break;
02638 }
02639
02640 case spacestress:
02641 switch (dir)
02642 {
02643 case 0:
02644 sig = "sig_e_x";
02645 break;
02646 case 1:
02647 sig = "sig_e_y";
02648 break;
02649 case 2:
02650 sig = "sig_e_z";
02651 break;
02652 case 3:
02653 sig = "tau_e_yz";
02654 break;
02655 case 4:
02656 sig = "tau_e_xz";
02657 break;
02658 case 5:
02659 sig = "tau_e_xy";
02660 break;
02661 default:
02662 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
02663 }
02664 break;
02665 default:
02666 print_err("unknown stress/strain state is required", __FILE__, __LINE__, __func__);
02667 }
02668 }
02669 if (scal == strain)
02670 {
02671 switch (Mm->ip[ipp].ssst)
02672 {
02673 case bar:
02674 switch (dir)
02675 {
02676 case 0:
02677 sig = "eps_e_1";
02678 break;
02679 case 1:
02680 sig = "eps_e_2";
02681 break;
02682 case 2:
02683 sig = "eps_e_3";
02684 break;
02685 case 3:
02686 sig = "eps_e_4";
02687 break;
02688 case 4:
02689 sig = "eps_e_5";
02690 break;
02691 case 5:
02692 sig = "eps_e_6";
02693 break;
02694 default:
02695 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
02696 }
02697 break;
02698 case plbeam:
02699 switch (dir)
02700 {
02701 case 0:
02702 sig = "u";
02703 break;
02704 case 1:
02705 sig = "w";
02706 break;
02707 case 2:
02708 sig = "phi";
02709 break;
02710 default:
02711 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
02712 }
02713 break;
02714 case planestress:
02715 case planestrain:
02716 switch (dir)
02717 {
02718 case 0:
02719 sig = "eps_e_x";
02720 break;
02721 case 1:
02722 sig = "eps_e_y";
02723 break;
02724 case 2:
02725 sig = "eps_e_xy";
02726 break;
02727 case 3:
02728 sig = "eps_e_z";
02729 break;
02730 default:
02731 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
02732 }
02733 break;
02734
02735 case plates:{
02736 switch (dir){
02737 case 0:{
02738 sig = "kappa_x";
02739 break;
02740 }
02741 case 1:{
02742 sig = "kappa_y";
02743 break;
02744 }
02745 case 2:{
02746 sig = "kappa_xy";
02747 break;
02748 }
02749 case 3:{
02750 sig = "gamma_x";
02751 break;
02752 }
02753 case 4:{
02754 sig = "gamma_y";
02755 break;
02756 }
02757
02758 default:{
02759 print_err("unknown direction of strain is required", __FILE__, __LINE__, __func__);
02760 }
02761 }
02762 break;
02763 }
02764
02765 case platek:{
02766 switch (dir){
02767 case 0:{
02768 sig = "kappa_x";
02769 break;
02770 }
02771 case 1:{
02772 sig = "kappa_y";
02773 break;
02774 }
02775 case 2:{
02776 sig = "kappa_xy";
02777 break;
02778 }
02779 default:{
02780 print_err("unknown direction of strain is required", __FILE__, __LINE__, __func__);
02781 }
02782 }
02783 break;
02784 }
02785
02786 case axisymm:
02787 switch (dir)
02788 {
02789 case 0:
02790 sig = "eps_e_r";
02791 break;
02792 case 1:
02793 sig = "eps_e_z";
02794 break;
02795 case 2:
02796 sig = "eps_e_phi";
02797 break;
02798 case 3:
02799 sig = "eps_e_rz";
02800 break;
02801 default:
02802 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
02803 }
02804 break;
02805
02806 case shell:{
02807 switch (dir){
02808 case 0:
02809 sig = "eps_x";
02810 break;
02811 case 1:
02812 sig = "eps_y";
02813 break;
02814 case 2:
02815 sig = "gamma_xy";
02816 break;
02817 case 3:
02818 sig = "kappa_x";
02819 break;
02820 case 4:
02821 sig = "kappa_y";
02822 break;
02823 case 5:
02824 sig = "kappa_xy";
02825 break;
02826 default:
02827 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
02828 }
02829 break;
02830 }
02831
02832 case spacestress:
02833 switch (dir)
02834 {
02835 case 0:
02836 sig = "eps_e_x";
02837 break;
02838 case 1:
02839 sig = "eps_e_y";
02840 break;
02841 case 2:
02842 sig = "eps_e_z";
02843 break;
02844 case 3:
02845 sig = "eps_e_yz";
02846 break;
02847 case 4:
02848 sig = "eps_e_xz";
02849 break;
02850 case 5:
02851 sig = "eps_e_xy";
02852 break;
02853 default:
02854 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
02855 }
02856 break;
02857 default:
02858 print_err("unknown stress/strain state is required", __FILE__, __LINE__, __func__);
02859 }
02860 }
02861 switch (te)
02862 {
02863 case bar2d:
02864 case bar3d:
02865 sprintf(gpname, "Lin_1D");
02866 break;
02867 case beam2d:
02868 case beam3d:
02869 case subsoilbeam:
02870 sprintf(gpname, "Beam_1D");
02871 break;
02872 case barq2d:
02873 case barq3d:
02874 sprintf(gpname, "Quad_1D");
02875 break;
02876 default:
02877 sprintf(gpname, "%d", te);
02878 }
02879 if (scal == other)
02880 fprintf(out, "\nResult \"other_e_%ld\" \"%ld\" %s Scalar OnGaussPoints \"%s\"\n", dir+1, lcid, desclcid, gpname);
02881 else
02882 fprintf(out, "\nResult \"%s\" \"%ld\" %s Scalar OnGaussPoints \"%s\"\n", sig, lcid, desclcid, gpname);
02883 fprintf(out, "Values\n");
02884 print_header = 0;
02885 }
02886 switch (scal)
02887 {
02888 case strain:
02889 fprintf(out, "%7ld", i+Outdm->ide1);
02890 for (j = 0; j < tnip; j++)
02891 {
02892 ipp = Mt->elements[i].ipp[0][0]+j;
02893 ncompstr = Mm->ip[ipp].ncompstr;
02894 if (j == 0)
02895 {
02896
02897 if (te == subsoilplatetr)
02898 continue;
02899
02900 if (Outdm->eog.transtra[ir] > 0)
02901 {
02902 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
02903 Outdm->lcs[Outdm->eog.transtra[ir]-1].evaluate_t(p, namevar, t);
02904 str.n = ncompstr;
02905 str.a = Mm->ip[ipp].strain+ncompstr*lcid;
02906 gl_comp_engvectortransf(str, aux, dir, t, Mm->ip[ipp].ssst, strain);
02907 str.a = NULL;
02908 fprintf(out, " % e\n", aux);
02909 }
02910 else
02911 fprintf(out, " % e\n", Mm->ip[ipp].strain[ncompstr*lcid+dir]);
02912 }
02913 else
02914 {
02915 if (Outdm->eog.transtra[ir] > 0)
02916 {
02917 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
02918 Outdm->lcs[Outdm->eog.transtra[ir]-1].evaluate_t(p, namevar, t);
02919 str.n = ncompstr;
02920 str.a = Mm->ip[ipp].strain+ncompstr*lcid;
02921 gl_comp_engvectortransf(str, aux, dir, t, Mm->ip[ipp].ssst, strain);
02922 str.a = NULL;
02923 fprintf(out, "%7c % e\n", ' ', aux);
02924 }
02925 else
02926 fprintf(out, "%7c % e\n", ' ', Mm->ip[ipp].strain[ncompstr*lcid+dir]);
02927 }
02928 }
02929 ipp = Mt->elements[i].ipp[0][0];
02930 ncompstr = Mm->ip[ipp].ncompstr;
02931 if (Outdm->eog.transtra[ir] > 0)
02932 {
02933 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
02934 Outdm->lcs[Outdm->eog.transtra[ir]-1].evaluate_t(p, namevar, t);
02935 str.n = ncompstr;
02936 str.a = Mm->ip[ipp].strain+ncompstr*lcid;
02937 gl_comp_engvectortransf(str, aux, dir, t, Mm->ip[ipp].ssst, strain);
02938 str.a = NULL;
02939 }
02940 for (j = 0; j < num_copy_ip; j++)
02941 {
02942 if (Outdm->eog.transtra[ir] > 0)
02943 fprintf(out, "%7c % e\n", ' ', aux);
02944 else
02945 fprintf(out, "%7c % e\n", ' ', Mm->ip[ipp].strain[ncompstr*lcid+dir]);
02946 }
02947 fprintf(out, "\n");
02948 break;
02949 case stress:
02950 fprintf(out, "%7ld", i+Outdm->ide1);
02951 for (j = 0; j < tnip; j++)
02952 {
02953 ipp = Mt->elements[i].ipp[0][0]+j;
02954 ncompstr = Mm->ip[ipp].ncompstr;
02955 if (j == 0)
02956 {
02957
02958 if (te == subsoilplatetr)
02959 continue;
02960
02961 if (Outdm->eog.transtre[ir] > 0)
02962 {
02963 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
02964 Outdm->lcs[Outdm->eog.transtre[ir]-1].evaluate_t(p, namevar, t);
02965 str.n = ncompstr;
02966 str.a = Mm->ip[ipp].stress+ncompstr*lcid;
02967 gl_comp_engvectortransf(str, aux, dir, t, Mm->ip[ipp].ssst, stress);
02968 str.a = NULL;
02969 fprintf(out, " % e\n", aux);
02970 }
02971 else
02972 {
02973 fprintf(out, " % e\n", Mm->ip[ipp].stress[ncompstr*lcid+dir]);
02974 }
02975 }
02976 else
02977 {
02978 if (Outdm->eog.transtre[ir] > 0)
02979 {
02980 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
02981 Outdm->lcs[Outdm->eog.transtre[ir]-1].evaluate_t(p, namevar, t);
02982 str.n = ncompstr;
02983 str.a = Mm->ip[ipp].stress+ncompstr*lcid;
02984 gl_comp_engvectortransf(str, aux, dir, t, Mm->ip[ipp].ssst, stress);
02985 str.a = NULL;
02986 fprintf(out, "%7c % e\n", ' ', aux);
02987 }
02988 else
02989 fprintf(out, "%7c % e\n", ' ', Mm->ip[ipp].stress[ncompstr*lcid+dir]);
02990 }
02991 }
02992 ipp = Mt->elements[i].ipp[0][0];
02993 ncompstr = Mm->ip[ipp].ncompstr;
02994 if (Outdm->eog.transtre[ir] > 0)
02995 {
02996 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
02997 Outdm->lcs[Outdm->eog.transtre[ir]-1].evaluate_t(p, namevar, t);
02998 str.n = ncompstr;
02999 str.a = Mm->ip[ipp].stress+ncompstr*lcid;
03000 gl_comp_engvectortransf(str, aux, dir, t, Mm->ip[ipp].ssst, stress);
03001 str.a = NULL;
03002 }
03003 for (j = 0; j < num_copy_ip; j++)
03004 {
03005 if (Outdm->eog.transtre[ir] > 0)
03006 fprintf(out, "%7c % e\n", ' ', aux);
03007 else
03008 fprintf(out, "%7c % e\n", ' ', Mm->ip[ipp].stress[ncompstr*lcid+dir]);
03009 }
03010 fprintf(out, "\n");
03011 break;
03012 case other:
03013 fprintf(out, "%7ld", i+Outdm->ide1);
03014 for (j = 0; j < tnip; j++)
03015 {
03016 ipp = Mt->elements[i].ipp[0][0]+j;
03017 ncompother = Mm->ip[ipp].ncompeqother;
03018 if (j == 0)
03019 {
03020
03021 if (te == subsoilplatetr)
03022 continue;
03023
03024 fprintf(out, " % e\n", Mm->ip[ipp].eqother[ncompother*lcid+dir]);
03025 }
03026 else
03027 fprintf(out, "%7c % e\n", ' ', Mm->ip[ipp].eqother[ncompother*lcid+dir]);
03028 }
03029 for (j = 0; j < num_copy_ip; j++)
03030 {
03031 ipp = Mt->elements[i].ipp[0][0];
03032 ncompother = Mm->ip[ipp].ncompeqother;
03033 fprintf(out, "%7c % e\n", ' ', Mm->ip[ipp].eqother[ncompother*lcid+dir]);
03034 }
03035
03036
03037
03038
03039
03040
03041
03042
03043
03044
03045
03046
03047
03048
03049
03050 fprintf(out, "\n");
03051 break;
03052 default:
03053 print_err("unknown value type is required",__FILE__, __LINE__, __func__);
03054 }
03055 }
03056 if (print_header == 0)
03057 fprintf(out, "End Values\n");
03058 }
03059
03060
03061
03062
03063
03064
03065
03066
03067
03068
03069
03070
03071
03072
03073
03074
03075
03076
03077 void write_gid_elemvector(FILE *out, strastre q, long lcid, long sid, const char *desclcid)
03078 {
03079 if (Bar2d)
03080 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, bar2d);
03081
03082 if (Bar3d)
03083 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, bar3d);
03084
03085 if (Barq2d)
03086 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, barq2d);
03087
03088 if (Barq3d)
03089 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, barq3d);
03090
03091 if (Beam2d)
03092 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, beam2d);
03093
03094 if (Beam3d)
03095 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, beam3d);
03096
03097 if (Sbeam)
03098 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, subsoilbeam);
03099
03100 if (Spring)
03101 {
03102 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, spring_1);
03103 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, spring_2);
03104 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, spring_3);
03105 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, spring_4);
03106 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, spring_5);
03107 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, spring_6);
03108 }
03109
03110 if (Pqcon)
03111 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, planequadcontact);
03112
03113 if (Pelt)
03114 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, planeelementlt);
03115
03116 if (Peqt)
03117 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, planeelementqt);
03118
03119 if (Perlt)
03120 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, planeelementrotlt);
03121
03122 if (Pelq)
03123 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, planeelementlq);
03124
03125 if (Peqq)
03126 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, planeelementqq);
03127
03128 if (Perlq)
03129 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, planeelementrotlq);
03130
03131 if (Pesqt)
03132 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, planeelementsubqt);
03133
03134 if (Cct)
03135 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, cctel);
03136
03137 if (Dkt)
03138 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, dktel);
03139
03140 if (Dst)
03141 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, dstel);
03142
03143 if (Q4pl)
03144 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, q4plateel);
03145
03146 if (Spltr)
03147 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, subsoilplatetr);
03148
03149 if (Splq)
03150 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, subsoilplateq);
03151
03152 if (Asymlq)
03153 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, axisymmlq);
03154
03155 if (Asymlt)
03156 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, axisymmlt);
03157
03158 if (Asymqq)
03159 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, axisymmqq);
03160
03161 if (Shtr)
03162 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, shelltrelem);
03163
03164 if (Shq)
03165 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, shellqelem);
03166
03167 if (Ltet)
03168 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, lineartet);
03169
03170 if (Qtet)
03171 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, quadrtet);
03172
03173 if (Lhex)
03174 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, linearhex);
03175
03176 if (Qhex)
03177 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, quadrhex);
03178
03179 if (Lwed)
03180 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, linearwed);
03181
03182 if (Qwed)
03183 write_gid_elem_type_vector(out, q, lcid, sid, desclcid, quadrwed);
03184 }
03185
03186
03187
03188
03189
03190
03191
03192
03193
03194
03195
03196
03197
03198
03199
03200
03201
03202
03203
03204 void write_gid_elem_type_vector(FILE *out, strastre q, long lcid, long sid, const char *desclcid, elemtype te)
03205 {
03206 char sig[70];
03207 char gpname[1000];
03208 long i, j, k, ipp,ncompstr,ncompother,tnip;
03209 long print_header = 1;
03210 long asid, q_id1, q_n;
03211 long num_copy_ip;
03212
03213
03214 switch (q)
03215 {
03216 case strain:
03217 q_id1 = Outdm->eog.selstra[sid].id1[0];
03218 q_n = Outdm->eog.selstra[sid].ncomp[0];
03219 sprintf(sig, "eps_e_v%ld-%ld_s%ld", q_id1+1, q_n, sid+1);
03220 break;
03221 case stress:
03222 q_id1 = Outdm->eog.selstre[sid].id1[0];
03223 q_n = Outdm->eog.selstre[sid].ncomp[0];
03224 sprintf(sig, "sig_e_v%ld-%ld_s%ld", q_id1+1, q_n, sid+1);
03225 break;
03226 case other:
03227 q_id1 = Outdm->eog.seloth[sid].id1[0];
03228 q_n = Outdm->eog.seloth[sid].ncomp[0];
03229 sprintf(sig, "other_e_v%ld-%ld_s%ld", q_id1+1, q_n, sid+1);
03230 break;
03231 default:
03232 print_err("unsupported type of quantity (%d) is required in function\n", __FILE__, __LINE__, __func__, q);
03233 abort();
03234 }
03235
03236 for (i = 0; i < Mt->ne; i++)
03237 {
03238 if (Gtm->leso[i]==0)
03239 continue;
03240
03241
03242 if (Mt->elements[i].te != te)
03243 continue;
03244
03245
03246 if (q == strain)
03247 {
03248 if (Outdm->eog.selestra.presence_id(i, sid, asid) == 0)
03249 continue;
03250 }
03251 if (q == stress)
03252 {
03253 if (Outdm->eog.selestre.presence_id(i, sid, asid) == 0)
03254 continue;
03255 }
03256 if (q == other)
03257 {
03258 if (Outdm->eog.seleoth.presence_id(i, sid, asid) == 0)
03259 continue;
03260 }
03261
03262 tnip = Mt->give_totnip(i);
03263 if (te == planeelementlq)
03264 tnip = 4;
03265 num_copy_ip = 0;
03266 switch (te)
03267 {
03268 case beam3d:
03269 case subsoilbeam:
03270 num_copy_ip = 1;
03271 break;
03272 default:
03273 break;
03274 }
03275
03276 if (print_header)
03277 {
03278 ipp = Mt->elements[i].ipp[0][0];
03279 switch (te)
03280 {
03281 case bar2d:
03282 case bar3d:
03283 case beam2d:
03284 case beam3d:
03285 case subsoilbeam:
03286 sprintf(gpname, "Lin_1D");
03287 break;
03288 case barq2d:
03289 case barq3d:
03290 sprintf(gpname, "Quad_1D");
03291 break;
03292 default:
03293 sprintf(gpname, "%d", te);
03294 }
03295 fprintf(out, "\nResult \"%s\" \"%ld\" %s Vector OnGaussPoints \"%s\"\n", sig, lcid, desclcid, gpname);
03296 fprintf(out, "Values\n");
03297 print_header = 0;
03298 }
03299 switch (q)
03300 {
03301 case stress:
03302 for (j = 0; j < tnip; j++)
03303 {
03304 ipp = Mt->elements[i].ipp[0][0]+j;
03305 ncompstr = Mm->ip[ipp].ncompstr;
03306 if (j == 0)
03307 fprintf(out, "%7ld", i+Outdm->ide1);
03308 else
03309 fprintf(out, "%7c", ' ');
03310 for (k = q_id1; k < q_id1+q_n; k++)
03311 fprintf(out, " % e", Mm->ip[ipp].stress[ncompstr*lcid+k]);
03312 for (k = 0; k < 3-q_n; k++)
03313 fprintf(out, " 0.0");
03314 fprintf(out, "\n");
03315 }
03316 for (j = 0; j < num_copy_ip; j++)
03317 {
03318 ipp = Mt->elements[i].ipp[0][0];
03319 ncompstr = Mm->ip[ipp].ncompstr;
03320 fprintf(out, "%7c", ' ');
03321 for (k = q_id1; k < q_id1+q_n; k++)
03322 fprintf(out, " % e", Mm->ip[ipp].stress[ncompstr*lcid+k]);
03323 for (k = 0; k < 3-q_n; k++)
03324 fprintf(out, " 0.0");
03325 fprintf(out, "\n");
03326 }
03327 fprintf(out, "\n");
03328 break;
03329 case strain:
03330 for (j = 0; j < tnip; j++)
03331 {
03332 ipp = Mt->elements[i].ipp[0][0]+j;
03333 ncompstr = Mm->ip[ipp].ncompstr;
03334 if (j == 0)
03335 fprintf(out, "%7ld", i+Outdm->ide1);
03336 else
03337 fprintf(out, "%7c", ' ');
03338 for (k = q_id1; k < q_id1+q_n; k++)
03339 fprintf(out, " % e", Mm->ip[ipp].strain[ncompstr*lcid+k]);
03340 for (k = 0; k < 3-q_n; k++)
03341 fprintf(out, " 0.0");
03342 fprintf(out, "\n");
03343 }
03344 for (j = 0; j < num_copy_ip; j++)
03345 {
03346 ipp = Mt->elements[i].ipp[0][0];
03347 ncompstr = Mm->ip[ipp].ncompstr;
03348 fprintf(out, "%7c", ' ');
03349 for (k = q_id1; k < q_id1+q_n; k++)
03350 fprintf(out, " % e", Mm->ip[ipp].strain[ncompstr*lcid+k]);
03351 for (k = 0; k < 3-q_n; k++)
03352 fprintf(out, " 0.0");
03353 fprintf(out, "\n");
03354 }
03355 fprintf(out, "\n");
03356 break;
03357 case other:
03358 for (j = 0; j < tnip; j++)
03359 {
03360 ipp = Mt->elements[i].ipp[0][0]+j;
03361 ncompother = Mm->ip[ipp].ncompeqother;
03362 if (j == 0)
03363 fprintf(out, "%7ld", i+Outdm->ide1);
03364 else
03365 fprintf(out, "%7c", ' ');
03366 for (k = q_id1; k < q_id1+q_n; k++)
03367 fprintf(out, " % e", Mm->ip[ipp].eqother[ncompother*lcid+k]);
03368 for (k = 0; k < 3-q_n; k++)
03369 fprintf(out, " 0.0");
03370 fprintf(out, "\n");
03371 }
03372 for (j = 0; j < num_copy_ip; j++)
03373 {
03374 ipp = Mt->elements[i].ipp[0][0];
03375 ncompother = Mm->ip[ipp].ncompeqother;
03376 fprintf(out, "%7c", ' ');
03377 for (k = q_id1; k < q_id1+q_n; k++)
03378 fprintf(out, " % e", Mm->ip[ipp].eqother[ncompother*lcid+k]);
03379 for (k = 0; k < 3-q_n; k++)
03380 fprintf(out, " 0.0");
03381 fprintf(out, "\n");
03382 }
03383 fprintf(out, "\n");
03384 break;
03385 default:
03386 print_err("unknown value type is rquired", __FILE__, __LINE__, __func__);
03387 }
03388 }
03389 if (print_header == 0)
03390 fprintf(out, "End Values\n");
03391 }
03392
03393
03394
03395
03396
03397
03398
03399
03400
03401
03402
03403
03404
03405
03406
03407
03408
03409
03410 void write_gid_elemtensor(FILE *out, strastre q, long lcid, long sid, const char *desclcid)
03411 {
03412 if (Bar2d)
03413 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, bar2d);
03414
03415 if (Bar3d)
03416 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, bar3d);
03417
03418 if (Barq2d)
03419 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, barq2d);
03420
03421 if (Barq3d)
03422 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, barq3d);
03423
03424 if (Beam2d)
03425 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, beam2d);
03426
03427 if (Beam3d)
03428 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, beam3d);
03429
03430 if (Sbeam)
03431 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, subsoilbeam);
03432
03433 if (Spring)
03434 {
03435 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, spring_1);
03436 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, spring_2);
03437 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, spring_3);
03438 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, spring_4);
03439 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, spring_5);
03440 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, spring_6);
03441 }
03442
03443 if (Pqcon)
03444 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, planequadcontact);
03445
03446 if (Pelt)
03447 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, planeelementlt);
03448
03449 if (Peqt)
03450 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, planeelementqt);
03451
03452 if (Perlt)
03453 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, planeelementrotlt);
03454
03455 if (Pelq)
03456 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, planeelementlq);
03457
03458 if (Peqq)
03459 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, planeelementqq);
03460
03461 if (Perlq)
03462 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, planeelementrotlq);
03463
03464 if (Pesqt)
03465 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, planeelementsubqt);
03466
03467 if (Cct)
03468 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, cctel);
03469
03470 if (Dkt)
03471 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, dktel);
03472
03473 if (Dst)
03474 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, dstel);
03475
03476 if (Q4pl)
03477 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, q4plateel);
03478
03479 if (Spltr)
03480 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, subsoilplatetr);
03481
03482 if (Splq)
03483 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, subsoilplateq);
03484
03485 if (Asymlq)
03486 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, axisymmlq);
03487
03488 if (Asymlt)
03489 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, axisymmlt);
03490
03491 if (Asymqq)
03492 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, axisymmqq);
03493
03494 if (Shtr)
03495 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, shelltrelem);
03496
03497 if (Shq)
03498 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, shellqelem);
03499
03500 if (Ltet)
03501 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, lineartet);
03502
03503 if (Qtet)
03504 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, quadrtet);
03505
03506 if (Lhex)
03507 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, linearhex);
03508
03509 if (Qhex)
03510 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, quadrhex);
03511
03512 if (Lwed)
03513 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, linearwed);
03514
03515 if (Qwed)
03516 write_gid_elem_type_tensor(out, q, lcid, sid, desclcid, quadrwed);
03517 }
03518
03519
03520
03521
03522
03523
03524
03525
03526
03527
03528
03529
03530
03531
03532
03533
03534
03535
03536
03537 void write_gid_elem_type_tensor(FILE *out, strastre q, long lcid, long sid, const char *desclcid, elemtype te)
03538 {
03539 char sig[70];
03540 char gpname[1000];
03541 long i, j, ipp, ncompstr, ncompother, tnip;
03542 long print_header = 1;
03543 long asid, q_id1, q_n;
03544 long num_copy_ip;
03545 vector v;
03546 matrix t(3,3);
03547 matrix tmat;
03548 vector str, p;
03549 char *namevar[] = {"x", "y", "z"};
03550
03551 if(Outdm->nlcs > 0)
03552 {
03553 allocv(3, p);
03554 allocm(3, 3, tmat);
03555 }
03556
03557
03558
03559 switch (q)
03560 {
03561 case strain:
03562 q_id1 = 0;
03563 q_n = 0;
03564 sprintf(sig, "eps_e_m_s%ld", sid+1);
03565 break;
03566 case stress:
03567 q_id1 = 0;
03568 q_n = 0;
03569 sprintf(sig, "sig_e_m_s%ld", sid+1);
03570 break;
03571 case other:
03572 q_id1 = Outdm->eog.seloth[sid].id1[0];
03573 q_n = Outdm->eog.seloth[sid].ncomp[0];
03574 sprintf(sig, "other_e_m%ld-%ld_s%ld", q_id1+1, q_n, sid+1);
03575 break;
03576 default:
03577 print_err("unsupported type of quantity (%d) is required", __FILE__, __LINE__, __func__, q);
03578 abort();
03579 }
03580
03581
03582 if (q != other)
03583 {
03584 for (i = 0; i < Mt->ne; i++)
03585 {
03586 if ((Gtm->leso[i]==0) || (Mt->elements[i].te != te))
03587 continue;
03588 ipp = Mt->elements[i].ipp[0][0];
03589 switch (q)
03590 {
03591 case strain:
03592 if (Outdm->eog.selestra.presence_id(i, sid, asid) == 0)
03593 continue;
03594 if (q_n < Mm->ip[ipp].ncompstr)
03595 q_n = Mm->ip[ipp].ncompstr;
03596 break;
03597 case stress:
03598 if (Outdm->eog.selestre.presence_id(i, sid, asid) == 0)
03599 continue;
03600 if (q_n < Mm->ip[ipp].ncompstr)
03601 q_n = Mm->ip[ipp].ncompstr;
03602 break;
03603 default:
03604 break;
03605 }
03606 }
03607 if (q_n == 0)
03608 return;
03609 }
03610 if ((q_n != 6) && (q_n != 4) && (q_n != 2) && (q_n != 1))
03611 {
03612 print_err("invalid number of required tensor components (%ld),\n"
03613 " only 1, 2, 4 or 6 components are allowed.", __FILE__, __LINE__, __func__, q_n);
03614 abort();
03615 }
03616
03617 for (i = 0; i < Mt->ne; i++)
03618 {
03619
03620 if (Gtm->leso[i]==0)
03621 continue;
03622
03623
03624 if (Mt->elements[i].te != te)
03625 continue;
03626
03627
03628 if ((q == strain) && (Outdm->eog.selestra.presence_id(i, sid, asid) == 0))
03629 continue;
03630 if ((q == stress) && (Outdm->eog.selestre.presence_id(i, sid, asid) == 0))
03631 continue;
03632 if ((q == other) && (Outdm->eog.seleoth.presence_id(i, sid, asid) == 0))
03633 continue;
03634
03635 tnip = Mt->give_totnip(i);
03636 if (te == planeelementlq)
03637 tnip = 4;
03638 num_copy_ip = 0;
03639 switch (te)
03640 {
03641 case beam3d:
03642 case subsoilbeam:
03643 num_copy_ip = 1;
03644 break;
03645 default:
03646 break;
03647 }
03648
03649 if (print_header)
03650 {
03651 ipp = Mt->elements[i].ipp[0][0];
03652 switch (te)
03653 {
03654 case bar2d:
03655 case bar3d:
03656 case beam2d:
03657 case beam3d:
03658 case subsoilbeam:
03659 sprintf(gpname, "Lin_1D");
03660 break;
03661 case barq2d:
03662 case barq3d:
03663 sprintf(gpname, "Quad_1D");
03664 break;
03665 default:
03666 sprintf(gpname, "%d", te);
03667 }
03668 if (q_n == 6)
03669 fprintf(out, "\nResult \"%s\" \"%ld\" %s Matrix OnGaussPoints \"%s\"\n", sig, lcid, desclcid, gpname);
03670 else
03671 fprintf(out, "\nResult \"%s\" \"%ld\" %s PlainDeformationMatrix OnGaussPoints \"%s\"\n", sig, lcid, desclcid, gpname);
03672 fprintf(out, "Values\n");
03673 print_header = 0;
03674 }
03675 switch (q)
03676 {
03677 case strain:
03678 for (j = 0; j < tnip; j++)
03679 {
03680 ipp = Mt->elements[i].ipp[0][0]+j;
03681 ncompstr = Mm->ip[ipp].ncompstr;
03682 if (j == 0)
03683 fprintf(out, "%7ld", i+Outdm->ide1);
03684 else
03685 fprintf(out, "%7c", ' ');
03686 allocv(ncompstr, v);
03687 copyv(Mm->ip[ipp].strain+ncompstr*lcid+q_id1, v);
03688 vector_tensor (v, t, Mm->ip[ipp].ssst, strain);
03689 if (Outdm->eog.transtra[sid] > 0)
03690 {
03691 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
03692 Outdm->lcs[Outdm->eog.transtra[sid]-1].evaluate_t(p, namevar, tmat);
03693 glmatrixtransf(t, tmat);
03694 }
03695 if (q_n == 6)
03696 fprintf(out, " % e % e % e % e % e % e\n", t[0][0], t[1][1], t[2][2], t[0][1], t[1][2], t[0][2]);
03697 else
03698 fprintf(out, " % e % e % e % e\n", t[0][0], t[1][1], t[0][1], t[2][2]);
03699 destrv(v);
03700 }
03701 for (j = 0; j < num_copy_ip; j++)
03702 {
03703 ipp = Mt->elements[i].ipp[0][0];
03704 ncompstr = Mm->ip[ipp].ncompstr;
03705 fprintf(out, "%7c", ' ');
03706 allocv(ncompstr, v);
03707 copyv(Mm->ip[ipp].strain+ncompstr*lcid+q_id1, v);
03708 vector_tensor (v, t, Mm->ip[ipp].ssst, strain);
03709 if (Outdm->eog.transtra[sid] > 0)
03710 {
03711 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
03712 Outdm->lcs[Outdm->eog.transtra[sid]-1].evaluate_t(p, namevar, tmat);
03713 glmatrixtransf(t, tmat);
03714 }
03715 if (q_n == 6)
03716 fprintf(out, " % e % e % e % e % e % e\n", t[0][0], t[1][1], t[2][2], t[0][1], t[1][2], t[0][2]);
03717 else
03718 fprintf(out, " % e % e % e % e\n", t[0][0], t[1][1], t[0][1], t[2][2]);
03719 destrv(v);
03720 }
03721 fprintf(out, "\n");
03722 break;
03723 case stress:
03724 for (j = 0; j < tnip; j++)
03725 {
03726 ipp = Mt->elements[i].ipp[0][0]+j;
03727 ncompstr = Mm->ip[ipp].ncompstr;
03728 if (j == 0)
03729 fprintf(out, "%7ld", i+Outdm->ide1);
03730 else
03731 fprintf(out, "%7c", ' ');
03732 allocv(ncompstr, v);
03733 copyv(Mm->ip[ipp].stress+ncompstr*lcid+q_id1, v);
03734 vector_tensor (v, t, Mm->ip[ipp].ssst, stress);
03735 if (Outdm->eog.transtre[sid] > 0)
03736 {
03737 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
03738 Outdm->lcs[Outdm->eog.transtre[sid]-1].evaluate_t(p, namevar, tmat);
03739 glmatrixtransf(t, tmat);
03740 }
03741 if (q_n == 6)
03742 fprintf(out, " % e % e % e % e % e % e\n", t[0][0], t[1][1], t[2][2], t[0][1], t[1][2], t[0][2]);
03743 else
03744 fprintf(out, " % e % e % e % e\n", t[0][0], t[1][1], t[0][1], t[2][2]);
03745 destrv(v);
03746 }
03747 for (j = 0; j < num_copy_ip; j++)
03748 {
03749 ipp = Mt->elements[i].ipp[0][0];
03750 ncompstr = Mm->ip[ipp].ncompstr;
03751 fprintf(out, "%7c", ' ');
03752 allocv(ncompstr, v);
03753 copyv(Mm->ip[ipp].stress+ncompstr*lcid+q_id1, v);
03754 vector_tensor (v, t, Mm->ip[ipp].ssst, stress);
03755 if (Outdm->eog.transtre[sid] > 0)
03756 {
03757 ipcoord(Mm->elip[ipp], ipp, 0, 0, p);
03758 Outdm->lcs[Outdm->eog.transtre[sid]-1].evaluate_t(p, namevar, tmat);
03759 glmatrixtransf(t, tmat);
03760 }
03761 if (q_n == 6)
03762 fprintf(out, " % e % e % e % e % e % e\n", t[0][0], t[1][1], t[2][2], t[0][1], t[1][2], t[0][2]);
03763 else
03764 fprintf(out, " % e % e % e % e\n", t[0][0], t[1][1], t[0][1], t[2][2]);
03765 destrv(v);
03766 }
03767 fprintf(out, "\n");
03768 break;
03769 case other:
03770 for (j = 0; j < tnip; j++)
03771 {
03772 ipp = Mt->elements[i].ipp[0][0]+j;
03773 ncompother = Mm->ip[ipp].ncompeqother;
03774 if (j == 0)
03775 fprintf(out, "%7ld", i+Outdm->ide1);
03776 else
03777 fprintf(out, "%7c", ' ');
03778 allocv(q_n, v);
03779 copyv(Mm->ip[ipp].eqother+ncompother*lcid+q_id1, v);
03780 vector_tensor (v, t, guess_ssst(q_n), stress);
03781 if (q_n == 6)
03782 fprintf(out, " % e % e % e % e % e % e\n", t[0][0], t[1][1], t[2][2], t[0][1], t[1][2], t[0][2]);
03783 else
03784 fprintf(out, " % e % e % e % e\n", t[0][0], t[1][1], t[0][1], t[2][2]);
03785 destrv(v);
03786 }
03787 for (j = 0; j < num_copy_ip; j++)
03788 {
03789 ipp = Mt->elements[i].ipp[0][0];
03790 ncompother = Mm->ip[ipp].ncompeqother;
03791 fprintf(out, "%7c", ' ');
03792 allocv(q_n, v);
03793 copyv(Mm->ip[ipp].eqother+ncompother*lcid+q_id1, v);
03794 vector_tensor (v, t, guess_ssst(q_n), stress);
03795 if (q_n == 6)
03796 fprintf(out, " % e % e % e % e % e % e\n", t[0][0], t[1][1], t[2][2], t[0][1], t[1][2], t[0][2]);
03797 else
03798 fprintf(out, " % e % e % e % e\n", t[0][0], t[1][1], t[0][1], t[2][2]);
03799 destrv(v);
03800 }
03801 fprintf(out, "\n");
03802 break;
03803 default:
03804 print_err("unknown value type is required", __FILE__, __LINE__, __func__);
03805 }
03806 }
03807 if (print_header == 0)
03808 fprintf(out, "End Values\n");
03809 }
03810
03811
03812
03813
03814
03815
03816
03817
03818
03819
03820
03821
03822 void export_femcad(FILE *out)
03823 {
03824 write_nodes(out);
03825 write_elements(out);
03826 }
03827
03828
03829
03830
03831
03832
03833
03834
03835
03836
03837
03838
03839 void write_nodes(FILE *out)
03840 {
03841 long i;
03842 fprintf(out, "\nNODES\n");
03843 for (i=0; i < Gtm->nn; i++)
03844 fprintf(out, "%ld %e %e %e\n", i, Gtm->gnodes[i].x, Gtm->gnodes[i].y, Gtm->gnodes[i].z);
03845 }
03846
03847
03848
03849
03850
03851
03852
03853
03854
03855
03856
03857
03858 void write_elements(FILE *out)
03859 {
03860 long i,j;
03861
03862 fprintf(out, "\nELEMENTS\n");
03863 for (i = 0; i < Mt->ne; i++)
03864 {
03865 if (Gtm->leso[i]==0)
03866 continue;
03867
03868 switch (Mt->elements[i].te)
03869 {
03870 case bar2d:
03871 case bar3d:
03872 case beam2d:
03873 case beam3d:
03874 case subsoilbeam:
03875
03876
03877
03878 break;
03879 case barq2d:
03880 case barq3d:
03881
03882
03883
03884 break;
03885 case spring_1:
03886 case spring_2:
03887 case spring_3:
03888 case spring_4:
03889 case spring_5:
03890 case spring_6:
03891 break;
03892 case cctel:
03893 case dktel:
03894 case dstel:
03895 case subsoilplatetr:
03896 case planeelementlt:
03897 case planeelementrotlt:
03898 fprintf(out, "%ld %d", i, TriangleLinear);
03899 for (j = 0; j < Gtm->gelements[i].nne; j++)
03900 fprintf(out, " %ld", Gtm->gelements[i].nodes[j]);
03901 break;
03902 case planeelementqt:
03903 fprintf(out, "%ld %d", i, TriangleQuadratic);
03904 for (j = 0; j < Gtm->gelements[i].nne; j++)
03905 fprintf(out, " %ld", Gtm->gelements[i].nodes[j]);
03906 break;
03907 case planeelementlq:
03908 case planeelementrotlq:
03909 case q4plateel:
03910 case subsoilplateq:
03911 fprintf(out, "%ld %d", i, ISOLinear2D);
03912 for (j = 0; j < Gtm->gelements[i].nne; j++)
03913 fprintf(out, " %ld", Gtm->gelements[i].nodes[j]);
03914 break;
03915 case planeelementqq:
03916 fprintf(out, "%ld %d", i, ISOQuadratic2D);
03917 for (j = 0; j < Gtm->gelements[i].nne; j++)
03918 fprintf(out, " %ld", Gtm->gelements[i].nodes[j]);
03919 break;
03920
03921
03922 case linearhex:
03923 fprintf(out, "%ld %d", i, ISOLinear2D);
03924 for (j = 0; j < 4; j++)
03925 fprintf(out, " %ld", Gtm->gelements[i].nodes[4+j]);
03926 break;
03927 default:
03928 print_err("unsupported element type is required",__FILE__, __LINE__, __func__);
03929 }
03930 fprintf(out, "\n");
03931 }
03932 }
03933
03934
03935
03936
03937
03938
03939
03940
03941
03942
03943
03944
03945
03946
03947 void write_nodes_prep(FILE *out, siftop *st)
03948 {
03949 long i;
03950 fprintf(out, "\nNODES\n");
03951 for (i=0; i < st->nn; i++)
03952 fprintf(out, "%ld %e %e %e\n", i, st->nodes[i].x, st->nodes[i].y, st->nodes[i].z);
03953 }
03954
03955
03956
03957
03958
03959
03960
03961
03962
03963
03964
03965
03966
03967
03968 void write_elements_prep(FILE *out, siftop *st)
03969 {
03970 long i,j;
03971
03972 fprintf(out, "\nELEMENTS\n");
03973 for (i = 0; i < st->ne; i++)
03974 {
03975 fprintf(out, "%ld %d", i, TriangleLinear);
03976 for (j = 0; j < st->elements[i].nne; j++)
03977 fprintf(out, " %ld", st->elements[i].nodes[j]-1);
03978 fprintf(out, "\n");
03979 }
03980 }
03981
03982
03983
03984
03985
03986
03987
03988
03989
03990
03991
03992
03993
03994
03995
03996 void write_displ(FILE *out, long lcid, const char *desclcid)
03997 {
03998 long i, j;
03999 double *r;
04000 fprintf (out,"\nVALUES");
04001 fprintf (out,"\nVECTOR3D NODES DEFORMATION %s\n", desclcid);
04002 for (i=0; i < Mt->nn; i++)
04003 {
04004 r = new double [Gtm->gnodes[i].ndofn];
04005 noddispl (lcid, r, i);
04006 fprintf(out, "%ld", i);
04007 for (j = 0; j < Gtm->gnodes[i].ndofn; j++)
04008 fprintf(out, " % e", r[j]);
04009 fprintf(out, " 0.0 \n");
04010 delete [] r;
04011 }
04012 }
04013
04014
04015
04016
04017
04018
04019
04020
04021
04022
04023
04024
04025
04026
04027
04028
04029 void write_nforces(FILE *out, long lcid, const char *desclcid, double *ifor)
04030 {
04031 long i, j, ii, ndof;
04032 vector f, g;
04033 matrix tm(3,3);
04034 long ndofn, transf;
04035
04036 fprintf (out,"\nVALUES");
04037 fprintf (out,"\nVECTOR3D NODES FORCES %s\n", desclcid);
04038 for (i=0; i < Mt->nn; i++)
04039 {
04040 ndof = Mt->give_ndofn(i);
04041 allocv(ndof, f);
04042 for (j=0;j<ndof;j++)
04043 {
04044 ii=Mt->give_dof(i,j);
04045 if (ii<1) f[j]=Mt->nodes[i].r[j];
04046 if (ii>0) f[j]=ifor[ii-1];
04047 }
04048
04049 transf = Mt->nodes[i].transf;
04050 if (transf>0){
04051
04052
04053
04054 ndofn = Gtm->gnodes[i].ndofn;
04055 allocv(ndofn, g);
04056 allocm(ndofn,ndofn,tm);
04057 fillm(0.0, tm);
04058
04059 tm[0][0]=Mt->nodes[i].e1[0];
04060 tm[1][0]=Mt->nodes[i].e1[1];
04061 if (transf == 3)
04062 tm[2][0]=Mt->nodes[i].e1[2];
04063
04064 tm[0][1]=Mt->nodes[i].e2[0];
04065 tm[1][1]=Mt->nodes[i].e2[1];
04066 if (transf == 3)
04067 tm[2][1]=Mt->nodes[i].e2[2];
04068
04069 if (transf == 3)
04070 {
04071 tm[0][2]=Mt->nodes[i].e3[0];
04072 tm[1][2]=Mt->nodes[i].e3[1];
04073 tm[2][2]=Mt->nodes[i].e3[2];
04074 }
04075 else
04076 {
04077 if (ndofn == 3)
04078 tm[2][2]=1.0;
04079 }
04080
04081 mxv(tm,f,g);
04082 copyv (g,f);
04083 destrv (g);
04084 destrm (tm);
04085 }
04086 fprintf(out, "%ld", i);
04087 if (ndof > 3) ndof = 3;
04088 for (j = 0; j < ndof; j++)
04089 fprintf(out, " % e", f[j]);
04090 for (j=ndof; j < 3; j++)
04091 fprintf(out, " 0.0");
04092 fprintf(out, "\n");
04093 destrv(f);
04094 }
04095 }
04096
04097
04098
04099
04100
04101
04102
04103
04104
04105
04106
04107
04108
04109
04110
04111
04112
04113 void write_nodscalar(FILE *out, strastre scal, long lcid, long dir, const char *desclcid)
04114 {
04115 const char *sig = "";
04116 long i, ncompstr, ncompother;
04117
04118 fprintf(out, "\nVALUES\nSCALAR NODES ");
04119 if (scal == stress)
04120 {
04121 if (Mt->nodes[0].ncompstr < 6)
04122 {
04123 switch (dir)
04124 {
04125 case 0:
04126 sig = "sig_x";
04127 break;
04128 case 1:
04129 sig = "sig_y";
04130 break;
04131 case 2:
04132 sig = "tau_xy";
04133 break;
04134 case 3:
04135 sig = "sig_z";
04136 break;
04137 default:
04138 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
04139 }
04140 }
04141 if (Mt->nodes[0].ncompstr == 6)
04142 {
04143 switch (dir)
04144 {
04145 case 0:
04146 sig = "sig_x";
04147 break;
04148 case 1:
04149 sig = "sig_y";
04150 break;
04151 case 2:
04152 sig = "sig_z";
04153 break;
04154 case 3:
04155 sig = "tau_yz";
04156 break;
04157 case 4:
04158 sig = "tau_xz";
04159 break;
04160 case 5:
04161 sig = "tau_xy";
04162 break;
04163 default:
04164 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
04165 }
04166 }
04167 }
04168 if (scal == strain)
04169 {
04170 if (Mt->nodes[0].ncompstr < 6)
04171 {
04172 switch (dir)
04173 {
04174 case 0:
04175 sig = "eps_x";
04176 break;
04177 case 1:
04178 sig = "eps_y";
04179 break;
04180 case 2:
04181 sig = "eps_xy";
04182 break;
04183 case 3:
04184 sig = "eps_z";
04185 break;
04186 default:
04187 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
04188 }
04189 }
04190 if (Mt->nodes[0].ncompstr == 6)
04191 {
04192 switch (dir)
04193 {
04194 case 0:
04195 sig = "eps_x";
04196 break;
04197 case 1:
04198 sig = "eps_y";
04199 break;
04200 case 2:
04201 sig = "eps_z";
04202 break;
04203 case 3:
04204 sig = "eps_yz";
04205 break;
04206 case 4:
04207 sig = "eps_xz";
04208 break;
04209 case 5:
04210 sig = "eps_xy";
04211 break;
04212 default:
04213 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
04214 }
04215 }
04216 }
04217 switch (scal)
04218 {
04219 case stress:
04220 fprintf(out, "%s %s\n", sig, desclcid);
04221 for (i = 0; i < Mt->nn; i++)
04222 {
04223 ncompstr = Mt->nodes[i].ncompstr;
04224 fprintf(out, "%ld % e\n", i, Mt->nodes[i].stress[ncompstr*lcid+dir]);
04225 }
04226 break;
04227 case strain:
04228 fprintf(out, "%s %s\n", sig, desclcid);
04229 for (i = 0; i < Mt->nn; i++)
04230 {
04231 ncompstr = Mt->nodes[i].ncompstr;
04232 fprintf(out, "%ld % e\n", i, Mt->nodes[i].strain[ncompstr*lcid+dir]);
04233 }
04234 break;
04235 case other:
04236 fprintf(out, "other_n%ld %s\n", dir, desclcid);
04237 for (i = 0; i < Mt->nn; i++)
04238 {
04239 ncompother = Mt->nodes[i].ncompother;
04240 fprintf(out, "%ld % e\n", i, Mt->nodes[i].other[ncompother*lcid+dir]);
04241 }
04242 break;
04243 default:
04244 print_err("unknown value type is required", __FILE__, __LINE__, __func__);
04245 }
04246 }
04247
04248
04249
04250
04251
04252
04253
04254
04255
04256
04257
04258
04259
04260
04261
04262
04263
04264 void write_elemscalar(FILE *out, strastre scal, long lcid, long dir, const char * desclcid)
04265 {
04266 const char *sig = "";
04267 long i, ipp,ncompstr,ncompother;
04268
04269 fprintf(out, "\nVALUES\nSCALAR ELEMENTS ");
04270 if (scal == stress)
04271 {
04272 if (Mm->ip[0].ncompstr < 6)
04273 {
04274 switch (dir)
04275 {
04276 case 0:
04277 sig = "sig_x";
04278 break;
04279 case 1:
04280 sig = "sig_y";
04281 break;
04282 case 2:
04283 sig = "tau_xy";
04284 break;
04285 case 3:
04286 sig = "sig_z";
04287 break;
04288 default:
04289 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
04290 }
04291 }
04292 if (Mm->ip[0].ncompstr == 6)
04293 {
04294 switch (dir)
04295 {
04296 case 0:
04297 sig = "sig_x";
04298 break;
04299 case 1:
04300 sig = "sig_y";
04301 break;
04302 case 2:
04303 sig = "sig_z";
04304 break;
04305 case 3:
04306 sig = "tau_yz";
04307 break;
04308 case 4:
04309 sig = "tau_xz";
04310 break;
04311 case 5:
04312 sig = "tau_xy";
04313 break;
04314 default:
04315 print_err("unknown direction of stresses is required", __FILE__, __LINE__, __func__);
04316 }
04317 }
04318 }
04319 if (scal == strain)
04320 {
04321 if (Mm->ip[0].ncompstr < 6)
04322 {
04323 switch (dir)
04324 {
04325 case 0:
04326 sig = "eps_x";
04327 break;
04328 case 1:
04329 sig = "eps_y";
04330 break;
04331 case 2:
04332 sig = "eps_xy";
04333 break;
04334 case 3:
04335 sig = "eps_z";
04336 break;
04337 default:
04338 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
04339 }
04340 }
04341 if (Mm->ip[0].ncompstr == 6)
04342 {
04343 switch (dir)
04344 {
04345 case 0:
04346 sig = "eps_x";
04347 break;
04348 case 1:
04349 sig = "eps_y";
04350 break;
04351 case 2:
04352 sig = "eps_z";
04353 break;
04354 case 3:
04355 sig = "eps_yz";
04356 break;
04357 case 4:
04358 sig = "eps_xz";
04359 break;
04360 case 5:
04361 sig = "eps_xy";
04362 break;
04363 default:
04364 print_err("unknown direction of strains is required", __FILE__, __LINE__, __func__);
04365 }
04366 }
04367 }
04368 switch (scal)
04369 {
04370 case stress:
04371 fprintf(out, "%s %s\n", sig, desclcid);
04372 for (i = 0; i < Mt->ne; i++)
04373 {
04374 if (Gtm->leso[i]==0)
04375 continue;
04376
04377 ipp = Mt->elements[i].ipp[0][0];
04378 ncompstr = Mm->ip[ipp].ncompstr;
04379 fprintf(out, "%ld % e\n", i, Mm->ip[ipp].stress[ncompstr*lcid+dir]);
04380 }
04381 break;
04382 case strain:
04383 fprintf(out, "%s %s\n", sig, desclcid);
04384 for (i = 0; i < Mt->nn; i++)
04385 {
04386 ipp = Mt->elements[i].ipp[0][0];
04387 ncompstr = Mm->ip[ipp].ncompstr;
04388 fprintf(out, "%ld % e\n", i, Mm->ip[ipp].strain[ncompstr*lcid+dir]);
04389 }
04390 break;
04391 case other:
04392 fprintf(out, "other_e%ld %s\n", dir, desclcid);
04393 for (i = 0; i < Mt->nn; i++)
04394 {
04395 ipp = Mt->elements[i].ipp[0][0];
04396 ncompother = Mm->ip[ipp].ncompother;
04397 fprintf(out, "%ld % e\n", i, Mm->ip[ipp].eqother[ncompother*lcid+dir]);
04398 }
04399 break;
04400 default:
04401 print_err("unknown value type is required", __FILE__, __LINE__, __func__);
04402 }
04403 }
04404
04405
04406
04407
04408
04409
04410
04411
04412
04413
04414
04415
04416
04417
04418
04419
04420 void write_nodscalar(FILE *out,double *val, const char *descr, long desclcid)
04421 {
04422 long i;
04423 fprintf(out, "\nVALUES\nSCALAR NODES ");
04424 fprintf(out, "%s %ld\n", descr, desclcid);
04425 for (i=0; i < Mt->nn; i++)
04426 fprintf(out, "%ld %e\n", i, val[i]);
04427 }
04428
04429
04430
04431
04432
04433
04434
04435
04436
04437
04438
04439
04440
04441
04442
04443
04444 void write_nodscalar(FILE *out, long dir, const char *desclcid)
04445 {
04446 long i;
04447 fprintf(out, "\nVALUES\nSCALAR NODES ");
04448 switch (dir)
04449 {
04450 case -1:
04451 fprintf(out, "%s %s\n", "peps_1", desclcid);
04452 for (i = 0; i < Mt->nn; i++)
04453 {
04454 if (Mt->nodes[i].pstra)
04455 fprintf(out, "%ld % e\n", i, Mt->nodes[i].pstra[-dir-1]);
04456 else
04457 fprintf(out, "%ld 0.0\n",i);
04458 }
04459 break;
04460 case -2:
04461 fprintf(out, "%s %s\n", "peps_2", desclcid);
04462 for (i = 0; i < Mt->nn; i++)
04463 {
04464 if (Mt->nodes[i].pstra)
04465 fprintf(out, "%ld % e\n", i, Mt->nodes[i].pstra[-dir-1]);
04466 else
04467 fprintf(out, "%ld 0.0\n",i);
04468 }
04469 break;
04470 case -3:
04471 fprintf(out, "%s %s\n", "peps_3", desclcid);
04472 for (i = 0; i < Mt->nn; i++)
04473 {
04474 if (Mt->nodes[i].pstra)
04475 fprintf(out, "%ld % e\n", i, Mt->nodes[i].pstra[-dir-1]);
04476 else
04477 fprintf(out, "%ld 0.0\n",i);
04478 }
04479 break;
04480 case 0:
04481 fprintf(out, "%s %s\n", "psig_1", desclcid);
04482 for (i = 0; i < Mt->nn; i++)
04483 {
04484 if (Mt->nodes[i].pstre)
04485 fprintf(out, "%ld % e\n", i, Mt->nodes[i].pstre[dir]);
04486 else
04487 fprintf(out, "%ld 0.0\n",i);
04488 }
04489 break;
04490 case 1:
04491 fprintf(out, "%s %s\n", "psig_2", desclcid);
04492 for (i = 0; i < Mt->nn; i++)
04493 {
04494 if (Mt->nodes[i].pstre)
04495 fprintf(out, "%ld % e\n", i, Mt->nodes[i].pstre[dir]);
04496 else
04497 fprintf(out, "%ld 0.0\n",i);
04498 }
04499 break;
04500 case 2:
04501 fprintf(out, "%s %s\n", "psig_3", desclcid);
04502 for (i = 0; i < Mt->nn; i++)
04503 {
04504 if (Mt->nodes[i].pstre)
04505 fprintf(out, "%ld % e\n", i, Mt->nodes[i].pstre[dir]);
04506 else
04507 fprintf(out, "%ld 0.0\n",i);
04508 }
04509 break;
04510 case 3 :
04511 fprintf(out, "%s %s\n", "tau_max", desclcid);
04512 for (i = 0; i < Mt->nn; i++)
04513 {
04514 if (Mt->nodes[i].pstre)
04515 fprintf(out, "%ld % e\n", i, (Mt->nodes[i].pstre[2]-Mt->nodes[i].pstre[0])/2);
04516 else
04517 fprintf(out, "%ld 0.0\n",i);
04518 }
04519 break;
04520 default:
04521 break;
04522 }
04523 }
04524
04525
04526
04527
04528
04529
04530
04531
04532
04533
04534
04535
04536
04537
04538
04539
04540 void write_elemscalar(FILE *out,double *val, const char *descr, const char *desclcid)
04541 {
04542 long i;
04543 fprintf(out, "\nVALUES\nSCALAR ELEMENTS ");
04544 fprintf(out, "%s %s\n", descr, desclcid);
04545 for (i=0; i < Mt->ne; i++){
04546 if (Gtm->leso[i]==0)
04547 continue;
04548
04549 fprintf(out, "%ld %e\n", i, val[i]);
04550 }
04551 }
04552
04553
04554
04555
04556
04557
04558
04559
04560
04561
04562
04563
04564
04565
04566
04567
04568 void write_deflection (FILE *out,long lcid,long dir,const char *desclcid)
04569 {
04570 long i;
04571 double *r;
04572
04573 fprintf(out, "\nVALUES\nSCALAR NODES ");
04574 fprintf(out, "w %s\n", desclcid);
04575 for (i=0; i < Mt->nn; i++){
04576 r = new double [Gtm->gnodes[i].ndofn];
04577 noddispl (lcid, r, i);
04578 fprintf(out, "%ld", i);
04579 fprintf(out, " % e\n", r[dir]);
04580 delete [] r;
04581 }
04582 }
04583
04584
04585
04586
04587
04588
04589
04590
04591 void print_displacements (FILE *out,long lcid)
04592 {
04593 long i,j,k,ndofn;
04594 double *r,*d;
04595
04596 r=Lsrs->lhs+lcid*Ndofm;
04597 switch (Mp->tprob){
04598 case linear_statics:{
04599 d=Mb->lc[lcid].pd;
04600 break;
04601 }
04602 case mat_nonlinear_statics:{
04603 d=Mb->lc[lcid].pd;
04604 break;
04605 }
04606 case geom_nonlinear_statics:{
04607 d=Mb->lc[lcid].pd;
04608 break;
04609 }
04610 case earth_pressure:{
04611 d=Mb->lc[lcid].pd;
04612 break;
04613 }
04614 case eigen_dynamics:{
04615 break;
04616 }
04617 case forced_dynamics:{
04618 d=NULL;
04619 break;
04620 }
04621 case growing_mech_structure:
04622 case mech_timedependent_prob:{
04623 d=NULL;
04624 break;
04625 }
04626 case layered_linear_statics:{
04627 d=Mb->lc[lcid].pd;
04628 break;
04629 }
04630 case lin_floating_subdomain:{
04631 d=Mb->lc[lcid].pd;
04632 break;
04633 }
04634 case nonlin_floating_subdomain:{
04635 d=Mb->lc[lcid].pd;
04636 break;
04637 }
04638 default:{
04639 print_err("unknown problem type is required",__FILE__,__LINE__, __func__);
04640 }
04641 }
04642
04643
04644 fprintf (out,"\n\n\n NODE DISPLACEMENTS \n");
04645 for (i=0;i<Mt->nn;i++){
04646 fprintf (out,"\n\n node %ld",i+1);
04647 ndofn=Mt->give_ndofn (i);
04648 for (j=0;j<ndofn;j++){
04649 k=Gtm->give_dof (i,j);
04650 if (k>0) fprintf (out,"\n gen. displ. %15.12e",r[k-1]);
04651 if (k==0) fprintf (out,"\n gen. displ. %15.12e",0.0);
04652 if (k<0) fprintf (out,"\n gen. displ. %15.12e",d[0-k-1]);
04653 }
04654 }
04655
04656
04657
04658
04659
04660
04661
04662
04663
04664
04665 }
04666
04667 void print_multipliers (FILE *out)
04668 {
04669 long i,j,k,l,nl,nmult,*cn;
04670 double *d;
04671
04672 d=Lsrs->lhs;
04673
04674 if (Mp->tprob!=layered_linear_statics){
04675 print_err("there are no multipliers in your problem.", __FILE__, __LINE__, __func__);
04676 abort ();
04677 }
04678
04679 if (Mp->tprob==layered_linear_statics){
04680 fprintf (out,"\n\n");
04681 for (i=0;i<Mt->nln;i++){
04682 nl=Gtm->lgnodes[i].nl;
04683 nmult=Gtm->give_nmult (i);
04684 cn = new long [nmult];
04685 fprintf (out,"\n\n layered node %ld",i);
04686 for (j=1;j<nl;j++){
04687 Gtm->give_mult_code_numbers (i,j,cn);
04688 for (k=0;k<nmult;k++){
04689 l=cn[k];
04690 if (l==0) fprintf (out,"\n mult. %15.12e",0.0);
04691 if (l>1) fprintf (out,"\n mult. %15.12e",d[l-1]);
04692 }
04693 }
04694 delete [] cn;
04695 }
04696 }
04697 }
04698
04699
04700
04701 void print_strains_old (FILE *out,long lcid)
04702 {
04703
04704
04705
04706
04707
04708
04709
04710
04711
04712
04713
04714
04715
04716
04717
04718
04719
04720
04721
04722
04723
04724
04725
04726
04727
04728
04729
04730
04731
04732
04733
04734
04735
04736
04737
04738
04739
04740
04741
04742
04743
04744
04745
04746
04747
04748
04749
04750
04751
04752
04753
04754
04755
04756
04757
04758
04759
04760
04761
04762
04763
04764
04765
04766
04767
04768
04769
04770
04771
04772
04773
04774
04775
04776
04777
04778
04779
04780
04781
04782
04783
04784
04785
04786
04787
04788
04789
04790
04791
04792
04793
04794
04795
04796
04797
04798
04799
04800
04801
04802
04803
04804
04805
04806
04807
04808
04809
04810
04811
04812
04813
04814
04815
04816
04817
04818
04819
04820
04821
04822
04823
04824
04825
04826
04827
04828
04829
04830
04831
04832
04833
04834
04835
04836
04837
04838
04839
04840
04841
04842
04843
04844
04845
04846
04847
04848
04849
04850
04851
04852
04853
04854
04855
04856
04857
04858
04859
04860
04861
04862
04863
04864 }
04865
04866
04867
04868
04869
04870
04871
04872
04873
04874
04875
04876
04877 void print_strains_nodes (FILE *out,long lcid)
04878 {
04879 long i,j,nn,ncomp;
04880 nn = Mt->nn;
04881
04882 fprintf (out,"\n\n\n\n STRAINS IN NODES IN LOAD CASE %ld\n",lcid+1);
04883 for (i=0;i<nn;i++){
04884 ncomp=Mt->nodes[i].ncompstr;
04885 fprintf (out,"\n node %4ld",i+1);
04886 for (j=0;j<ncomp;j++){
04887 fprintf (out," %e",Mt->nodes[i].strain[lcid*ncomp+j]);
04888 }
04889 }
04890 }
04891
04892
04893
04894
04895
04896
04897
04898
04899
04900 void print_strains_udp (FILE *out,long eid)
04901 {
04902 long i,j,nape,ncomp;
04903
04904 nape=Mm->stra.nape[eid];
04905 ncomp=Mt->give_tncomp(eid);
04906 fprintf (out,"\n\n strains on element number %ld",eid);
04907 for (i=0;i<nape;i++){
04908 fprintf (out,"\n point number %ld\n",i);
04909 for (j=0;j<ncomp;j++){
04910 fprintf (out," %20.10f",Mm->stra.ev[eid][i][j]);
04911 }
04912 }
04913 }
04914
04915
04916
04917
04918
04919
04920
04921
04922
04923
04924
04925
04926
04927
04928
04929
04930
04931
04932
04933
04934
04935
04936
04937
04938
04939
04940
04941
04942
04943
04944
04945
04946
04947
04948
04949
04950
04951
04952
04953
04954
04955
04956
04957
04958
04959
04960
04961
04962
04963
04964
04965
04966
04967
04968
04969
04970
04971
04972 void print_stresses_nodes (FILE *out,long lcid)
04973 {
04974 long i,j,nn,ncomp;
04975 nn = Mt->nn;
04976
04977 fprintf (out,"\n\n\n\n STRESSES IN NODES IN LOAD CASE %ld\n",lcid+1);
04978 for (i=0;i<nn;i++){
04979 ncomp=Mt->nodes[i].ncompstr;
04980 fprintf (out,"\n node %4ld",i+1);
04981 for (j=0;j<ncomp;j++){
04982 fprintf (out," %e",Mt->nodes[i].stress[lcid*ncomp+j]);
04983 }
04984 }
04985 }
04986
04987
04988
04989
04990
04991
04992
04993
04994
04995 void print_stresses_udp (FILE *out,long eid)
04996 {
04997 long i,j,nape,ncomp;
04998
04999 nape=Mm->stra.nape[eid];
05000 ncomp=Mt->give_tncomp(eid);
05001 fprintf (out,"\n\n stresses on element number %ld",eid);
05002 for (i=0;i<nape;i++){
05003 fprintf (out,"\n point number %ld\n",i);
05004 for (j=0;j<ncomp;j++){
05005 fprintf (out," %20.10f",Mm->stre.ev[eid][i][j]);
05006 }
05007 }
05008 }
05009
05010
05011
05012
05013
05014
05015
05016
05017
05018
05019
05020
05021
05022
05023
05024
05025
05026
05027
05028
05029
05030
05031
05032
05033
05034
05035
05036
05037
05038
05039
05040
05041
05042
05043
05044
05045
05046
05047
05048
05049
05050
05051
05052
05053
05054
05055
05056
05057
05058
05059
05060
05061
05062
05063
05064
05065
05066
05067
05068
05069
05070 void print_stresses_old (FILE *out,long lcid)
05071 {
05072
05073
05074
05075
05076
05077
05078
05079
05080
05081
05082
05083
05084
05085
05086
05087
05088
05089
05090
05091
05092
05093
05094
05095
05096
05097
05098
05099
05100
05101
05102
05103
05104
05105
05106
05107
05108
05109
05110
05111
05112
05113
05114
05115
05116
05117
05118
05119
05120
05121
05122
05123
05124
05125
05126
05127
05128
05129
05130
05131
05132
05133
05134
05135
05136
05137
05138
05139
05140
05141
05142
05143
05144
05145
05146
05147
05148
05149
05150
05151
05152
05153
05154
05155
05156
05157
05158
05159
05160
05161
05162
05163
05164
05165
05166
05167
05168
05169
05170
05171
05172
05173
05174
05175
05176
05177
05178
05179
05180
05181
05182
05183
05184
05185
05186
05187
05188
05189
05190
05191
05192
05193
05194
05195
05196
05197
05198
05199
05200
05201
05202
05203
05204
05205
05206
05207
05208
05209
05210
05211
05212
05213
05214
05215
05216
05217
05218
05219
05220
05221
05222
05223
05224
05225
05226
05227
05228
05229
05230
05231
05232
05233 }
05234
05235 void print_other (FILE *out,long lcid)
05236
05237 {
05238 long i,j;
05239
05240 fprintf (out,"\n\n\n OTHER \n\n");
05241 for (i=0;i<Mm->tnip;i++){
05242 fprintf (out,"\n");
05243 for (j=0;j<Mm->ip[i].ncompother;j++){
05244 fprintf (out," %f",Mm->ip[i].eqother[j]);
05245 }
05246 }
05247
05248
05249
05250
05251
05252
05253
05254
05255
05256
05257
05258
05259
05260
05261
05262
05263
05264
05265
05266
05267
05268
05269
05270
05271
05272
05273
05274
05275
05276
05277
05278
05279
05280
05281
05282
05283
05284
05285
05286
05287
05288
05289
05290
05291
05292
05293
05294
05295
05296
05297
05298
05299
05300
05301
05302
05303
05304
05305
05306
05307
05308
05309
05310
05311
05312
05313
05314
05315
05316
05317
05318
05319
05320
05321
05322
05323
05324
05325
05326
05327
05328
05329
05330
05331 }
05332
05333 void print_intforces (FILE *out,double *fi)
05334 {
05335 long i,j,k,ndofn;
05336
05337 fprintf (out,"\n\n\n NODE INTERNAL FORCES \n");
05338 for (i=0;i<Mt->nn;i++){
05339 fprintf (out,"\n\n node %ld ",i+1);
05340 ndofn=Mt->give_ndofn (i);
05341 for (j=0;j<ndofn;j++){
05342 k=Gtm->give_dof (i,j);
05343 if (k>0) fprintf (out," %15.12e",fi[k-1]);
05344 else fprintf (out,"%15s", "reakce");
05345 }
05346 fprintf (out,"\n");
05347 }
05348
05349 }
05350
05351 void print_reactions (FILE *out,long lcid)
05352 {
05353 long i,j,nn,ndofn;
05354 nn=Mt->nn;
05355
05356 fprintf (out,"\n\n\n REACTIONS \n");
05357 for (i=0;i<nn;i++){
05358 if (Mt->nodes[i].react==1){
05359 ndofn=Mt->give_ndofn (i);
05360 fprintf (out,"\n%4ld",i+1);
05361 for (j=0;j<ndofn;j++){
05362 fprintf (out," %20.10f",Mt->nodes[i].r[j]);
05363 }
05364 }
05365 }
05366 }
05367
05368 void print_eigenvalues (double *w)
05369 {
05370 long i,n;
05371
05372
05373 n=Mp->eigsol.neigv;
05374
05375 fprintf (Out,"\n\n\n EIGENVALUES \n");
05376
05377 for (i=0;i<n;i++){
05378 fprintf (Out,"\n omega %4ld %30.20f %30.20f",i+1,w[i],sqrt(w[i]));
05379
05380 }
05381
05382 }
05383
05384 void print_eigenvectors ()
05385 {
05386 long i,j,k,n,m;
05387 double *v;
05388
05389 n=Mp->eigsol.neigv;
05390 m=Ndofm;
05391
05392 v=Lsrs->give_lhs (0);
05393
05394
05395 fprintf (Out,"\n\n\n EIGENVECTORS \n");
05396 fprintf (Out,"\n\n");
05397 k=0;
05398 for (i=0;i<n;i++){
05399 fprintf (Out,"\n\n Eigenvector number %ld",i+1);
05400 for (j=0;j<m;j++){
05401 fprintf (Out,"\n %e",v[k]);
05402 k++;
05403 }
05404 }
05405
05406
05407
05408
05409
05410
05411
05412
05413
05414
05415
05416
05417
05418
05419
05420
05421
05422
05423
05424
05425
05426
05427
05428
05429
05430
05431
05432
05433
05434
05435
05436
05437
05438
05439
05440
05441
05442
05443
05444
05445
05446
05447
05448
05449
05450
05451
05452
05453
05454
05455
05456
05457
05458
05459
05460
05461
05462
05463
05464
05465
05466
05467
05468
05469
05470
05471
05472
05473
05474
05475
05476
05477
05478
05479
05480
05481
05482
05483
05484
05485
05486
05487
05488
05489
05490
05491
05492
05493
05494
05495
05496
05497
05498
05499
05500
05501
05502
05503
05504
05505
05506 }
05507
05508
05509 void print_eigenvect_martin (FILE *out)
05510 {
05511 long i,j,n,ndofn;
05512 double *v;
05513
05514
05515
05516
05517
05518
05519
05520
05521
05522
05523
05524
05525
05526
05527
05528
05529
05530
05531
05532
05533
05534
05535
05536
05537
05538
05539
05540
05541
05542 n=Mp->eigsol.neigv;
05543
05544 for (i=0;i<n;i++){
05545 v=Lsrs->give_lhs (i);
05546 fprintf (out,"\n\n Eigenvector number %ld",i+1);
05547 for (j=0;j<Mt->nn;j++){
05548 fprintf (out,"\n%ld",j);
05549 ndofn=Mt->give_ndofn (j);
05550 if (ndofn==2){
05551 if (Gtm->give_dof (j,0)>0) fprintf (out," %e",v[Gtm->give_dof (j,0)-1]);
05552 else fprintf (out," 0.0");
05553 if (Gtm->give_dof (j,1)>0) fprintf (out," %e",v[Gtm->give_dof (j,1)-1]);
05554 else fprintf (out," 0.0");
05555 fprintf (out," 0.0");
05556 }
05557 if (ndofn==3){
05558 if (Gtm->give_dof (j,0)>0) fprintf (out," %e",v[Gtm->give_dof (j,0)-1]);
05559 else fprintf (out," 0.0");
05560 if (Gtm->give_dof (j,1)>0) fprintf (out," %e",v[Gtm->give_dof (j,1)-1]);
05561 else fprintf (out," 0.0");
05562 if (Gtm->give_dof (j,2)>0) fprintf (out," %e",v[Gtm->give_dof (j,2)-1]);
05563 else fprintf (out," 0.0");
05564 }
05565 }
05566 }
05567
05568 }
05569
05570
05571
05572
05573
05574
05575
05576
05577
05578
05579
05580
05581
05582
05583
05584
05585
05586
05587
05588
05589
05590
05591
05592
05593
05594
05595
05596
05597
05598
05599
05600
05601
05602
05603
05604
05605
05606
05607
05608
05609
05610
05611
05612
05613
05614
05615
05616
05617
05618
05619
05620 void print_valel (gtopology *gt, probdesc *mp, const char *file, char *caption, double *valel,char flag)
05621 {
05622 long i,j;
05623 double *valnod;
05624 char *fff;
05625 char *tcaption;
05626
05627
05628 tcaption = new char[strlen(caption)];
05629 i = 0;
05630 while (caption[i] != '\0')
05631 {
05632 tcaption[i] = caption[i];
05633 i++;
05634 }
05635 tcaption[i] = '\0';
05636
05637 fff = NULL;
05638 valnod = new double [gt->nn];
05639
05640 if (file == NULL){
05641 fff = new char [255];
05642 if (flag == 'e') sprintf (fff,"%s%s.ex",mp->path,tcaption);
05643 if (flag == 'd') sprintf (fff,"%s%s.dx",mp->path,tcaption);
05644
05645 file = fff;
05646 }
05647
05648 if (gt->nadjelnod == NULL)
05649 gt->adjacelem (Out);
05650
05651 for (i=0;i<gt->nn;i++){
05652 valnod[i] = 0.0;
05653 for (j=0;j<gt->nadjelnod[i];j++){
05654 valnod[i] += valel[gt->adjelnod[i][j]];
05655 }
05656 valnod[i] /= (double)gt->nadjelnod[i];
05657 }
05658
05659 if (flag == 'e')
05660 print_ex (gt,file,valnod,valel);
05661
05662 if (flag == 'd'){
05663 long dimindex = 1;
05664 print_dx (gt,file,&valnod,&valel,'n',&dimindex,&tcaption,1);
05665 }
05666
05667 if (fff != NULL){ file = NULL; delete [] fff; }
05668 delete [] valnod;
05669 delete [] tcaption;
05670 }
05671
05672
05673
05674
05675
05676
05677
05678
05679
05680
05681
05682
05683
05684
05685
05686
05687
05688
05689
05690
05691
05692
05693
05694
05695
05696
05697
05698
05699
05700
05701
05702
05703
05704
05705
05706
05707
05708
05709
05710
05711
05712
05713
05714
05715
05716
05717
05718
05719
05720
05721
05722
05723
05724
05725
05726
05727
05728
05729
05730
05731
05732
05733
05734
05735
05736
05737
05738
05739
05740
05741
05742
05743
05744
05745
05746
05747
05748
05749
05750
05751
05752
05753
05754 void print_default_dx (gtopology *gt,probdesc *mp,mechtop *mt,mechmat *mm,long lcid,const char *file)
05755 {
05756 long i,j,nn,ne,vali,capi,ndofn,ncomp,dd;
05757 double **valnod,**valel,*r;
05758 char **caption;
05759 long d[6],*dimindex;
05760
05761 dd = 1;
05762
05763 if (!mp->detnodstrain){
05764
05765
05766 }
05767
05768 switch (mt->elements[0].te){
05769 case planeelementlt:
05770 case planeelementlq:
05771 case planeelementqt: { ndofn = 2; d[0]= 1; d[1]= 1; d[2]= 0; d[3]= 0; d[4]= 0; d[5]= 1; break;}
05772 case planeelementqq: { ndofn = 2; d[0]= 1; d[1]= 1; d[2]= 0; d[3]= 0; d[4]= 0; d[5]= 1; break;}
05773 case lineartet:
05774 case quadrtet:
05775 case linearhex: { ndofn = 3; d[0]= 1; d[1]= 1; d[2]= 1; d[3]= 1; d[4]= 1; d[5]= 1; break;}
05776 default:{
05777 print_err("unknown element type is required", __FILE__,__LINE__, __func__);
05778 }
05779 }
05780
05781 vali=0;
05782 nn = gt->nn;
05783 ne = gt->ne;
05784 ncomp = mt->give_tncomp (0);
05785
05786 r = new double [ndofn];
05787 valel = new double* [ne];
05788 valnod = new double* [(1+dd)*ndofn+2*ncomp];
05789 caption = new char* [ 1+dd *ndofn+2*ncomp];
05790 dimindex = new long [ 1+dd *ndofn+2*ncomp];
05791
05792
05793
05794 for (i=0;i<ndofn;i++){
05795 valnod[vali+i] = new double [nn];
05796 if (dd) valnod[vali+i+ndofn] = new double [nn];
05797 }
05798
05799 for (i=0;i<nn;i++){
05800 noddispl (lcid,r,i);
05801 for (j=0;j<ndofn;j++){
05802 valnod[vali+j][i] = r[j];
05803 if (dd) valnod[vali+j+ndofn][i] = r[j];
05804 }
05805 }
05806 vali+=(1+dd)*ndofn;
05807
05808
05809 for (i=0;i<ncomp;i++){
05810 valnod[vali] = new double [nn];
05811 for (j=0;j<nn;j++)
05812 valnod[vali][j] = mt->nodes[j].strain[lcid*ncomp+i];
05813
05814 vali++;
05815 }
05816
05817
05818 for (i=0;i<ncomp;i++){
05819 valnod[vali] = new double [nn];
05820 for (j=0;j<nn;j++)
05821 valnod[vali][j] = mt->nodes[j].stress[lcid*ncomp+i];
05822
05823 vali++;
05824 }
05825
05826
05827 for (i=0;i<ne;i++)
05828 if (planeelementqq == mt->elements[i].te){
05829 valel[i] = new double[vali];
05830 Peqq->midpoints (i,dd,valel[i]);
05831 }
05832 else
05833 valel[i] = NULL;
05834
05835
05836 capi=0; dimindex[capi]=ndofn; caption[capi] = new char[4]; sprintf (caption[capi],"def"); capi++;
05837
05838 if (dd && ndofn>0) { dimindex[capi]=1; caption[capi] = new char[5]; sprintf (caption[capi],"defx"); capi++; }
05839 if (dd && ndofn>1) { dimindex[capi]=1; caption[capi] = new char[5]; sprintf (caption[capi],"defy"); capi++; }
05840 if (dd && ndofn>2) { dimindex[capi]=1; caption[capi] = new char[5]; sprintf (caption[capi],"defz"); capi++; }
05841
05842 if (d[0]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"epsxx"); capi++; }
05843 if (d[1]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"epsyy"); capi++; }
05844 if (d[2]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"epszz"); capi++; }
05845 if (d[3]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"epsyz"); capi++; }
05846 if (d[4]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"epsxz"); capi++; }
05847 if (d[5]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"epsxy"); capi++; }
05848
05849 if (d[0]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"sigxx"); capi++; }
05850 if (d[1]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"sigyy"); capi++; }
05851 if (d[2]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"sigzz"); capi++; }
05852 if (d[3]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"sigyz"); capi++; }
05853 if (d[4]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"sigxz"); capi++; }
05854 if (d[5]) { dimindex[capi]=1; caption[capi] = new char[6]; sprintf (caption[capi],"sigxy"); capi++; }
05855
05856 print_dx (gt,file,valnod,valel,'t',dimindex,caption,capi);
05857
05858
05859 delete [] r;
05860 delete [] dimindex;
05861 for (i=0;i<vali;i++) delete [] valnod[i]; delete [] valnod;
05862 for (i=0;i<capi;i++) delete [] caption[i]; delete [] caption;
05863 for (i=0;i<ne;i++) if (valel[i]!=NULL) delete [] valnod[i];
05864 delete [] valel;
05865 }
05866
05867
05868
05869
05870
05871
05872
05873
05874
05875
05876
05877
05878
05879
05880
05881
05882
05883
05884
05885
05886
05887
05888
05889
05890
05891
05892
05893
05894
05895
05896
05897
05898
05899
05900
05901
05902
05903
05904
05905
05906
05907
05908
05909
05910
05911
05912
05913
05914
05915
05916
05917
05918
05919
05920
05921
05922
05923
05924
05925
05926
05927
05928
05929
05930
05931
05932
05933
05934
05935
05936
05937
05938
05939
05940
05941
05942
05943
05944
05945
05946
05947
05948 void aux_mech_nonlin_print (FILE *aux,double *r,double l)
05949 {
05950
05951
05952
05953
05954
05955
05956
05957
05958
05959
05960
05961
05962
05963
05964
05965 }
05966
05967 void aux_mech_time_print (FILE *aux,double *r,double l)
05968 {
05969
05970
05971
05972
05973
05974
05975
05976
05977
05978
05979
05980
05981
05982
05983
05984
05985 }