00001 #include <math.h>
00002 #include <string.h>
00003 #include "mechmat.h"
00004 #include "global.h"
00005 #include "globmat.h"
00006 #include "elemhead.h"
00007 #include "vecttens.h"
00008 #include "intpoints.h"
00009 #include "element.h"
00010 #include "node.h"
00011 #include "elastisomat.h"
00012 #include "elastgmat3d.h"
00013 #include "splas1d.h"
00014 #include "j2flow.h"
00015 #include "j2flow2.h"
00016 #include "microM4.h"
00017 #include "microSIM.h"
00018 #include "microfiber.h"
00019 #include "mohrc.h"
00020 #include "mohrcparab.h"
00021 #include "boermat.h"
00022 #include "drprag.h"
00023 #include "drprag2.h"
00024 #include "camclay.h"
00025 #include "camclaycoup.h"
00026 #include "shefplast.h"
00027 #include "hissplas.h"
00028 #include "glasgmech.h"
00029 #include "glasgowdam.h"
00030 #include "creepdam.h"
00031 #include "timeswmat.h"
00032 #include "effstress.h"
00033 #include "simviscous.h"
00034 #include "lemaitre.h"
00035 #include "ortodam.h"
00036 #include "ortodamrot.h"
00037 #include "anisodam.h"
00038 #include "anisodamrot.h"
00039 #include "scaldam.h"
00040 #include "scaldamcc.h"
00041 #include "aci78.h"
00042 #include "cebfip78.h"
00043 #include "graphmat.h"
00044 #include "varelastisomat.h"
00045 #include "geoelast.h"
00046 #include "creep.h"
00047 #include "creep_b3.h"
00048 #include "creep_rspec.h"
00049 #include "creepb.h"
00050 #include "creep_dpl.h"
00051 #include "creep_effym.h"
00052 #include "shrinkmat.h"
00053 #include "consol.h"
00054 #include "winpast.h"
00055 #include "nonlocmicroM4.h"
00056 #include "therisomat.h"
00057 #include "nonlocplast.h"
00058 #include "nonlocdamg.h"
00059 #include "damplast.h"
00060 #include "tensor.h"
00061 #include "visplast.h"
00062 #include "chen.h"
00063 #include "lenjonesmat.h"
00064 #include "contactmat.h"
00065 #include "cebfipcontactmat.h"
00066 #include "elemswitch.h"
00067 #include "relaxeuroc.h"
00068 #include "therisomat.h"
00069 #include "shrinkmat.h"
00070 #include "hypoplast.h"
00071
00072
00073
00074
00075
00076
00077
00078
00079
00080 mechmat::mechmat ()
00081 {
00082
00083 nmt=0;
00084
00085 tnip=0;
00086
00087
00088 plast=0;
00089
00090
00091 max_nncompo=0;
00092
00093 max_encompo=0;
00094
00095 max_ncompstrn=0;
00096
00097 max_ncompstre=0;
00098
00099
00100 ip=NULL;
00101
00102 ic=NULL;
00103
00104 elip = NULL;
00105
00106
00107 eliso=NULL; elgm3d=NULL;
00108 spl1d=NULL; j2f=NULL; j2f2=NULL; mcoul=NULL; mcpar=NULL;
00109 boerm=NULL; drprm=NULL; drprm2=NULL; cclay=NULL; cclayc=NULL; shpl=NULL; hisspl=NULL;
00110 glasgmat=NULL; glasgdam=NULL; crdam=NULL; tswmat=NULL; effstr=NULL;
00111 mpM4=NULL; mpSIM=NULL; mpfib=NULL; nmpM4=NULL;
00112 svis=NULL; svipl=NULL; lmtr=NULL;
00113 scdam=NULL; anidam=NULL; anidamrot=NULL; ortdam=NULL;
00114 ortdam=NULL; ortdamrot=NULL;
00115 aci78mod=NULL; cebfip78mod=NULL;
00116 grmat=NULL; geoel=NULL; veliso=NULL;
00117 crb3=NULL;
00118 crrs=NULL;
00119 crdpl=NULL;
00120 crbaz=NULL; csol=NULL; wpast=NULL;
00121 tidilat=NULL; relaxec=NULL;
00122 nlplast=NULL; nldamg=NULL;
00123
00124 dampl = NULL; visplas = NULL;
00125 chplast = NULL;
00126 lenjon = NULL;
00127 conmat=NULL;
00128 cebfipconmat=NULL;
00129
00130
00131 shmat=NULL;
00132
00133
00134 ipv = NULL;
00135
00136
00137 eigstrid = NULL;
00138
00139 eigstrains = NULL;
00140
00141 eigstresses = NULL;
00142
00143 tempstrains = NULL;
00144
00145 nonmechq = NULL;
00146 nnmq = 0;
00147 nmqo = NULL;
00148 for (long i=0; i<tnknmq; i++)
00149 nmqid[i] = -1;
00150
00151 mtype = NULL;
00152 numtype = NULL;
00153 }
00154
00155
00156
00157
00158
00159
00160
00161
00162
00163 mechmat::~mechmat ()
00164 {
00165 long i;
00166
00167 delete [] ip;
00168 delete [] elip;
00169 delete [] eliso; delete [] elgm3d;
00170
00171 delete [] spl1d; delete [] j2f; delete [] j2f2;
00172 delete [] mcoul; delete [] mcpar;
00173
00174 delete [] boerm; delete [] drprm; delete [] drprm2; delete [] cclay; delete [] shpl; delete [] hisspl;
00175
00176 delete [] glasgmat; delete [] glasgdam; delete [] crdam; delete [] tswmat; delete [] effstr;
00177
00178 delete [] mpM4; delete [] mpSIM; delete [] mpfib; delete [] nmpM4;
00179
00180 delete [] svis;
00181 delete [] lmtr;
00182
00183 delete [] scdam; delete [] anidam; delete [] anidamrot;
00184 delete [] ortdam; delete [] ortdamrot;
00185
00186 delete [] aci78mod; delete [] cebfip78mod;
00187 delete [] grmat; delete [] geoel; delete [] veliso;
00188 delete [] crb3; delete [] crrs; delete [] crdpl;
00189 delete [] crbaz; delete [] csol; delete [] wpast;
00190 delete [] tidilat; delete [] relaxec;
00191 delete [] nlplast; delete [] nldamg;
00192 delete [] dampl; delete [] visplas;
00193 delete [] chplast;
00194 delete [] lenjon; delete [] conmat; delete [] cebfipconmat;
00195
00196
00197 delete [] shmat;
00198
00199 if (ic != NULL)
00200 {
00201 for (i = 0; i < tnip; i++)
00202 delete [] ic[i];
00203 }
00204 delete [] ic;
00205
00206 delete [] ipv;
00207
00208 if (nonmechq != NULL){
00209 for (i=0;i<nnmq;i++){
00210 delete [] nonmechq[i];
00211 }
00212 delete [] nonmechq;
00213 }
00214
00215 delete [] nmqo;
00216
00217 if (eigstrains != NULL){
00218 for (i=0;i<tnip;i++){
00219 delete [] eigstrains[i];
00220 }
00221 delete [] eigstrains;
00222 }
00223
00224 if (eigstrid != NULL){
00225 for (i=0;i<tnip;i++){
00226 delete [] eigstrid[i];
00227 }
00228 delete [] eigstrid;
00229 }
00230
00231 if (eigstresses != NULL){
00232 for (i=0;i<tnip;i++){
00233 delete [] eigstresses[i];
00234 }
00235 delete [] eigstresses;
00236 }
00237 if (tempstrains != NULL){
00238 for (i=0;i<tnip;i++){
00239 delete [] tempstrains[i];
00240 }
00241 delete [] tempstrains;
00242 }
00243
00244 delete [] mtype;
00245 delete [] numtype;
00246 }
00247
00248
00249
00250
00251
00252
00253
00254
00255
00256
00257
00258 long mechmat::intpnum (void)
00259 {
00260 long i,j,k,n,nb;
00261
00262 n=0;
00263 for (i=0;i<Mt->ne;i++){
00264 nb=Mt->give_nb (i);
00265 Mt->elements[i].ipp = new long* [nb];
00266 for (j=0;j<nb;j++){
00267 Mt->elements[i].ipp[j] = new long [nb];
00268 }
00269 for (j=0;j<nb;j++){
00270 for (k=0;k<nb;k++){
00271 Mt->elements[i].ipp[j][k]=n;
00272 n+=Mt->give_nip (i,j,k);
00273 }
00274 }
00275 }
00276 return n;
00277 }
00278
00279
00280
00281
00282
00283
00284
00285
00286
00287
00288 void mechmat::intpointalloc ()
00289 {
00290 ip = new intpoints [tnip];
00291 }
00292
00293
00294
00295
00296
00297
00298
00299
00300
00301
00302 void mechmat::alloc_ic ()
00303 {
00304 long i,n;
00305
00306 ic = new double* [tnip];
00307 for (i=0;i<tnip;i++){
00308 n = ip[i].ncompstr;
00309 ic[i] = new double[3+n];
00310 }
00311 }
00312
00313
00314
00315
00316
00317
00318
00319
00320
00321
00322
00323
00324
00325
00326 void mechmat::alloc_nonmechq(long n)
00327 {
00328 long i, j;
00329
00330 nnmq = n;
00331 nonmechq = new double*[nnmq];
00332
00333 for (i=0; i<nnmq; i++)
00334 {
00335 nonmechq[i] = new double[tnip];
00336
00337 memset(nonmechq[i], 0, sizeof(*nonmechq[i])*tnip);
00338 }
00339 if (nmqo == NULL)
00340 {
00341 nmqo = new nonmechquant[nnmq];
00342 j=0;
00343 for (i=0; i<tnknmq; i++)
00344 {
00345 if (nmqid[i] > -1)
00346 {
00347 nmqo[j] = nonmechquant(i+1);
00348 j++;
00349 }
00350 }
00351 if (j != nnmq)
00352 {
00353 print_err("indices of used non-mechanical qunatities must be initialized before allocation", __FILE__, __LINE__, __func__);
00354 abort();
00355 }
00356 }
00357 return;
00358 }
00359
00360
00361
00362
00363
00364
00365
00366
00367
00368
00369
00370
00371
00372
00373
00374 void mechmat::readip (XFILE *in)
00375 {
00376 long i;
00377 for (i=0;i<tnip;i++){
00378 ip[i].read (in);
00379 }
00380 }
00381
00382
00383
00384
00385
00386
00387
00388
00389
00390
00391
00392 void mechmat::intpointinit ()
00393 {
00394 long i, j, ii, jj, k, nb, nip, ipp;
00395
00396
00397 elip = new long[tnip];
00398
00399 for (i=0;i<Mt->ne;i++){
00400 nb=Mt->give_nb (i);
00401 for (ii=0;ii<nb;ii++){
00402 for (jj=0;jj<nb;jj++){
00403 nip=Mt->give_nip(i,ii,jj);
00404 ipp=Mt->elements[i].ipp[ii][jj];
00405 for (j=0;j<nip;j++){
00406 if (elip != NULL)
00407 elip[ipp] = i;
00408 ip[ipp].nm = Mt->elements[i].nm;
00409 ip[ipp].tm = new mattype[ip[ipp].nm];
00410 ip[ipp].idm = new long[ip[ipp].nm];
00411 for (k = 0; k < ip[ipp].nm; k++)
00412 {
00413 ip[ipp].tm[k] = Mt->elements[i].tm[k];
00414 ip[ipp].idm[k] = Mt->elements[i].idm[k];
00415 switch (ip[ipp].tm[k])
00416 {
00417 case therisodilat:
00418 ip[ipp].hmt ^= 1;
00419 break;
00420 case nonlocplastmat:
00421 ip[ipp].hmt ^= 2;
00422 break;
00423 case nonlocdamgmat:
00424 ip[ipp].hmt ^= 2;
00425 break;
00426 case nonlocalmicroM4:
00427 ip[ipp].hmt ^= 2;
00428 break;
00429 case relaxationeuro:
00430 ip[ipp].hmt ^= 4;
00431 break;
00432 case effective_stress:
00433 ip[ipp].hmt ^= 8;
00434 break;
00435 default:
00436 break;
00437 }
00438 }
00439 ipp++;
00440 }
00441 }
00442 }
00443 }
00444 }
00445
00446
00447
00448
00449
00450
00451
00452
00453
00454
00455
00456 void mechmat::init_ip_1 (void)
00457 {
00458 long i,j,ne,ii,jj,ipp,nb,nip,ncomp;
00459 strastrestate ssst;
00460
00461 ne=Mt->ne;
00462
00463 for (i=0;i<ne;i++){
00464 nb=Mt->give_nb (i);
00465 for (ii=0;ii<nb;ii++){
00466
00467 ncomp=Mt->give_ncomp (i,ii);
00468
00469 ssst=Mt->give_ssst (i,ii);
00470 for (jj=0;jj<nb;jj++){
00471
00472 ipp=Mt->elements[i].ipp[ii][jj];
00473
00474 nip=Mt->give_nip (i,ii,jj);
00475 for (j=0;j<nip;j++){
00476 ip[ipp].ncompstr=ncomp;
00477 ip[ipp].ssst=ssst;
00478
00479 if ((ssst == planestrain) || (ssst == planestress))
00480 ip[ipp].ncompstr=4;
00481
00482 ipp++;
00483 }
00484 }
00485 }
00486 }
00487 }
00488
00489
00490
00491
00492
00493
00494
00495
00496
00497
00498
00499
00500
00501
00502
00503 void mechmat::init_ip_2 (void)
00504 {
00505 long i,n,ii,nlmid,ncompnl;
00506
00507
00508
00509
00510 if ((Mp->temperature == 1) || (Mp->pore_press == 1))
00511 {
00512
00513 n = search_reqnmq(nmqo);
00514 if (n)
00515 {
00516 for (i=0; i<n; i++)
00517 nmqid[nmqo[i]-1] = i;
00518
00519 alloc_nonmechq(n);
00520 }
00521 }
00522 else
00523 {
00524
00525 }
00526
00527
00528 for (i=0;i<tnip;i++)
00529 {
00530 ncompnl = 0;
00531 ii = ip[i].hmt & 2;
00532 if (ii > 0)
00533
00534 {
00535 nlmid = givenonlocid(i);
00536 ncompnl = give_num_averq(i, nlmid);
00537 }
00538
00539 ip[i].alloc(Mb->nlc,i,ncompnl);
00540 }
00541 }
00542
00543
00544
00545
00546
00547
00548
00549
00550
00551
00552
00553
00554
00555 void mechmat::read (XFILE *in)
00556 {
00557
00558
00559 readmatchar (in);
00560
00561
00562 tnip = intpnum ();
00563 if (Mespr==1){
00564 fprintf (stdout,"\n number of integration points %ld",tnip);
00565 }
00566
00567
00568 intpointalloc ();
00569
00570
00571
00572 intpointinit ();
00573
00574
00575
00576
00577
00578
00579
00580
00581
00582
00583 }
00584
00585
00586
00587
00588
00589
00590
00591
00592
00593
00594
00595
00596
00597
00598 void mechmat::readmatchar (XFILE *in)
00599 {
00600 long i;
00601
00602 xfscanf (in, "%k%ld", "num_mat_types", &nmt);
00603
00604 numtype = new long [nmt];
00605 mtype = new mattype [nmt];
00606
00607 if (nmt<1)
00608 print_err("wrong number of material types",__FILE__,__LINE__,__func__);
00609
00610 if (Mespr==1) fprintf (stdout,"\n number of different types of materials %ld",nmt);
00611
00612 for (i=0;i<nmt;i++)
00613 {
00614 xfscanf (in,"%k%m %k%ld", "mattype", &mattype_kwdset, &mtype[i], "num_inst", &numtype[i]);
00615 if (numtype[i]<1)
00616 print_err("wrong number of material characteristics",__FILE__,__LINE__,__func__);
00617
00618 readmattype(in, mtype[i], numtype[i]);
00619 }
00620 }
00621
00622
00623
00624
00625
00626
00627
00628
00629
00630
00631
00632
00633
00634
00635 void mechmat::readmattype(XFILE *in, mattype mtype, long numt)
00636 {
00637 long j,k;
00638
00639 switch (mtype)
00640 {
00641 case elisomat:{
00642 if (Mespr==1) fprintf (stdout,"\n number of elastic isotropic materials %ld",numt);
00643 eliso = new elastisomat [numt];
00644 for (j=0;j<numt;j++){
00645 k=numt+1;
00646 xfscanf (in,"%ld",&k);
00647 if (k>numt || k<1)
00648 print_err("wrong number of elastic isotropic material",__FILE__,__LINE__,__func__);
00649 eliso[k-1].read (in);
00650 }
00651 break;
00652 }
00653 case elgmat3d:{
00654 if (Mespr==1) fprintf (stdout,"\n number of elastic anisotropic materials %ld",numt);
00655 elgm3d = new elastgmat3d [numt];
00656 for (j=0;j<numt;j++){
00657 k=numt+1;
00658 xfscanf (in,"%ld",&k);
00659 if (k>numt || k<1)
00660 print_err("wrong number of elastic anisotropic material",__FILE__,__LINE__,__func__);
00661 elgm3d[k-1].read (in);
00662 }
00663 break;
00664 }
00665
00666 case simplas1d:{
00667 if (Mespr==1) fprintf (stdout,"\n number of simple plastic 1D materials %ld",numt);
00668 spl1d = new splas1d [numt];
00669 for (j=0;j<numt;j++){
00670 k=numt+1;
00671 xfscanf (in,"%ld",&k);
00672 if (k>numt || k<1)
00673 print_err("wrong number of simple plastic 1D material",__FILE__,__LINE__,__func__);
00674 spl1d[k-1].read (in);
00675 }
00676 break;
00677 }
00678 case jflow:{
00679 if (Mespr==1) fprintf (stdout,"\n number of J2 flow materials %ld",numt);
00680 j2f = new j2flow [numt];
00681 for (j=0;j<numt;j++){
00682 k=numt+1;
00683 xfscanf (in,"%ld",&k);
00684 if (k>numt || k<1)
00685 print_err("wrong number of J2 flow material",__FILE__,__LINE__,__func__);
00686 j2f[k-1].read (in);
00687 }
00688 break;
00689 }
00690 case jflow2:{
00691 if (Mespr==1) fprintf (stdout,"\n number of J2 flow materials %ld",numt);
00692 j2f2 = new j2flow2 [numt];
00693 for (j=0;j<numt;j++){
00694 k=numt+1;
00695 xfscanf (in,"%ld",&k);
00696 if (k>numt || k<1)
00697 print_err("wrong number of J2 flow material",__FILE__,__LINE__,__func__);
00698 j2f2[k-1].read (in);
00699 }
00700 break;
00701 }
00702 case microplaneM4:{
00703 if (Mespr==1) fprintf (stdout,"\n number of microplane materials %ld",numt);
00704 mpM4 = new microM4 [numt];
00705 for (j=0;j<numt;j++){
00706 k=numt+1;
00707 xfscanf (in,"%ld",&k);
00708 if (k>numt || k<1)
00709 print_err("wrong number of microplane material",__FILE__,__LINE__,__func__);
00710 mpM4[k-1].read (in);
00711 }
00712 break;
00713 }
00714 case nonlocalmicroM4:{
00715 if (Mespr==1) fprintf (stdout,"\n number of nonlocal microplane materials %ld",numt);
00716 nmpM4 = new nonlocmicroM4 [numt];
00717 for (j=0;j<numt;j++){
00718 k=numt+1;
00719 xfscanf (in,"%ld",&k);
00720 if (k>numt || k<1)
00721 print_err("wrong number of nonlocal microplane material",__FILE__,__LINE__,__func__);
00722 nmpM4[k-1].read (in);
00723 }
00724 break;
00725 }
00726 case microsimp:{
00727 if (Mespr==1) fprintf (stdout,"\n number of microplane materials %ld",numt);
00728 mpSIM = new microSIM [numt];
00729 for (j=0;j<numt;j++){
00730 k=numt+1;
00731 xfscanf (in,"%ld",&k);
00732 if (k>numt || k<1)
00733 print_err("wrong number of microplane material",__FILE__,__LINE__,__func__);
00734 mpSIM[k-1].read (in);
00735 }
00736 break;
00737 }
00738 case microfibro:{
00739 if (Mespr==1) fprintf (stdout,"\n number of microplane materials %ld",numt);
00740 mpfib = new microfiber [numt];
00741 for (j=0;j<numt;j++){
00742 k=numt+1;
00743 xfscanf (in,"%ld",&k);
00744 if (k>numt || k<1)
00745 print_err("wrong number of microplane material",__FILE__,__LINE__,__func__);
00746 mpfib[k-1].read (in);
00747 }
00748 break;
00749 }
00750 case mohrcoul:{
00751 if (Mespr==1) fprintf (stdout,"\n number of Mohr-Coulomb materials %ld",numt);
00752 mcoul = new mohrcoulomb [numt];
00753 for (j=0;j<numt;j++){
00754 k=numt+1;
00755 xfscanf (in,"%ld",&k);
00756 if (k>numt || k<1)
00757 print_err("wrong number of Mohr-Coulomb material",__FILE__,__LINE__,__func__);
00758 mcoul[k-1].read (in);
00759 }
00760 break;
00761 }
00762 case mohrcoulparab:{
00763 if (Mespr==1) fprintf (stdout,"\n number of parabolic Mohr-Coulomb materials %ld",numt);
00764 mcpar = new mohrcoulombpar [numt];
00765 for (j=0;j<numt;j++){
00766 k=numt+1;
00767 xfscanf (in,"%ld",&k);
00768 if (k>numt || k<1)
00769 print_err("wrong number of parabolic Mohr-Coulomb material",__FILE__,__LINE__,__func__);
00770 mcpar[k-1].read (in);
00771 }
00772 break;
00773 }
00774 case boermaterial:{
00775 if (Mespr==1) fprintf (stdout,"\n number of Boer materials %ld",numt);
00776 boerm = new boermat [numt];
00777 for (j=0;j<numt;j++){
00778 k=numt+1;
00779 xfscanf (in,"%ld",&k);
00780 if (k>numt || k<1)
00781 print_err("wrong number of Boer material",__FILE__,__LINE__,__func__);
00782 boerm[k-1].read (in);
00783 }
00784 break;
00785 }
00786 case druckerprager:{
00787 if (Mespr==1) fprintf (stdout,"\n number of Drucker-Prager materials %ld",numt);
00788 drprm = new drprag [numt];
00789 for (j=0;j<numt;j++){
00790 k=numt+1;
00791 xfscanf (in,"%ld",&k);
00792 if (k>numt || k<1)
00793 print_err("wrong number of Drucker-Prager material",__FILE__,__LINE__,__func__);
00794 drprm[k-1].read (in);
00795 }
00796 break;
00797 }
00798 case druckerprager2:{
00799 if (Mespr==1) fprintf (stdout,"\n number of Drucker-Prager2 materials %ld",numt);
00800 drprm2 = new drprag2 [numt];
00801 for (j=0;j<numt;j++){
00802 k=numt+1;
00803 xfscanf (in,"%ld",&k);
00804 if (k>numt || k<1)
00805 print_err("wrong number of Drucker-Prager material",__FILE__,__LINE__,__func__);
00806 drprm2[k-1].read (in);
00807 }
00808 break;
00809 }
00810 case modcamclaymat:{
00811 if (Mespr==1) fprintf (stdout,"\n number of modified cam-clay materials %ld",numt);
00812 cclay = new camclay [numt];
00813 for (j=0;j<numt;j++){
00814 k=numt+1;
00815 xfscanf (in,"%ld",&k);
00816 if (k>numt || k<1)
00817 print_err("wrong number of modified cam-clay material",__FILE__,__LINE__,__func__);
00818 cclay[k-1].read (in);
00819 }
00820 break;
00821 }
00822 case modcamclaycoupmat:{
00823 if (Mespr==1) fprintf (stdout,"\n number of modified cam-clay materials for coupling %ld",numt);
00824 cclayc = new camclaycoup [numt];
00825 for (j=0;j<numt;j++){
00826 k=numt+1;
00827 xfscanf (in,"%ld",&k);
00828 if (k>numt || k<1)
00829 print_err("wrong number of modified cam-clay material for coupling",__FILE__,__LINE__,__func__);
00830 cclayc[k-1].read (in);
00831 }
00832 break;
00833 }
00834 case hypoplastmat:{
00835 if (Mespr==1) fprintf (stdout,"\n number of modified hypoplast materials for coupling %ld",numt);
00836 hypopl = new hypoplast [numt];
00837 for (j=0;j<numt;j++){
00838 k=numt+1;
00839 xfscanf (in,"%ld",&k);
00840 if (k>numt || k<1)
00841 print_err("wrong number of hypoplast material for coupling",__FILE__,__LINE__,__func__);
00842 hypopl[k-1].read (in);
00843 }
00844 break;
00845 }
00846 case shefpl:{
00847 if (Mespr==1) fprintf (stdout,"\n number of Sheffield plasticity materials %ld",numt);
00848 shpl = new shefplast [numt];
00849 for (j=0;j<numt;j++){
00850 k=numt+1;
00851 xfscanf (in,"%ld",&k);
00852 if (k>numt || k<1)
00853 print_err("wrong number of Sheffiled material",__FILE__,__LINE__,__func__);
00854 shpl[k-1].read (in);
00855 }
00856 break;
00857 }
00858 case hissplasmat:{
00859 if (Mespr==1) fprintf (stdout,"\n number of HISS plasticity materials %ld",numt);
00860 hisspl = new hissplas [numt];
00861 for (j=0;j<numt;j++){
00862 k=numt+1;
00863 xfscanf (in,"%ld",&k);
00864 if (k>numt || k<1)
00865 print_err("wrong number of HISS plasticity material",__FILE__,__LINE__,__func__);
00866 hisspl[k-1].read (in);
00867 }
00868 break;
00869 }
00870 case chenplast:{
00871 if (Mespr==1) fprintf (stdout,"\n number of Chen plasticity materials %ld",numt);
00872 chplast = new chen [numt];
00873 for (j=0;j<numt;j++){
00874 k=numt+1;
00875 xfscanf (in,"%ld",&k);
00876 if (k>numt || k<1)
00877 print_err("wrong number of Chen material",__FILE__,__LINE__,__func__);
00878 chplast[k-1].read (in);
00879 }
00880 break;
00881 }
00882 case glasgowmechmat:{
00883 if (Mespr==1) fprintf (stdout,"\n number of Glasgow plasticity materials %ld",numt);
00884 glasgmat = new glasgmech [numt];
00885 for (j=0;j<numt;j++){
00886 k=numt+1;
00887 xfscanf (in,"%ld",&k);
00888 if (k>numt || k<1)
00889 print_err("wrong number of Glasgow material",__FILE__,__LINE__,__func__);
00890 glasgmat[k-1].read (in);
00891 }
00892 break;
00893 }
00894 case glasgowdamage:{
00895 if (Mespr==1) fprintf (stdout,"\n number of Glasgow damage materials %ld",numt);
00896 glasgdam = new glasgowdam [numt];
00897 for (j=0;j<numt;j++){
00898 k=numt+1;
00899 xfscanf (in,"%ld",&k);
00900 if (k>numt || k<1)
00901 print_err("wrong number of Glasgow damage material",__FILE__,__LINE__,__func__);
00902 glasgdam[k-1].read (in);
00903 }
00904 break;
00905 }
00906 case creep_damage:{
00907 if (Mespr==1) fprintf (stdout,"\n number of artificial material creep-damage %ld",numt);
00908 crdam = new creepdam [numt];
00909 for (j=0;j<numt;j++){
00910 k=numt+1;
00911 xfscanf (in,"%ld",&k);
00912 if (k>numt || k<1)
00913 print_err("wrong number of artificial material creep-damge",__FILE__,__LINE__,__func__);
00914 }
00915 break;
00916 }
00917 case time_switchmat:{
00918 if (Mespr==1) fprintf (stdout,"\n number of artificial time-switched material %ld",numt);
00919 tswmat = new timeswmat [numt];
00920 for (j=0;j<numt;j++){
00921 k=numt+1;
00922 xfscanf (in,"%ld",&k);
00923 if (k>numt || k<1)
00924 print_err("wrong number of artificial time-switched material",__FILE__,__LINE__,__func__);
00925 }
00926 break;
00927 }
00928 case effective_stress:{
00929 if (Mespr==1) fprintf (stdout,"\n number of artificial material for effective stresses %ld",numt);
00930 effstr = new effstress [numt];
00931 for (j=0;j<numt;j++){
00932 k=numt+1;
00933 xfscanf (in,"%ld",&k);
00934 if (k>numt || k<1)
00935 print_err("wrong number of artificial material for effective stresses",__FILE__,__LINE__,__func__);
00936 effstr[k-1].read (in);
00937 }
00938 break;
00939 }
00940 case simvisc:{
00941 if (Mespr==1) fprintf (stdout,"\n number of simple viscous materials %ld",numt);
00942 svis = new simviscous [numt];
00943 for (j=0;j<numt;j++){
00944 k=numt+1;
00945 xfscanf (in,"%ld",&k);
00946 if (k>numt || k<1)
00947 print_err("wrong number of simple viscous material in function",__FILE__,__LINE__,__func__);
00948 svis[k-1].read (in);
00949 }
00950 break;
00951 }
00952 case lemaitr:{
00953 if (Mespr==1) fprintf (stdout,"\n number of Lemaitre visco materials %ld",numt);
00954 lmtr = new lemaitre [numt];
00955 for (j=0;j<numt;j++){
00956 k=numt+1;
00957 xfscanf (in,"%ld",&k);
00958 if (k>numt || k<1)
00959 print_err("wrong number of Lemaitre visco material",__FILE__,__LINE__,__func__);
00960 lmtr[k-1].read (in);
00961 }
00962 break;
00963 }
00964 case scaldamage:{
00965 if (Mespr==1) fprintf (stdout,"\n number of simple damage materials %ld",numt);
00966 scdam = new scaldam [numt];
00967 for (j=0;j<numt;j++){
00968 k=numt+1;
00969 xfscanf (in,"%ld",&k);
00970 if (k>numt || k<1)
00971 print_err("wrong number of simple damage material",__FILE__,__LINE__,__func__);
00972 scdam[k-1].read (in);
00973 }
00974 break;
00975 }
00976 case scaldamagecc:{
00977 if (Mespr==1) fprintf (stdout,"\n number of simple damage materials with crack closure %ld",numt);
00978 scdamcc = new scaldamcc [numt];
00979 for (j=0;j<numt;j++){
00980 k=numt+1;
00981 xfscanf (in,"%ld",&k);
00982 if (k>numt || k<1)
00983 print_err("wrong number of simple damage material with crack closure",__FILE__,__LINE__,__func__);
00984 scdamcc[k-1].read (in);
00985 }
00986 break;
00987 }
00988 case ortodamage:{
00989 if (Mespr==1) fprintf (stdout,"\n number of ortotropic damage materials %ld",numt);
00990 ortdam = new ortodam [numt];
00991 for (j=0;j<numt;j++){
00992 k=numt+1;
00993 xfscanf (in,"%ld",&k);
00994 if (k>numt || k<1)
00995 print_err("wrong number of ortotropic damage material",__FILE__,__LINE__,__func__);
00996 ortdam[k-1].read (in);
00997 }
00998 break;
00999 }
01000 case ortodamagerot:{
01001 if (Mespr==1) fprintf (stdout,"\n number of orthotropic damage materials with crack rotation %ld",numt);
01002 ortdamrot = new ortodamrot [numt];
01003 for (j=0;j<numt;j++){
01004 k=numt+1;
01005 xfscanf (in,"%ld",&k);
01006 if (k>numt || k<1)
01007 print_err("wrong number of orthotropic damage material with crack rotation",__FILE__,__LINE__,__func__);
01008 ortdamrot[k-1].read (in);
01009 }
01010 break;
01011 }
01012 case anisodamage:{
01013 if (Mespr==1) fprintf (stdout,"\n number of anisotropic damage materials %ld",numt);
01014 anidam = new anisodam [numt];
01015 for (j=0;j<numt;j++){
01016 k=numt+1;
01017 xfscanf (in,"%ld",&k);
01018 if (k>numt || k<1)
01019 print_err("wrong number of anisotropic damage material",__FILE__,__LINE__,__func__);
01020 anidam[k-1].read (in);
01021 }
01022 break;
01023 }
01024 case anisodamagerot:{
01025 if (Mespr==1) fprintf (stdout,"\n number of rotational anisotropic damage materials with crack rotation %ld",numt);
01026 anidamrot = new anisodamrot [numt];
01027 for (j=0;j<numt;j++){
01028 k=numt+1;
01029 xfscanf (in,"%ld",&k);
01030 if (k>numt || k<1)
01031 print_err("wrong number of anisotropic damage material with crack rotation",__FILE__,__LINE__,__func__);
01032 anidamrot[k-1].read (in);
01033 }
01034 break;
01035 }
01036
01037 case aci:{
01038 if (Mespr==1) fprintf (stdout,"\n number of aci78 materials %ld",numt);
01039 aci78mod = new aci78 [numt];
01040 for (j=0;j<numt;j++){
01041 k=numt+1;
01042 xfscanf (in,"%ld",&k);
01043 if (k>numt || k<1)
01044 print_err("wrong number of aci78 material",__FILE__,__LINE__,__func__);
01045 aci78mod[k-1].read (in);
01046 }
01047 break;
01048 }
01049 case cebfip:{
01050 if (Mespr==1) fprintf (stdout,"\n number of cebfip78 materials %ld",numt);
01051 cebfip78mod = new cebfip78 [numt];
01052 for (j=0;j<numt;j++){
01053 k=numt+1;
01054 xfscanf (in,"%ld",&k);
01055 if (k>numt || k<1)
01056 print_err("wrong number of cebfip78 material in function",__FILE__,__LINE__,__func__);
01057 cebfip78mod[k-1].read (in);
01058 }
01059 break;
01060 }
01061 case graphm:{
01062 if (Mespr==1) fprintf (stdout,"\n number of graphic material materials %ld",numt);
01063 grmat = new graphmat [numt];
01064 for (j=0;j<numt;j++){
01065 k=numt+1;
01066 xfscanf (in,"%ld",&k);
01067 if (k>numt || k<1)
01068 print_err("wrong number of graphic material",__FILE__,__LINE__,__func__);
01069 grmat[k-1].read (in);
01070 }
01071 break;
01072 }
01073 case varelisomat:{
01074 if (Mespr==1) fprintf (stdout,"\n number of variable elastic isotropic materials %ld",numt);
01075 veliso = new varelastisomat [numt];
01076 for (j=0;j<numt;j++){
01077 k=numt+1;
01078 xfscanf (in,"%ld",&k);
01079 if (k>numt || k<1)
01080 print_err("wrong number of variable elastic isotropic material",__FILE__,__LINE__,__func__);
01081 veliso[k-1].read (in);
01082 }
01083 break;
01084 }
01085 case geoelast:{
01086 if (Mespr==1) fprintf (stdout,"\n number of geoelast material materials %ld",numt);
01087 geoel = new geoelastmat [numt];
01088 for (j=0;j<numt;j++){
01089 k=numt+1;
01090 xfscanf (in,"%ld",&k);
01091 if (k>numt || k<1)
01092 print_err("wrong number of graphic material",__FILE__,__LINE__,__func__);
01093 geoel[k-1].read (in);
01094 }
01095 break;
01096 }
01097 case creepb3:{
01098 if (Mespr==1) fprintf (stdout,"\n number of b3 creep material materials %ld",numt);
01099 crb3 = new b3mat [numt];
01100 for (j=0;j<numt;j++){
01101 k=numt+1;
01102 xfscanf (in,"%ld",&k);
01103 if (k>numt || k<1)
01104 print_err("wrong number of b3 creep material",__FILE__,__LINE__,__func__);
01105 crb3[k-1].read (in);
01106 }
01107 break;
01108 }
01109 case creeprs:{
01110 if (Mespr==1) fprintf (stdout,"\n number of b3 creep material materials %ld",numt);
01111 crrs = new rspecmat [numt];
01112 for (j=0;j<numt;j++){
01113 k=numt+1;
01114 xfscanf (in,"%ld",&k);
01115 if (k>numt || k<1)
01116 print_err("wrong number of b3 creep material",__FILE__,__LINE__,__func__);
01117 crrs[k-1].read (in);
01118 }
01119 break;
01120 }
01121 case creepdpl:{
01122 if (Mespr==1) fprintf (stdout,"\n number of double-power law creep material materials %ld",numt);
01123 crdpl = new dplmat [numt];
01124 for (j=0;j<numt;j++){
01125 k=numt+1;
01126 xfscanf (in,"%ld",&k);
01127 if (k>numt || k<1)
01128 print_err("wrong number of double-power law creep material",__FILE__,__LINE__,__func__);
01129 crdpl[k-1].read (in);
01130 }
01131 break;
01132 }
01133 case creepbaz:{
01134 if (Mespr==1) fprintf (stdout,"\n number of creep material materials %ld",numt);
01135 crbaz = new creepb [numt];
01136 for (j=0;j<numt;j++){
01137 k=numt+1;
01138 xfscanf (in,"%ld",&k);
01139 if (k>numt || k<1)
01140 print_err("wrong number of creep material in function",__FILE__,__LINE__,__func__);
01141 crbaz[k-1].read (in);
01142 }
01143 break;
01144 }
01145 case creepeffym:{
01146 if (Mespr==1) fprintf (stdout,"\n number of creep_effym materials %ld",numt);
01147 creffym = new creep_effym [numt];
01148 for (j=0;j<numt;j++){
01149 k=numt+1;
01150 xfscanf (in,"%ld",&k);
01151 if (k>numt || k<1)
01152 print_err ("wrong number of creep_effym material in function",__FILE__,__LINE__,__func__);
01153 }
01154 break;
01155 }
01156 case consolidation:{
01157 if (Mespr==1) fprintf (stdout,"\n number of consolidation material materials %ld",numt);
01158 csol = new consol [numt];
01159 for (j=0;j<numt;j++){
01160 k=numt+1;
01161 xfscanf (in,"%ld",&k);
01162 if (k>numt || k<1)
01163 print_err("wrong number of consolidation material model",__FILE__,__LINE__,__func__);
01164 csol[k-1].read (in);
01165 }
01166 break;
01167 }
01168 case winklerpasternak:{
01169 if (Mespr==1) fprintf (stdout,"\n number of Winkler-Pasternak material materials %ld",numt);
01170 wpast = new winpast [numt];
01171 for (j=0;j<numt;j++){
01172 k=numt+1;
01173 xfscanf (in,"%ld",&k);
01174 if (k>numt || k<1)
01175 print_err("wrong number of Winkler-Pasternak material model",__FILE__,__LINE__,__func__);
01176 wpast[k-1].read (in);
01177 }
01178 break;
01179 }
01180 case therisodilat:{
01181 if (Mespr==1) fprintf (stdout,"\n number of isotropic thermal dilatancy materials %ld",numt);
01182 tidilat = new therisomat [numt];
01183 for (j=0;j<numt;j++){
01184 k=numt+1;
01185 xfscanf (in,"%ld",&k);
01186 if (k>numt || k<1)
01187 print_err("wrong number of isotropic thermal dilatancy material",__FILE__,__LINE__,__func__);
01188 tidilat[k-1].read (in);
01189 }
01190 break;
01191 }
01192 case relaxationeuro:{
01193 if (Mespr==1) fprintf (stdout,"\n number of models of stress relaxation based on Eurocode2 %ld",numt);
01194 relaxec = new relaxeuroc [numt];
01195 for (j=0;j<numt;j++){
01196 k=numt+1;
01197 xfscanf (in,"%ld",&k);
01198 if (k>numt || k<1)
01199 print_err("wrong number of models of stress relaxation based on Eurocode2",__FILE__,__LINE__,__func__);
01200 relaxec[k-1].read (in);
01201 }
01202 break;
01203 }
01204
01205 case nonlocplastmat:{
01206 if (Mespr==1) fprintf (stdout,"\n number of nonlocal plasticity materials %ld",numt);
01207 nlplast = new nonlocplast [numt];
01208 for (j=0;j<numt;j++){
01209 k=numt+1;
01210 xfscanf (in,"%ld",&k);
01211 if (k>numt || k<1)
01212 print_err("wrong number of nonlocal plasticity material",__FILE__,__LINE__,__func__);
01213 nlplast[k-1].read (in);
01214 }
01215
01216 Mp->matmodel = nonlocal;
01217 break;
01218 }
01219 case nonlocdamgmat:{
01220 if (Mespr==1) fprintf (stdout,"\n number of nonlocal damage materials %ld",numt);
01221 nldamg = new nonlocdamg [numt];
01222 for (j=0;j<numt;j++){
01223 k=numt+1;
01224 xfscanf (in,"%ld",&k);
01225 if (k>numt || k<1)
01226 print_err("wrong number of nonlocal damage material", __FILE__,__LINE__,__func__);
01227 nldamg[k-1].read (in);
01228 }
01229
01230 Mp->matmodel = nonlocal;
01231 break;
01232 }
01233
01234 case damage_plasticity:{
01235 if (Mespr==1) fprintf (stdout,"\n number of artificial material damage-plasticity %ld",numt);
01236 dampl = new damplast [numt];
01237 for (j=0;j<numt;j++){
01238 k=numt+1;
01239 xfscanf (in,"%ld",&k);
01240 if (k>numt || k<1)
01241 print_err("wrong number of artificial material damge-plasticity",__FILE__,__LINE__,__func__);
01242 }
01243 break;
01244 }
01245
01246 case viscoplasticity:{
01247 if (Mespr==1) fprintf (stdout,"\n number of artificial material of viscoplasticity %ld",numt);
01248 visplas = new visplast [numt];
01249 for (j=0;j<numt;j++){
01250 k=numt+1;
01251 xfscanf (in,"%ld",&k);
01252 if (k>numt || k<1)
01253 print_err("wrong number of artificial material of viscoplasticity",__FILE__,__LINE__,__func__);
01254 }
01255 break;
01256 }
01257
01258 case lenjonespot:{
01259 if (Mespr==1) fprintf (stdout,"\n number of Lennard-Jones interatomic potentials %ld",numt);
01260 lenjon = new lenjonesmat [numt];
01261 for (j=0;j<numt;j++){
01262 k=numt+1;
01263 xfscanf (in,"%ld",&k);
01264 if (k>numt || k<1)
01265 print_err("wrong number of Lennard-Jones interatomic potentials",__FILE__,__LINE__,__func__);
01266 lenjon[k-1].read (in);
01267 }
01268 break;
01269 }
01270
01271 case contmat:{
01272 if (Mespr==1) fprintf (stdout,"\n number of materials for contact problems %ld",numt);
01273 conmat = new contactmat [numt];
01274 for (j=0;j<numt;j++){
01275 k=numt+1;
01276 xfscanf (in,"%ld",&k);
01277 if (k>numt || k<1)
01278 print_err("wrong number of materials for contact problems",__FILE__,__LINE__,__func__);
01279 conmat[k-1].read (in);
01280 }
01281 break;
01282 }
01283
01284 case cebfipcontmat:{
01285 if (Mespr==1) fprintf (stdout,"\n number of materials for contact problems %ld",numt);
01286 cebfipconmat = new cebfipcontactmat [numt];
01287 for (j=0;j<numt;j++){
01288 k=numt+1;
01289 xfscanf (in,"%ld",&k);
01290 if (k>numt || k<1)
01291 print_err("wrong number of materials for contact problems",__FILE__,__LINE__,__func__);
01292 cebfipconmat[k-1].read (in);
01293 }
01294 break;
01295 }
01296
01297 case shrinkagemat:{
01298 if (Mespr==1) fprintf (stdout,"\n number of materials for shrinkage description %ld",numt);
01299 shmat = new shrinkmat [numt];
01300 for (j=0;j<numt;j++){
01301 k=numt+1;
01302 xfscanf (in,"%ld",&k);
01303 if (k>numt || k<1)
01304 print_err("wrong number of materials for shrinkage description",__FILE__,__LINE__,__func__);
01305 shmat[k-1].read (in);
01306 }
01307 break;
01308 }
01309
01310 default:
01311 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
01312 }
01313 }
01314
01315
01316
01317
01318
01319
01320
01321
01322
01323
01324
01325
01326
01327 void mechmat::print (FILE *out)
01328 {
01329 long i, j;
01330
01331 fprintf(out,"\n## materials:\n");
01332 fprintf(out,"\n# number of material types\n");
01333 fprintf (out,"\n%ld\n",nmt);
01334 for (i=0;i<nmt;i++)
01335 {
01336 fprintf (out,"\n\n %d %ld", mtype[i], numtype[i]);
01337 for (j=0;j<numtype[i];j++)
01338 {
01339 fprintf (out,"\n%ld ",j+1);
01340 printmatchar (out, mtype[i], j);
01341 fprintf(out, "\n");
01342 }
01343 }
01344 }
01345
01346
01347
01348
01349
01350
01351
01352
01353
01354
01355
01356
01357
01358
01359
01360 void mechmat::printmatchar (FILE *out, mattype mt, long numinst)
01361 {
01362 switch (mt){
01363 case elisomat:{
01364 eliso[numinst].print (out);
01365 break;
01366 }
01367
01368 case scaldamage:{
01369 scdam[numinst].print (out);
01370 break;
01371 }
01372
01373 case therisodilat:{
01374 tidilat[numinst].print (out);
01375 break;
01376 }
01377
01378 case relaxationeuro:{
01379 relaxec[numinst].print (out);
01380 break;
01381 }
01382
01383 case shrinkagemat:{
01384 shmat[numinst].print (out);
01385 break;
01386 }
01387
01388 default:
01389 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
01390
01391 }
01392 }
01393
01394
01395
01396
01397
01398
01399
01400
01401
01402
01403 void mechmat::alloceigstrains ()
01404 {
01405 long i,j,nc;
01406
01407 if (eigstrains==NULL){
01408 eigstrains = new double* [tnip];
01409 for (i=0;i<tnip;i++){
01410 nc=ip[i].ncompstr;
01411 eigstrains[i] = new double [nc];
01412 for (j=0;j<nc;j++){
01413 eigstrains[i][j]=0.0;
01414 }
01415 }
01416 }
01417
01418 }
01419
01420
01421
01422
01423
01424
01425
01426
01427
01428
01429
01430 void mechmat::nulleigstrains ()
01431 {
01432 long i,j,nc;
01433
01434 if (eigstrains != NULL){
01435 for (i=0;i<tnip;i++){
01436 nc=ip[i].ncompstr;
01437 for (j=0;j<nc;j++){
01438 eigstrains[i][j]=0.0;
01439 }
01440 }
01441 }
01442 }
01443
01444
01445
01446
01447
01448
01449
01450
01451
01452
01453 void mechmat::alloctempstrains ()
01454 {
01455 long i,j,nc;
01456
01457 if (tempstrains == NULL){
01458 tempstrains = new double* [tnip];
01459 for (i=0;i<tnip;i++){
01460 nc=ip[i].ncompstr;
01461 tempstrains[i] = new double [nc];
01462 for (j=0;j<nc;j++){
01463 tempstrains[i][j]=0.0;
01464 }
01465 }
01466 }
01467
01468 }
01469
01470
01471
01472
01473
01474
01475
01476
01477
01478
01479 void mechmat::nulltempstrains ()
01480 {
01481 long i,j,nc;
01482
01483 if (tempstrains != NULL){
01484 for (i=0;i<tnip;i++){
01485 nc=ip[i].ncompstr;
01486 for (j=0;j<nc;j++){
01487 tempstrains[i][j]=0.0;
01488 }
01489 }
01490 }
01491 }
01492
01493
01494
01495
01496
01497
01498
01499
01500
01501
01502 void mechmat::alloceigstresses ()
01503 {
01504 long i,j,nc;
01505
01506 if (eigstresses == NULL){
01507 eigstresses = new double* [tnip];
01508 for (i=0;i<tnip;i++){
01509 nc=ip[i].ncompstr;
01510 eigstresses[i] = new double [nc];
01511 for (j=0;j<nc;j++){
01512 eigstresses[i][j]=0.0;
01513 }
01514 }
01515 }
01516
01517 }
01518
01519
01520
01521
01522
01523
01524
01525
01526
01527
01528
01529
01530 void mechmat::initmaterialmodels (long lcid)
01531 {
01532 long i,j,nip,ipp,nm;
01533 elemtype te;
01534
01535 for (i=0;i<Mt->ne;i++){
01536 if (Gtm->leso[i]==1){
01537 te = Mt->give_elem_type (i);
01538 nip = Mt->give_tnip (i);
01539
01540 ipp=Mt->elements[i].ipp[0][0];
01541 for(j=0; j<nip; j++){
01542 nm = ip[ipp+j].nm;
01543 initvalues(lcid, ipp+j, 0, 0);
01544 }
01545 }
01546 }
01547 }
01548
01549
01550
01551
01552
01553
01554
01555
01556
01557
01558
01559
01560
01561
01562
01563
01564 void mechmat::initvalues (long lcid, long ipp,long im,long ido)
01565 {
01566 switch (ip[ipp].tm[im])
01567 {
01568 case elisomat:
01569 case simplas1d:
01570 case jflow:
01571 case jflow2:
01572 case microplaneM4:
01573 case nonlocalmicroM4:
01574 case microsimp:
01575 case microfibro:
01576 case mohrcoul:
01577 case mohrcoulparab:
01578 case boermaterial:
01579 case druckerprager:
01580 case druckerprager2:
01581 break;
01582 case modcamclaymat:
01583 if (Mp->strainstate == 0)
01584 {
01585 compute_ipstrains(lcid);
01586 Mp->strainstate = 1;
01587 }
01588 cclay[ip[ipp].idm[im]].initval(ipp, ido);
01589 break;
01590 case modcamclaycoupmat:
01591 if (Mp->strainstate == 0)
01592 {
01593 compute_ipstrains(lcid);
01594 Mp->strainstate = 1;
01595 }
01596 cclayc[ip[ipp].idm[im]].initval(ipp, ido);
01597 break;
01598 case hypoplastmat:
01599 if (Mp->strainstate == 0)
01600 {
01601 compute_ipstrains(lcid);
01602 Mp->strainstate = 1;
01603 }
01604 hypopl[ip[ipp].idm[im]].initval(ipp, ido);
01605 break;
01606 case shefpl:
01607 case chenplast:
01608 case glasgowmechmat:
01609 case glasgowdamage:
01610 break;
01611 case creep_damage:
01612 crdam[ip[ipp].idm[im]].initvalues(lcid, ipp, im, ido);
01613 break;
01614 case shrinkagemat:
01615 shmat[ip[ipp].idm[im]].initvalues(lcid, ipp, im, ido);
01616 break;
01617 case time_switchmat:
01618 tswmat[ip[ipp].idm[im]].initvalues(lcid, ipp, im, ido);
01619 break;
01620 case effective_stress:
01621 effstr[ip[ipp].idm[im]].initvalues(ipp);
01622 initvalues(lcid,ipp, im+1, ido);
01623 break;
01624 case simvisc:
01625
01626 case lemaitr:
01627 case scaldamage:
01628 case scaldamagecc:
01629 break;
01630 case aci:
01631 case cebfip:
01632 case graphm:
01633 case varelisomat:
01634 case geoelast:
01635 case ortodamage:
01636 ortdam[ip[ipp].idm[im]].initvalues(ipp, ido);
01637 break;
01638 case ortodamagerot:
01639 ortdamrot[ip[ipp].idm[im]].initvalues(ipp, ido);
01640 break;
01641 case anisodamage:
01642 anidam[ip[ipp].idm[im]].initvalues(ipp, ido);
01643 break;
01644 case anisodamagerot:
01645 anidamrot[ip[ipp].idm[im]].initvalues(ipp, ido);
01646 break;
01647 case contmat:
01648 case cebfipcontmat:
01649 break;
01650 case creeprs:
01651 case creepb3:
01652 case creepdpl:{
01653 creep_initmaterialmodel(ipp,im,ido);
01654 }
01655 case creepbaz:
01656 case consolidation:
01657 case winklerpasternak:
01658 case therisodilat:
01659 case relaxationeuro:
01660 case nonlocplastmat:
01661 case nonlocdamgmat:
01662 case damage_plasticity:
01663 case viscoplasticity:
01664 case creepeffym:
01665 case rcmatmodnorm:
01666 break;
01667 default:
01668 print_err("unknown material type is required", __FILE__,__LINE__,__func__);
01669 }
01670 }
01671
01672
01673
01674
01675
01676
01677
01678
01679
01680
01681
01682
01683
01684
01685
01686
01687
01688
01689
01690
01691
01692
01693
01694
01695
01696 void mechmat::storestrain (long lcid,long ipp,vector &eps)
01697 {
01698 long i,j,m,n;
01699 m=lcid*ip[ipp].ncompstr;
01700 n=eps.n;
01701 j=0;
01702 for (i=m;i<m+n;i++){
01703 ip[ipp].strain[i]=eps[j];
01704 j++;
01705 }
01706 }
01707
01708
01709
01710
01711
01712
01713
01714
01715
01716
01717
01718
01719
01720
01721
01722 void mechmat::storestrain (long lcid,long ipp,long fi,vector &eps)
01723 {
01724 long i,j,m,n;
01725 m=lcid*ip[ipp].ncompstr;
01726 n=eps.n;
01727 j=0;
01728 for (i=m+fi;i<m+fi+n;i++){
01729 ip[ipp].strain[i]=eps[j];
01730 j++;
01731 }
01732 }
01733
01734
01735
01736
01737
01738
01739
01740
01741
01742
01743
01744
01745
01746
01747
01748
01749
01750 void mechmat::storestrain (long lcid,long ipp,long fi,long ncomp,vector &eps)
01751 {
01752 long i,j,m;
01753 m=lcid*ip[ipp].ncompstr;
01754 j=fi;
01755 for (i=m+fi;i<m+fi+ncomp;i++){
01756 ip[ipp].strain[i]=eps[j];
01757 j++;
01758 }
01759 }
01760
01761
01762
01763
01764
01765
01766
01767
01768
01769
01770
01771
01772
01773
01774
01775
01776
01777 void mechmat::storestrain (long lcid, long ipp, long fi, long ncomp, double *eps)
01778 {
01779 long i,j,m;
01780 m=lcid*ip[ipp].ncompstr;
01781 j=fi;
01782 for (i=m+fi;i<m+fi+ncomp;i++){
01783 ip[ipp].strain[i]=eps[j];
01784 j++;
01785 }
01786 }
01787
01788
01789
01790
01791
01792
01793
01794
01795
01796
01797
01798
01799
01800
01801 void mechmat::givestrain (long lcid,long ipp,vector &eps)
01802 {
01803 long i,j,m,n;
01804 m=lcid*ip[ipp].ncompstr;
01805 n=eps.n;
01806 j=0;
01807 for (i=m;i<n+m;i++){
01808 eps[j]=ip[ipp].strain[i];
01809 j++;
01810 }
01811 }
01812
01813
01814
01815
01816
01817
01818
01819
01820
01821
01822
01823
01824
01825
01826
01827 void mechmat::givestrain (long lcid,long ipp,long fi,vector &eps)
01828 {
01829 long i,j,m,n;
01830 m=lcid*ip[ipp].ncompstr;
01831 n=eps.n;
01832 j=0;
01833 for (i=m+fi;i<m+fi+n;i++){
01834 eps[j]=ip[ipp].strain[i];
01835 j++;
01836 }
01837 }
01838
01839
01840
01841
01842
01843
01844
01845
01846
01847
01848
01849
01850
01851
01852
01853
01854
01855
01856 void mechmat::givestrain (long lcid,long ipp,long fi,long ncomp,vector &eps)
01857 {
01858 long i,j,m;
01859 m=lcid*ip[ipp].ncompstr;
01860 j=fi;
01861 for (i=m+fi;i<m+fi+ncomp;i++){
01862 eps[j]=ip[ipp].strain[i];
01863 j++;
01864 }
01865 }
01866
01867
01868
01869
01870
01871
01872
01873
01874
01875
01876
01877 void mechmat::cleanstrain ()
01878 {
01879 long i;
01880
01881 for (i=0;i<tnip;i++){
01882
01883 ip[i].clean_strains (1);
01884 }
01885 }
01886
01887
01888
01889
01890
01891
01892
01893
01894
01895
01896
01897
01898
01899
01900 void mechmat::storestress (long lcid,long ipp,vector &sig)
01901 {
01902 long i,j,m,n;
01903 m=lcid*ip[ipp].ncompstr;
01904 n=sig.n;
01905 j=0;
01906 for (i=m;i<m+n;i++){
01907 ip[ipp].stress[i]=sig[j];
01908 j++;
01909 }
01910 }
01911
01912
01913
01914
01915
01916
01917
01918
01919
01920
01921
01922
01923
01924
01925
01926 void mechmat::storestress (long lcid,long ipp,long fi,vector &sig)
01927 {
01928 long i,j,m,n;
01929 m=lcid*ip[ipp].ncompstr;
01930 n=sig.n;
01931 j=0;
01932 for (i=m+fi;i<m+fi+n;i++){
01933 ip[ipp].stress[i]=sig[j];
01934 j++;
01935 }
01936 }
01937
01938
01939
01940
01941
01942
01943
01944
01945
01946
01947
01948
01949
01950
01951
01952
01953
01954
01955 void mechmat::storestress (long lcid,long ipp,long fi,long ncomp,vector &sig)
01956 {
01957 long i,j,m;
01958 m=lcid*ip[ipp].ncompstr;
01959 j=fi;
01960 for (i=m+fi;i<m+fi+ncomp;i++){
01961 ip[ipp].stress[i]=sig[j];
01962 j++;
01963 }
01964 }
01965
01966
01967
01968
01969
01970
01971
01972
01973
01974
01975
01976
01977
01978
01979
01980
01981
01982
01983 void mechmat::storestress (long lcid, long ipp, long fi, long ncomp, double *sig)
01984 {
01985 long i,j,m;
01986 m=lcid*ip[ipp].ncompstr;
01987 j=fi;
01988 for (i=m+fi;i<m+fi+ncomp;i++){
01989 ip[ipp].stress[i]=sig[j];
01990 j++;
01991 }
01992 }
01993
01994
01995
01996
01997
01998
01999
02000
02001
02002
02003
02004
02005
02006
02007 void mechmat::givestress (long lcid,long ipp,vector &sig)
02008 {
02009 long i,j,m,n;
02010 m=lcid*ip[ipp].ncompstr;
02011 n=sig.n;
02012 j=0;
02013 for (i=m;i<m+n;i++){
02014 sig[j]=ip[ipp].stress[i];
02015 j++;
02016 }
02017 }
02018
02019
02020
02021
02022
02023
02024
02025
02026
02027
02028
02029
02030
02031
02032
02033 void mechmat::givestress (long lcid,long ipp,long fi,vector &sig)
02034 {
02035 long i,j,m,n;
02036 m=lcid*ip[ipp].ncompstr;
02037 n=sig.n;
02038 j=0;
02039 for (i=m+fi;i<m+fi+n;i++){
02040 sig[j]=ip[ipp].stress[i];
02041 j++;
02042 }
02043 }
02044
02045
02046
02047
02048
02049
02050
02051
02052
02053
02054
02055
02056
02057
02058
02059
02060
02061
02062 void mechmat::givestress (long lcid,long ipp,long fi,long ncomp,vector &sig)
02063 {
02064 long i,j,m;
02065 m=lcid*ip[ipp].ncompstr;
02066 j=fi;
02067 for (i=m+fi;i<m+fi+ncomp;i++){
02068 sig[j]=ip[ipp].stress[i];
02069 j++;
02070 }
02071 }
02072
02073
02074
02075
02076
02077
02078
02079
02080
02081
02082
02083
02084
02085 void mechmat::storeeigstrain (long ipp,vector &eps)
02086 {
02087 long i,nc;
02088
02089 nc=ip[ipp].ncompstr;
02090 if (est==eigstrain){
02091 for (i=0;i<nc;i++){
02092 eigstrains[ipp][i]=eps[i];
02093 }
02094 }
02095 if (est==tempstrain){
02096 for (i=0;i<nc;i++){
02097 tempstrains[ipp][i]=eps[i];
02098 }
02099 }
02100 }
02101
02102
02103
02104
02105
02106
02107
02108
02109
02110
02111
02112
02113
02114
02115
02116 void mechmat::storeeigstrain (long ipp,long fi,long ncomp,vector &eps)
02117 {
02118 long i,j;
02119
02120 if (est==eigstrain){
02121 j=0;
02122 for (i=fi;i<fi+ncomp;i++){
02123 eigstrains[ipp][i] = eps[j];
02124 j++;
02125 }
02126 }
02127 if (est==tempstrain){
02128 j=0;
02129 for (i=fi;i<fi+ncomp;i++){
02130 tempstrains[ipp][i] = eps[j];
02131 j++;
02132 }
02133 }
02134 }
02135
02136
02137
02138
02139
02140
02141
02142
02143
02144
02145
02146
02147
02148
02149 void mechmat::giveeigstrain (long ipp,vector &eps)
02150 {
02151 long i,n;
02152
02153 n=eps.n;
02154 if (est==eigstrain){
02155 for (i=0;i<n;i++){
02156 eps[i] = eigstrains[ipp][i];
02157 }
02158 }
02159 if (est==tempstrain){
02160 for (i=0;i<n;i++){
02161 eps[i] = tempstrains[ipp][i];
02162 }
02163 }
02164 }
02165
02166
02167
02168
02169
02170
02171
02172
02173
02174
02175
02176
02177
02178
02179
02180
02181 void mechmat::giveeigstrain (long ipp,long fi,long ncomp,vector &eps)
02182 {
02183 long i,j;
02184
02185 if (est==eigstrain){
02186 j=0;
02187 for (i=fi;i<fi+ncomp;i++){
02188 eps[j] = eigstrains[ipp][i];
02189 j++;
02190 }
02191 }
02192 if (est==tempstrain){
02193 j=0;
02194 for (i=fi;i<fi+ncomp;i++){
02195 eps[j] = tempstrains[ipp][i];
02196 j++;
02197 }
02198 }
02199 }
02200
02201
02202
02203
02204
02205
02206
02207
02208
02209
02210
02211
02212
02213 void mechmat::storeeigstress (long ipp,vector &sig)
02214 {
02215 long i,nc;
02216
02217 nc=sig.n;
02218 for (i=0;i<nc;i++){
02219 eigstresses[ipp][i]=sig[i];
02220 }
02221
02222 }
02223
02224
02225
02226
02227
02228
02229
02230
02231
02232
02233
02234
02235
02236
02237
02238
02239 void mechmat::storeeigstress (long ipp,long fi,long ncomp,vector &sig)
02240 {
02241 long i,j;
02242
02243 j=0;
02244 for (i=fi;i<fi+ncomp;i++){
02245 eigstresses[ipp][i] = sig[j];
02246 j++;
02247 }
02248 }
02249
02250
02251
02252
02253
02254
02255
02256
02257
02258
02259
02260
02261
02262 void mechmat::giveeigstress (long ipp,vector &sig)
02263 {
02264 long i,n;
02265
02266 n=sig.n;
02267 for (i=0;i<n;i++){
02268 sig[i] = eigstresses[ipp][i];
02269 }
02270 }
02271
02272
02273
02274
02275
02276
02277
02278
02279
02280
02281
02282
02283
02284
02285 void mechmat::giveeigstress (long ipp,long fi,vector &sig)
02286 {
02287 long i,j,ncomp;
02288 ncomp=sig.n;
02289 j=0;
02290 for (i=fi;i<fi+ncomp;i++){
02291 sig[j] = eigstresses[ipp][i];
02292 j++;
02293 }
02294 }
02295
02296
02297
02298
02299
02300
02301
02302
02303
02304
02305
02306
02307
02308
02309
02310 void mechmat::giveeigstress (long ipp,long fi,long ncomp,vector &sig)
02311 {
02312 long i,j;
02313 j=fi;
02314 for (i=fi;i<fi+ncomp;i++){
02315 sig[j] = eigstresses[ipp][i];
02316 j++;
02317 }
02318 }
02319
02320
02321
02322
02323
02324
02325
02326
02327
02328
02329
02330
02331
02332
02333
02334 void mechmat::storeeqother (long ipp,long fi,long ncomp,double *comp)
02335 {
02336 long i,j;
02337 j=0;
02338 for (i=fi;i<fi+ncomp;i++){
02339 ip[ipp].eqother[i]=comp[j];
02340 j++;
02341 }
02342 }
02343
02344
02345
02346
02347
02348
02349
02350
02351
02352
02353
02354
02355
02356
02357
02358 void mechmat::giveeqother (long ipp,long fi,long ncomp,double *comp)
02359 {
02360 long i,j;
02361 j=0;
02362 for (i=fi;i<fi+ncomp;i++){
02363 comp[j]=ip[ipp].eqother[i];
02364 j++;
02365 }
02366 }
02367
02368
02369
02370
02371
02372
02373
02374
02375
02376
02377
02378
02379
02380
02381
02382 void mechmat::storeother (long ipp,long fi,long ncomp,double *comp)
02383 {
02384 long i,j;
02385 j=0;
02386 for (i=fi;i<fi+ncomp;i++){
02387 ip[ipp].other[i]=comp[j];
02388 j++;
02389 }
02390 }
02391
02392
02393
02394
02395
02396
02397
02398
02399
02400
02401
02402
02403
02404
02405
02406 long mechmat::givencompother (long ipp,long im)
02407 {
02408 long ncompo = 0;
02409 long ncompstr = ip[ipp].ncompstr;
02410
02411 switch (ip[ipp].tm[im])
02412 {
02413 case elisomat:
02414 case elgmat3d:
02415 case winklerpasternak:
02416 case graphm:
02417 case varelisomat:
02418 case geoelast:
02419 case simplas1d:
02420 case jflow:
02421 case jflow2:
02422 case microplaneM4:
02423 case microsimp:
02424 case microfibro:
02425 case mohrcoul:
02426 case mohrcoulparab:
02427 case boermaterial:
02428 case druckerprager:
02429 case druckerprager2:
02430 case chenplast:
02431 case modcamclaymat:
02432 case modcamclaycoupmat:
02433 case hypoplastmat:
02434 case shefpl:
02435 case scaldamage:
02436 case scaldamagecc:
02437 case ortodamage:
02438 case ortodamagerot:
02439 case anisodamage:
02440 case anisodamagerot:
02441 case creepbaz:
02442 case creepeffym:
02443 case contmat:
02444 case cebfipcontmat:
02445 case consolidation:
02446 case nonlocplastmat:
02447 case nonlocdamgmat:
02448 case nonlocalmicroM4:
02449 case glasgowdamage:
02450 case damage_plasticity:
02451 case therisodilat:
02452 case relaxationeuro:
02453 case simvisc:
02454 case lenjonespot:{
02455 ncompo = givencompeqother (ipp,im);
02456 break;
02457 }
02458 case creeprs:
02459 case creepb3:
02460 case creepdpl:{
02461 ncompo=creep_ncompo(ipp,im);
02462 break;
02463 }
02464 case creep_damage:
02465 ncompo += givencompother (ipp, im+1);
02466 ncompo += givencompother (ipp, im+2);
02467 break;
02468 case shrinkagemat:
02469 ncompo += 1 + 2*ncompstr;
02470 ncompo += givencompother (ipp, im+1);
02471 break;
02472 case effective_stress:
02473 ncompo += givencompother (ipp, im+1);
02474 break;
02475
02476 case lemaitr:
02477 case glasgowmechmat:
02478 case viscoplasticity:
02479 ncompo=0;
02480 break;
02481 case time_switchmat:
02482 ncompo = tswmat[ip[ipp].idm[im]].givencompother (ipp,im);
02483 break;
02484
02485 default:
02486 print_err("unknown material type is required", __FILE__, __LINE__, __func__);
02487 }
02488
02489 return ncompo;
02490 }
02491
02492
02493
02494
02495
02496
02497
02498
02499
02500
02501
02502
02503
02504
02505
02506 long mechmat::givencompeqother (long ipp,long im)
02507 {
02508 long ncompo = 0;
02509 long ncompstr = ip[ipp].ncompstr;
02510 long n;
02511
02512 switch (ip[ipp].tm[im])
02513 {
02514 case elisomat:
02515 break;
02516 case elgmat3d:
02517 break;
02518
02519
02520
02521 case winklerpasternak:
02522 break;
02523 case graphm:
02524 break;
02525 case varelisomat:
02526 break;
02527 case geoelast:
02528 ncompo = 1;
02529 ncompo += givencompeqother(ipp, im+1);
02530 break;
02531 case simplas1d:
02532 ncompo=3;
02533 break;
02534 case jflow:
02535 ncompo=ncompstr+2;
02536 break;
02537 case jflow2:
02538 ncompo=ncompstr+2;
02539 break;
02540 case microplaneM4:
02541 n=mpM4[ip[ipp].idm[im]].numberOfMicroplanes;
02542 ncompo=3*n+1+ncompstr;
02543 break;
02544 case microsimp:
02545 n=mpSIM[ip[ipp].idm[im]].numberOfMicroplanes;
02546 ncompo=n+1+6;
02547 break;
02548 case microfibro:
02549 n=mpfib[ip[ipp].idm[im]].numberOfMicroplanes;
02550 ncompo=5*n+1+6;
02551 break;
02552 case mohrcoul:
02553 ncompo=ncompstr+1;
02554 break;
02555 case mohrcoulparab:
02556 ncompo=ncompstr+1;
02557 break;
02558 case boermaterial:
02559 ncompo=ncompstr+1;
02560 break;
02561 case druckerprager:
02562 ncompo=ncompstr+4;
02563 break;
02564 case druckerprager2:
02565 ncompo=ncompstr+4;
02566 break;
02567 case chenplast:
02568 ncompo=ncompstr+3;
02569 break;
02570 case modcamclaymat:
02571 ncompo=ncompstr+2+3+4;
02572 break;
02573 case modcamclaycoupmat:
02574 ncompo=ncompstr+3+2+4;
02575 case hypoplastmat:
02576 ncompo=2*ncompstr+hypopl[ip[ipp].idm[im]].nstatv+1;
02577 break;
02578 case shefpl:
02579 ncompo=ncompstr+2;
02580 break;
02581 case glasgowmechmat:
02582 ncompo=ncompstr*4+5;
02583 break;
02584 case glasgowdamage:
02585 ncompo=2;
02586 break;
02587 case creep_damage:
02588 ncompo += givencompeqother (ipp, im+1);
02589 ncompo += givencompeqother (ipp, im+2);
02590 break;
02591 case shrinkagemat:
02592 ncompo += 1 + 2*ncompstr;
02593 ncompo += givencompeqother (ipp, im+1);
02594 break;
02595 case time_switchmat:
02596 ncompo = tswmat[ip[ipp].idm[im]].givencompeqother (ipp,im);
02597 break;
02598 case effective_stress:
02599 ncompo += givencompeqother (ipp, im+1);
02600 break;
02601 case simvisc:
02602 ncompo=3*ncompstr;
02603 break;
02604
02605
02606
02607 case lemaitr:
02608 ncompo=3*ncompstr+2;
02609 break;
02610 case scaldamage:
02611 ncompo=3;
02612 break;
02613 case scaldamagecc:
02614 ncompo=2;
02615 break;
02616 case ortodamage:
02617 ncompo=19;
02618 break;
02619 case ortodamagerot:
02620 ncompo=23;
02621 break;
02622 case anisodamage:
02623
02624 ncompo=26;
02625 break;
02626 case anisodamagerot:
02627
02628 ncompo=22;
02629 break;
02630 case cebfipcontmat:
02631 case contmat:
02632 ncompo = 2;
02633 break;
02634 case creeprs:
02635 case creepb3:
02636 case creepdpl:{
02637 ncompo=creep_ncompo(ipp,im);
02638 break;
02639 }
02640 case creepbaz:
02641
02642 n=crbaz[ip[ipp].idm[im]].numberOfCreepb();
02643 ncompo=ncompstr*(n+3)+3+1;
02644 break;
02645 case creepeffym:
02646 ncompo=ncompstr;
02647 break;
02648 case consolidation:
02649 n=csol[ip[ipp].idm[im]].numberOfConsol();
02650 if(ncompstr==6) ncompo=ncompstr+2*n+1;
02651 if(ncompstr==5) ncompo=ncompstr+3*n+1;
02652 if(ncompstr==3) ncompo=ncompstr+3*n+1;
02653
02654
02655 break;
02656
02657 case nonlocplastmat:
02658 ncompo = givencompeqother (ipp, im+1);
02659 break;
02660 case nonlocdamgmat:
02661 ncompo = givencompeqother (ipp, im+1);
02662 break;
02663 case nonlocalmicroM4:
02664 n=mpM4[ip[ipp].idm[im]].numberOfMicroplanes;
02665 ncompo=3*n+1+ncompstr+ncompstr;
02666 break;
02667 case damage_plasticity:
02668 ncompo = givencompeqother (ipp, im+1);
02669 ncompo += givencompeqother (ipp, im+2);
02670 break;
02671 case viscoplasticity:
02672 ncompo = givencompeqother (ipp, im+1);
02673 ncompo += givencompeqother (ipp, im+2);
02674 break;
02675 case therisodilat:{
02676 break;
02677 }
02678 case relaxationeuro:{
02679 break;
02680 }
02681 case lenjonespot:{
02682 ncompo = 0;
02683 break;
02684 }
02685 default:{
02686 print_err("unknown material type is required", __FILE__, __LINE__, __func__);
02687 }
02688 }
02689
02690 if (((ip[ipp].hmt & 1) && (im < ip[ipp].nm-1)) && (im == 0))
02691
02692
02693
02694 {
02695 switch (ip[ipp].tm[ip[ipp].nm-1])
02696 {
02697 case therisodilat:
02698 break;
02699 default:{
02700 print_err("unknown thermal material type is required", __FILE__, __LINE__, __func__);
02701 }
02702 }
02703 }
02704 return ncompo;
02705 }
02706
02707
02708
02709
02710
02711
02712
02713
02714
02715
02716
02717
02718
02719
02720
02721 void mechmat::storenonloc (long ipp,long fi,long ncomp,double *comp)
02722 {
02723 long i,j;
02724 j=0;
02725 for (i=fi;i<fi+ncomp;i++){
02726 ip[ipp].nonloc[i]=comp[j];
02727 j++;
02728 }
02729 }
02730
02731
02732
02733
02734
02735
02736
02737
02738
02739
02740
02741
02742
02743
02744
02745 void mechmat::givenonloc (long ipp,long fi,long ncomp,double *comp)
02746 {
02747 long i,j;
02748 j=0;
02749 for (i=fi;i<fi+ncomp;i++){
02750 comp[j]=ip[ipp].nonloc[i];
02751 j++;
02752 }
02753 }
02754
02755
02756
02757
02758
02759
02760
02761
02762
02763
02764
02765
02766
02767
02768
02769
02770
02771
02772
02773
02774
02775
02776
02777
02778
02779
02780
02781
02782
02783 void mechmat::givequantity (integratedquant iq,long lcid,long ipp,long fi,vector &ipv)
02784 {
02785 switch (iq){
02786 case locstress:{
02787
02788 givestress (lcid,ipp,fi,ipv);
02789 break;
02790 }
02791 case nonlocstress:{
02792
02793 givestress (lcid,ipp,fi,ipv);
02794 break;
02795 }
02796 case stressincr:{
02797
02798 givestressincr (lcid,ipp,0,0,fi,ipv);
02799 break;
02800 }
02801 case eigstress:{
02802
02803 giveeigstress (ipp,fi,ipv);
02804 break;
02805 }
02806 default:{
02807 print_err("unknown type of quantity is required",__FILE__,__LINE__,__func__);
02808 }
02809 }
02810 }
02811
02812
02813
02814
02815
02816
02817
02818
02819
02820
02821
02822
02823
02824
02825
02826 double mechmat::give_actual_ym(long ipp, long im, long ido)
02827 {
02828 long ncompo;
02829 double e = 0.0;
02830
02831 switch (ip[ipp].tm[im])
02832 {
02833 case elisomat:
02834 e = eliso[ip[ipp].idm[im]].e;
02835 break;
02836 case modcamclaymat:
02837 e = cclay[ip[ipp].idm[im]].give_actual_ym(ipp, im, ido);
02838 break;
02839 case winklerpasternak:
02840 case simplas1d:
02841 case jflow:
02842 case jflow2:
02843 case microplaneM4:
02844 case microsimp:
02845 case microfibro:
02846 case mohrcoul:
02847 case mohrcoulparab:
02848 case boermaterial:
02849 case druckerprager:
02850 case druckerprager2:
02851 case chenplast:
02852 case shefpl:
02853 case glasgowmechmat:
02854 case glasgowdamage:
02855 case simvisc:
02856
02857 case lemaitr:
02858 case scaldamage:
02859 case scaldamagecc:
02860 case ortodamage:
02861 case ortodamagerot:
02862 case anisodamage:
02863 case anisodamagerot:
02864 case creepdpl:
02865 case creepbaz:
02866 case creepeffym:
02867 case consolidation:
02868 ncompo = givencompeqother(ipp, im);
02869 e = give_actual_ym(ipp, im+1, ido+ncompo);
02870 break;
02871 case damage_plasticity:
02872 ncompo = givencompeqother(ipp, im+1);
02873 e= give_actual_ym(ipp, im+2, ido+ncompo);
02874 break;
02875 case viscoplasticity:
02876 ncompo = givencompeqother(ipp, im+1);
02877 e = give_actual_ym(ipp, im+2, ido+ncompo);
02878 break;
02879 case nonlocplastmat:
02880 case nonlocdamgmat:
02881 case nonlocalmicroM4:
02882 case creep_damage:
02883 e = give_actual_ym(ipp, im+1, ido);
02884 break;
02885 case shrinkagemat:
02886 ncompo = givencompeqother(ipp, im);
02887 e = give_actual_ym (ipp, im+1, ido+ncompo);
02888 break;
02889 case creepb3:
02890 case creeprs:
02891 e = creep_give_actual_ym(ipp, im, ido);
02892 break;
02893 case time_switchmat:
02894 e = tswmat[ip[ipp].idm[im]].give_actual_ym(ipp, im, ido);
02895 break;
02896 case effective_stress:{
02897 ncompo = givencompeqother(ipp, im);
02898 e = give_actual_ym(ipp, im+1, ido+ncompo);
02899 break;
02900 }
02901 default:
02902 print_err("unknown elastic material type is required",__FILE__,__LINE__,__func__);
02903 }
02904
02905 return e;
02906 }
02907
02908
02909
02910
02911
02912
02913
02914
02915
02916
02917
02918
02919
02920
02921
02922 double mechmat::give_initial_ym(long ipp, long im, long ido)
02923 {
02924 long ncompo;
02925 double e = 0.0;
02926
02927 switch (ip[ipp].tm[im])
02928 {
02929 case elisomat:
02930 e = eliso[ip[ipp].idm[im]].e;
02931 break;
02932 case varelisomat:
02933 e = veliso[ip[ipp].idm[im]].e;
02934 break;
02935 case winklerpasternak:
02936 case simplas1d:
02937 case jflow:
02938 case jflow2:
02939 case microplaneM4:
02940 case microsimp:
02941 case microfibro:
02942 case mohrcoul:
02943 case mohrcoulparab:
02944 case boermaterial:
02945 case druckerprager:
02946 case druckerprager2:
02947 case chenplast:
02948 case modcamclaymat:
02949 case shefpl:
02950 case glasgowmechmat:
02951 case glasgowdamage:
02952 case simvisc:
02953
02954 case lemaitr:
02955 case scaldamage:
02956 case scaldamagecc:
02957 case ortodamage:
02958 case ortodamagerot:
02959 case anisodamage:
02960 case anisodamagerot:
02961 case creepdpl:
02962 case creepbaz:
02963 case creepeffym:
02964 case consolidation:
02965 ncompo = givencompeqother(ipp, im);
02966 e = give_initial_ym(ipp, im+1, ido+ncompo);
02967 break;
02968 case damage_plasticity:
02969 ncompo = givencompeqother(ipp, im+1);
02970 e= give_initial_ym(ipp, im+2, ido+ncompo);
02971 break;
02972 case viscoplasticity:
02973 ncompo = givencompeqother(ipp, im+1);
02974 e = give_initial_ym(ipp, im+2, ido+ncompo);
02975 break;
02976 case nonlocplastmat:
02977 case nonlocdamgmat:
02978 case nonlocalmicroM4:
02979 case creep_damage:
02980 case effective_stress:
02981 e = give_initial_ym(ipp, im+1, ido);
02982 break;
02983 case shrinkagemat:
02984 ncompo = givencompeqother(ipp, im);
02985 e = give_initial_ym (ipp, im+1, ido+ncompo);
02986 break;
02987 case creepb3:
02988 case creeprs:
02989 e = creep_compute_inital_ym(ipp, im, ido);
02990 break;
02991 case time_switchmat:
02992 e = tswmat[ip[ipp].idm[im]].give_initial_ym(ipp, im, ido);
02993 break;
02994 default:
02995 print_err("unknown elastic material type is required",__FILE__,__LINE__,__func__);
02996 }
02997 return e;
02998 }
02999
03000
03001
03002
03003
03004
03005
03006
03007
03008
03009
03010
03011
03012
03013
03014 double mechmat::give_actual_nu(long ipp, long im, long ido)
03015 {
03016 long ncompo;
03017 double nu = 0.0;
03018
03019 switch (ip[ipp].tm[im])
03020 {
03021 case elisomat:
03022 nu = eliso[ip[ipp].idm[im]].nu;
03023 break;
03024 case modcamclaymat:
03025 nu = cclay[ip[ipp].idm[im]].give_actual_nu(ipp, im, ido);
03026 break;
03027 case winklerpasternak:
03028 case simplas1d:
03029 case jflow:
03030 case jflow2:
03031 case microplaneM4:
03032 case microsimp:
03033 case microfibro:
03034 case mohrcoul:
03035 case mohrcoulparab:
03036 case boermaterial:
03037 case druckerprager:
03038 case druckerprager2:
03039 case chenplast:
03040 case shefpl:
03041 case glasgowmechmat:
03042 case glasgowdamage:
03043 case simvisc:
03044
03045 case lemaitr:
03046 case scaldamage:
03047 case scaldamagecc:
03048 case ortodamage:
03049 case ortodamagerot:
03050 case anisodamage:
03051 case anisodamagerot:
03052 case creepdpl:
03053 case creepbaz:
03054 case creepb3:
03055 case creeprs:
03056 case creepeffym:
03057 case consolidation:
03058 ncompo = givencompeqother(ipp, im);
03059 nu = give_actual_nu(ipp, im+1, ido+ncompo);
03060 break;
03061 case damage_plasticity:
03062 ncompo = givencompeqother(ipp, im+1);
03063 nu = give_actual_nu(ipp, im+2, ido+ncompo);
03064 break;
03065 case viscoplasticity:
03066 ncompo = givencompeqother(ipp, im+1);
03067 nu = give_actual_nu(ipp, im+2, ido+ncompo);
03068 break;
03069 case creep_damage:
03070 ncompo = givencompeqother(ipp, im+1);
03071 nu = give_actual_nu(ipp, im+2, ido+ncompo);
03072 break;
03073 case shrinkagemat:
03074 ncompo = givencompeqother(ipp, im);
03075 nu = give_actual_nu(ipp, im+1, ido+ncompo);
03076 break;
03077 case nonlocplastmat:
03078 case nonlocdamgmat:
03079 case nonlocalmicroM4:
03080 case effective_stress:
03081 nu = give_actual_nu(ipp, im+1, ido);
03082 break;
03083 case time_switchmat:
03084 nu = tswmat[ip[ipp].idm[im]].give_actual_nu(ipp, im, ido);
03085 break;
03086 default:
03087 print_err("unknown elastic material type is required",__FILE__,__LINE__,__func__);
03088 }
03089 return nu;
03090 }
03091
03092
03093
03094
03095
03096
03097
03098
03099
03100
03101
03102
03103
03104
03105
03106 double mechmat::give_initial_nu(long ipp, long im, long ido)
03107 {
03108 long ncompo;
03109 double nu = 0.0;
03110
03111 switch (ip[ipp].tm[im])
03112 {
03113 case elisomat:
03114 nu = eliso[ip[ipp].idm[im]].nu;
03115 break;
03116 case varelisomat:
03117 nu = veliso[ip[ipp].idm[im]].nu;
03118 break;
03119 case winklerpasternak:
03120 case simplas1d:
03121 case jflow:
03122 case jflow2:
03123 case microplaneM4:
03124 case microsimp:
03125 case microfibro:
03126 case mohrcoul:
03127 case mohrcoulparab:
03128 case boermaterial:
03129 case druckerprager:
03130 case druckerprager2:
03131 case chenplast:
03132 case modcamclaymat:
03133 case shefpl:
03134 case glasgowmechmat:
03135 case glasgowdamage:
03136 case simvisc:
03137
03138 case lemaitr:
03139 case scaldamage:
03140 case scaldamagecc:
03141 case ortodamage:
03142 case ortodamagerot:
03143 case anisodamage:
03144 case anisodamagerot:
03145 case creepdpl:
03146 case creepbaz:
03147 case creepb3:
03148 case creeprs:
03149 case creepeffym:
03150 case consolidation:
03151 ncompo = givencompeqother(ipp, im);
03152 nu = give_actual_nu(ipp, im+1, ido+ncompo);
03153 break;
03154 case damage_plasticity:
03155 ncompo = givencompeqother(ipp, im+1);
03156 nu = give_actual_nu(ipp, im+2, ido+ncompo);
03157 break;
03158 case viscoplasticity:
03159 ncompo = givencompeqother(ipp, im+1);
03160 nu = give_actual_nu(ipp, im+2, ido+ncompo);
03161 break;
03162 case creep_damage:
03163 ncompo = givencompeqother(ipp, im+1);
03164 nu = give_actual_nu(ipp, im+2, ido+ncompo);
03165 break;
03166 case shrinkagemat:
03167 ncompo = givencompeqother(ipp, im);
03168 nu = give_actual_nu(ipp, im+1, ido+ncompo);
03169 break;
03170 case nonlocplastmat:
03171 case nonlocdamgmat:
03172 case nonlocalmicroM4:
03173 case effective_stress:
03174 nu = give_actual_nu(ipp, im+1, ido);
03175 break;
03176 case time_switchmat:
03177 nu = tswmat[ip[ipp].idm[im]].give_actual_nu(ipp, im, ido);
03178 break;
03179 default:
03180 print_err("unknown elastic material type is required",__FILE__,__LINE__,__func__);
03181 }
03182 return nu;
03183 }
03184
03185
03186
03187
03188
03189
03190
03191
03192
03193
03194
03195
03196
03197
03198
03199 double mechmat::give_actual_ft(long ipp, long im, long ido)
03200 {
03201 double ft = 0.0;
03202
03203 switch (ip[ipp].tm[im])
03204 {
03205 case scaldamage:
03206 ft = scdam[ip[ipp].idm[im]].give_actual_ft(ipp, im, ido);
03207 break;
03208 case ortodamage:
03209 ft = ortdam[ip[ipp].idm[im]].give_actual_ft(ipp);
03210 break;
03211 case ortodamagerot:
03212 ft = ortdamrot[ip[ipp].idm[im]].give_actual_ft(ipp);
03213 break;
03214 case anisodamage:
03215 ft = anidam[ip[ipp].idm[im]].give_actual_ft(ipp);
03216 break;
03217 case anisodamagerot:
03218 ft = anidamrot[ip[ipp].idm[im]].give_actual_ft(ipp);
03219 break;
03220 case nonlocdamgmat:
03221 ft = give_actual_ft(ipp,im+1,ido);
03222 break;
03223 case damage_plasticity:
03224 ft = dampl[ip[ipp].idm[im]].give_actual_ft(ipp, im, ido);
03225 break;
03226 case creep_damage:
03227 ft = crdam[ip[ipp].idm[im]].give_actual_ft(ipp, im, ido);
03228 break;
03229 case shrinkagemat:
03230 ft = shmat[ip[ipp].idm[im]].give_actual_ft(ipp, im, ido);
03231 break;
03232 case effective_stress:
03233 ft = give_actual_ft(ipp,im+1,ido);
03234 break;
03235 case creepb3:
03236 case creeprs:
03237 ft = creep_give_actual_ft(ipp, im, ido);
03238 break;
03239 case creepeffym:
03240 ft = creffym[ip[ipp].idm[im]].give_actual_ft(ipp, im, ido);
03241 break;
03242 case time_switchmat:
03243 ft = tswmat[ip[ipp].idm[im]].give_actual_ft(ipp, im, ido);
03244 break;
03245 default:
03246 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
03247 }
03248 return ft;
03249 }
03250
03251
03252
03253
03254
03255
03256
03257
03258
03259
03260
03261
03262
03263
03264
03265 double mechmat::give_actual_fc(long ipp, long im, long ido)
03266 {
03267 double fc = 0.0;
03268
03269 switch (ip[ipp].tm[im])
03270 {
03271 case ortodamage:
03272 fc = ortdam[ip[ipp].idm[im]].give_actual_fc(ipp);
03273 break;
03274 case ortodamagerot:
03275 fc = ortdamrot[ip[ipp].idm[im]].give_actual_fc(ipp);
03276 break;
03277 case anisodamage:
03278 fc = anidam[ip[ipp].idm[im]].give_actual_fc(ipp);
03279 break;
03280 case anisodamagerot:
03281 fc = anidamrot[ip[ipp].idm[im]].give_actual_fc(ipp);
03282 break;
03283 case nonlocdamgmat:
03284 fc = give_actual_fc(ipp,im+1,ido);
03285 break;
03286 case creep_damage:
03287 fc = crdam[ip[ipp].idm[im]].give_actual_fc(ipp, im, ido);
03288 break;
03289 case shrinkagemat:
03290 fc = shmat[ip[ipp].idm[im]].give_actual_fc(ipp, im, ido);
03291 break;
03292 default:
03293 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
03294 }
03295 return fc;
03296 }
03297
03298
03299
03300
03301
03302
03303
03304
03305
03306
03307
03308
03309
03310
03311 double mechmat::givemechq(nontransquant qt, long ipp)
03312 {
03313 double ret=0.0;
03314
03315 switch (qt){
03316 case precons_press:{
03317
03318 ret = give_preconspress(ipp);
03319 break;
03320 }
03321 case mean_stress_eff:{
03322 vector sig(Mm->ip[ipp].ncompstr);
03323 matrix sigt(3,3);
03324 givestress(0,ipp,sig);
03325 vector_tensor (sig,sigt,Mm->ip[ipp].ssst,stress);
03326
03327 ret = first_invar(sigt)/3.0;
03328 if (Mm->givestatusnmq(pore_pressure))
03329 {
03330
03331 ret -= Mm->givenonmechq(pore_pressure, ipp);
03332 }
03333 break;
03334 }
03335 case virgin_porosity:{
03336
03337 ret = give_virgporosity(ipp);
03338 break;
03339 }
03340 case init_porosity:{
03341
03342 ret = give_iniporosity(ipp);
03343 break;
03344 }
03345 case scal_iso_damage:{
03346
03347 ret = give_dampar (ipp);
03348 break;
03349 }
03350 case proc_zone_length:{
03351
03352 ret = give_proczonelength (ipp);
03353 break;
03354 }
03355 case crack_width:{
03356
03357 ret = give_crackwidth (ipp);
03358 break;
03359 }
03360 default:{
03361 print_err("unknown type of quantity is required",__FILE__,__LINE__,__func__);
03362 }
03363 }
03364
03365 return ret;
03366 }
03367
03368
03369
03370
03371
03372
03373
03374
03375
03376
03377
03378
03379
03380
03381
03382 double mechmat::give_preconspress(long ipp, long im, long ido)
03383 {
03384 long i;
03385 double pc;
03386
03387 i = ip[ipp].idm[im];
03388
03389 switch (ip[ipp].tm[im]){
03390 case modcamclaymat:{
03391 pc = cclay[i].give_preconspress(ipp,ido);
03392 break;
03393 }
03394 case modcamclaycoupmat:{
03395 pc = cclayc[i].give_preconspress(ipp,ido);
03396 break;
03397 }
03398 case effective_stress:
03399 pc = give_preconspress(ipp, im+1, ido);
03400 break;
03401 case damage_plasticity:
03402 pc = give_preconspress(ipp, im+1, ido);
03403 break;
03404 default:{
03405 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
03406 }
03407 }
03408
03409 return pc;
03410 }
03411
03412
03413
03414
03415
03416
03417
03418
03419
03420
03421
03422
03423
03424
03425
03426 double mechmat::give_virgporosity(long ipp, long im, long ido)
03427 {
03428 long i;
03429 double e_lambda1;
03430
03431 i = ip[ipp].idm[im];
03432
03433 switch (ip[ipp].tm[im]){
03434 case modcamclaymat:{
03435 e_lambda1 = cclay[i].give_virgporosity(ipp,ido);
03436 break;
03437 }
03438 case modcamclaycoupmat:{
03439 e_lambda1 = cclayc[i].give_virgporosity(ipp,ido);
03440 break;
03441 }
03442 case hypoplastmat:{
03443 e_lambda1 = hypopl[i].give_virgporosity(ipp,ido);
03444 break;
03445 }
03446 case effective_stress:
03447 e_lambda1 = give_virgporosity(ipp, im+1, ido);
03448 break;
03449 case damage_plasticity:
03450 e_lambda1 = give_virgporosity(ipp, im+1, ido);
03451 break;
03452 default:{
03453 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
03454 }
03455 }
03456
03457 return e_lambda1;
03458 }
03459
03460
03461
03462
03463
03464
03465
03466
03467
03468
03469
03470
03471
03472
03473
03474 double mechmat::give_iniporosity(long ipp, long im, long ido)
03475 {
03476 long i;
03477 double e_ini;
03478
03479 i = ip[ipp].idm[im];
03480
03481 switch (ip[ipp].tm[im]){
03482 case modcamclaymat:{
03483 e_ini = cclay[i].give_iniporosity(ipp,ido);
03484 break;
03485 }
03486 case modcamclaycoupmat:{
03487 e_ini = cclayc[i].give_iniporosity(ipp,ido);
03488 break;
03489 }
03490 case hypoplastmat:{
03491 e_ini = hypopl[i].give_iniporosity(ipp,ido);
03492 break;
03493 }
03494 case effective_stress:
03495 e_ini = give_iniporosity(ipp, im+1, ido);
03496 break;
03497 case damage_plasticity:
03498 e_ini = give_iniporosity(ipp, im+1, ido);
03499 break;
03500 default:{
03501 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
03502 }
03503 }
03504
03505 return e_ini;
03506 }
03507
03508
03509
03510
03511
03512
03513
03514
03515
03516
03517
03518
03519
03520
03521 double mechmat::give_consparam (long ipp, long im, long ido)
03522 {
03523 long i;
03524 double cp;
03525
03526 i = ip[ipp].idm[im];
03527
03528 switch (ip[ipp].tm[im]){
03529 case simplas1d:{
03530 cp = spl1d[i].give_consparam (ipp,ido);
03531 break;
03532 }
03533 case jflow:{
03534 cp = j2f[i].give_consparam (ipp,ido);
03535 break;
03536 }
03537 case mohrcoul:{
03538 cp = mcoul[i].give_consparam (ipp,ido);
03539 break;
03540 }
03541 case mohrcoulparab:{
03542 cp = mcpar[i].give_consparam (ipp,ido);
03543 break;
03544 }
03545 case boermaterial:{
03546 cp = boerm[i].give_consparam (ipp,ido);
03547 break;
03548 }
03549 case druckerprager:{
03550 cp = drprm[i].give_consparam (ipp,ido);
03551 break;
03552 }
03553 case druckerprager2:{
03554 cp = drprm2[i].give_consparam (ipp,ido);
03555 break;
03556 }
03557 case modcamclaymat:{
03558 cp = cclay[i].give_consparam (ipp,ido);
03559 break;
03560 }
03561 case modcamclaycoupmat:{
03562 cp = cclayc[i].give_consparam (ipp,ido);
03563 break;
03564 }
03565 case hissplasmat:{
03566 cp = hisspl[i].give_consparam (ipp,ido);
03567 break;
03568 }
03569 case damage_plasticity:
03570 cp = give_consparam(ipp, im+1, ido);
03571 break;
03572 case effective_stress:
03573 cp = give_consparam(ipp, im+1, ido);
03574 break;
03575 default:{
03576 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
03577 }
03578 }
03579
03580 return cp;
03581 }
03582
03583
03584
03585
03586
03587
03588
03589
03590
03591
03592
03593
03594
03595
03596 double mechmat::give_dampar(long ipp, long im, long ido)
03597 {
03598 double dam = 0.0;
03599 switch (ip[ipp].tm[im])
03600 {
03601 case scaldamage:{
03602 dam = scdam[ip[ipp].idm[im]].givedamage(ipp, ido);
03603 break;
03604 }
03605 case nonlocdamgmat:{
03606 dam = give_dampar(ipp, im+1, ido);
03607 break;
03608 }
03609 case effective_stress:{
03610 dam = give_dampar (ipp, im+1, ido);
03611 break;
03612 }
03613 default:{
03614 print_err("unknown material model is required",__FILE__,__LINE__,__func__);
03615 }
03616 }
03617 return dam;
03618 }
03619
03620
03621
03622
03623
03624
03625
03626
03627
03628
03629
03630
03631 double mechmat::give_proczonelength (long ipp, long im, long ido)
03632 {
03633 double dam = 0.0;
03634 switch (ip[ipp].tm[im])
03635 {
03636 case scaldamage:{
03637 dam = scdam[ip[ipp].idm[im]].give_proczonelength(ipp, ido);
03638 break;
03639 }
03640 case effective_stress:{
03641 dam = give_proczonelength (ipp, im+1, ido);
03642 break;
03643 }
03644 case nonlocdamgmat:
03645 dam = give_proczonelength (ipp, im+1, ido);
03646 break;
03647 default:{
03648 print_err("unknown material model is required",__FILE__,__LINE__,__func__);
03649 }
03650 }
03651 return dam;
03652 }
03653
03654
03655
03656
03657
03658
03659
03660
03661
03662
03663
03664
03665
03666 double mechmat::give_crackwidth (long ipp, long im, long ido)
03667 {
03668 double dam = 0.0;
03669 switch (ip[ipp].tm[im])
03670 {
03671 case scaldamage:{
03672 dam = scdam[ip[ipp].idm[im]].give_crackwidth (ipp, ido);
03673 break;
03674 }
03675 case effective_stress:{
03676 dam = give_crackwidth (ipp, im+1, ido);
03677 break;
03678 }
03679 case nonlocdamgmat:
03680 dam = give_crackwidth (ipp, im+1, ido);
03681 break;
03682 default:{
03683 print_err("unknown material model is required",__FILE__,__LINE__,__func__);
03684 }
03685 }
03686 return dam;
03687 }
03688
03689
03690
03691
03692
03693
03694
03695
03696
03697
03698
03699
03700
03701
03702
03703 void mechmat::give_kappa(long ipp,long im,long ido, vector &kappa)
03704 {
03705 switch (ip[ipp].tm[im])
03706 {
03707 case scaldamage:
03708 kappa[0] = ip[ipp].other[ido];
03709 break;
03710 case scaldamagecc:
03711 kappa[0] = ip[ipp].other[ido];
03712
03713 break;
03714 default:
03715 print_err ("unknown type of damage model is required",__FILE__,__LINE__,__func__);
03716 }
03717 return;
03718 }
03719
03720
03721
03722
03723
03724
03725
03726
03727
03728
03729
03730
03731
03732
03733
03734 void mechmat::giveirrstrains(long ipp, long im, long ido, vector &epsirr)
03735 {
03736 fillv(0.0, epsirr);
03737 switch (ip[ipp].tm[im])
03738 {
03739 case simplas1d:{
03740 spl1d[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03741 break;
03742 }
03743 case jflow:{
03744 j2f[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03745 break;
03746 }
03747 case jflow2:{
03748 j2f2[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03749 break;
03750 }
03751 case mohrcoul:{
03752 mcoul[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03753 break;
03754 }
03755 case mohrcoulparab:{
03756 mcpar[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03757 break;
03758 }
03759 case boermaterial:{
03760 boerm[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03761 break;
03762 }
03763 case druckerprager:{
03764 drprm[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03765 break;
03766 }
03767 case druckerprager2:{
03768 drprm2[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03769 break;
03770 }
03771 case chenplast:{
03772
03773 break;
03774 }
03775 case modcamclaymat:{
03776 cclay[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03777 break;
03778 }
03779 case modcamclaycoupmat:{
03780 cclayc[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03781 break;
03782 }
03783 case hissplasmat:{
03784 hisspl[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03785 break;
03786 }
03787 case effective_stress:{
03788 giveirrstrains(ipp, im+1, ido, epsirr);
03789 break;
03790 }
03791 case creepbaz:{
03792 crbaz[ip[ipp].idm[im]].giveirrstrains(ipp, ido, epsirr);
03793 break;
03794 }
03795 case creepeffym:{
03796 giveirrstrains(ipp, im+1, ido, epsirr);
03797 break;
03798 }
03799 case creepb3:
03800 case creeprs:
03801 case creepdpl:{
03802 creep_giveirrstrains(ipp, im, ido, epsirr);
03803 break;
03804 }
03805 case shrinkagemat:
03806 shmat[ip[ipp].idm[im]].giveirrstrains(ipp, im, ido, epsirr);
03807 default:
03808 print_err("unknown material model is required",__FILE__,__LINE__,__func__);
03809 }
03810 }
03811
03812
03813
03814
03815
03816
03817
03818
03819
03820
03821
03822
03823
03824
03825
03826
03827
03828
03829
03830
03831
03832
03833
03834
03835
03836
03837
03838
03839
03840
03841 long mechmat::search_reqnmq(nonmechquant* &rnmq)
03842 {
03843 long i, j, n;
03844 long anmq[tnknmq];
03845 memset(anmq, 0, sizeof(*anmq)*tnknmq);
03846
03847
03848 for (i=0; i<tnip; i++)
03849 {
03850 for (j=0; j<ip[i].nm; j++)
03851 give_reqnmq(i, j, anmq);
03852 }
03853
03854
03855 n = 0;
03856 for (i=0; i<tnknmq; i++)
03857 {
03858 if (anmq[i] == 1)
03859 n++;
03860 }
03861
03862
03863 if (n == 0)
03864 {
03865 rnmq = NULL;
03866 return n;
03867 }
03868 else
03869 {
03870
03871 rnmq = new nonmechquant[n];
03872 j = 0;
03873 for (i=0; i<tnknmq; i++)
03874 {
03875 if (anmq[i] == 1){
03876 rnmq[j] = nonmechquant(i+1);
03877 j++;
03878 }
03879 }
03880 }
03881
03882 return n;
03883 }
03884
03885
03886
03887
03888
03889
03890
03891
03892
03893
03894
03895
03896
03897
03898
03899
03900
03901
03902
03903 void mechmat::give_reqnmq(long ipp, long im, long *anmq)
03904 {
03905 switch(ip[ipp].tm[im])
03906 {
03907 case elisomat:
03908 case elgmat3d:
03909 case winklerpasternak:
03910 case graphm:
03911 case varelisomat:
03912 case geoelast:
03913 case simplas1d:
03914 case jflow:
03915 case jflow2:
03916 case microplaneM4:
03917 case microsimp:
03918 case microfibro:
03919 case mohrcoul:
03920 case mohrcoulparab:
03921 case boermaterial:
03922 case druckerprager:
03923 case druckerprager2:
03924 case chenplast:
03925 case modcamclaymat:
03926 case shefpl:
03927 case scaldamage:
03928 case scaldamagecc:
03929 case ortodamage:
03930 case ortodamagerot:
03931 case anisodamage:
03932 case anisodamagerot:
03933 case creepbaz:
03934 case creepeffym:
03935 case contmat:
03936 case cebfipcontmat:
03937 case consolidation:
03938 case nonlocplastmat:
03939 case nonlocdamgmat:
03940 case nonlocalmicroM4:
03941 case glasgowdamage:
03942 case relaxationeuro:
03943 case simvisc:
03944 case lemaitr:
03945 case glasgowmechmat:
03946 case viscoplasticity:
03947 case lenjonespot:
03948 case creep_damage:
03949 case creepdpl:
03950 case damage_plasticity:
03951 break;
03952
03953
03954
03955
03956 case creeprs:
03957 crrs[im].give_reqnmq(anmq);
03958 break;
03959 case creepb3:
03960 crb3[im].give_reqnmq(anmq);
03961 break;
03962 case modcamclaycoupmat:
03963 cclayc[im].give_reqnmq(anmq);
03964 break;
03965 case hypoplastmat:
03966 hypopl[im].give_reqnmq(anmq);
03967 break;
03968 case therisodilat:
03969 tidilat[im].give_reqnmq(anmq);
03970 break;
03971 case effective_stress:
03972 effstr[im].give_reqnmq(anmq);
03973 break;
03974 case shrinkagemat:
03975 shmat[im].give_reqnmq(anmq);
03976 break;
03977 case time_switchmat:
03978
03979 break;
03980
03981 default:
03982 print_err("unknown material type is required", __FILE__, __LINE__, __func__);
03983 }
03984
03985 return;
03986 }
03987
03988
03989
03990
03991
03992
03993
03994
03995
03996
03997
03998
03999
04000
04001
04002 double mechmat::givenonmechq(nonmechquant qt, long ipp)
04003 {
04004 long id = nmqid[qt-1];
04005 if (id < 0)
04006 {
04007 print_err("Required quantity %s is not used in the problem solved",
04008 __FILE__, __LINE__, __func__, nonmechquantstr[qt-1].alias);
04009 abort();
04010 }
04011 return nonmechq[id][ipp];
04012 }
04013
04014
04015
04016
04017
04018
04019
04020
04021
04022
04023
04024
04025
04026
04027
04028
04029 void mechmat::storenonmechq(nonmechquant qt, long ipp, double val)
04030 {
04031 long id = nmqid[qt-1];
04032 if (id < 0)
04033 {
04034 print_err("Required quantity %s is not used in the problem solved",
04035 __FILE__, __LINE__, __func__, nonmechquantstr[qt-1].alias);
04036 abort();
04037 }
04038 nonmechq[id][ipp] = val;
04039 }
04040
04041
04042
04043
04044
04045
04046
04047
04048
04049
04050
04051
04052
04053
04054 long mechmat::givestatusnmq (nonmechquant qt)
04055 {
04056 if (nmqid[qt-1] < 0)
04057 return 0;
04058
04059 return 1;
04060 }
04061
04062
04063
04064
04065
04066
04067
04068
04069
04070
04071
04072
04073
04074
04075 long mechmat::givenonmechqid (nonmechquant qt)
04076 {
04077 return nmqid[qt-1];
04078 }
04079
04080
04081
04082
04083
04084
04085
04086
04087
04088
04089
04090
04091
04092
04093
04094
04095
04096
04097
04098
04099
04100
04101
04102
04103
04104
04105
04106
04107
04108 void mechmat::matstiff (matrix &d,long ipp,long im,long ido)
04109 {
04110 long i;
04111 switch (ip[ipp].tm[im]){
04112 case elisomat:{
04113 i=ip[ipp].idm[im];
04114 eliso[i].matstiff (d,ip[ipp].ssst);
04115 break;
04116 }
04117 case elgmat3d:{
04118 i=ip[ipp].idm[im];
04119 elgm3d[i].matstiff (d);
04120 break;
04121 }
04122 case simplas1d:{
04123 i=ip[ipp].idm[im];
04124 spl1d[i].matstiff (d,ipp,ido);
04125 break;
04126 }
04127 case jflow:{
04128 i=ip[ipp].idm[im];
04129 j2f[i].matstiff (d,ipp,ido);
04130 break;
04131 }
04132 case jflow2:{
04133 i=ip[ipp].idm[im];
04134 j2f2[i].matstiff (d,ipp,ido);
04135 break;
04136 }
04137 case microplaneM4:{
04138 i=ip[ipp].idm[im];
04139 mpM4[i].matstiff (d,ipp,ido);
04140 break;
04141 }
04142 case microsimp:{
04143 i=ip[ipp].idm[im];
04144 mpSIM[i].matstiff (d,ipp,ido);
04145 break;
04146 }
04147 case microfibro:{
04148 i=ip[ipp].idm[im];
04149 mpfib[i].matstiff (d,ipp,ido);
04150 break;
04151 }
04152 case mohrcoul:{
04153 i=ip[ipp].idm[im];
04154 mcoul[i].matstiff (d,ipp,ido);
04155 break;
04156 }
04157 case mohrcoulparab:{
04158 i=ip[ipp].idm[im];
04159 mcpar[i].matstiff (d,ipp,ido);
04160 break;
04161 }
04162 case boermaterial:{
04163 i=ip[ipp].idm[im];
04164 boerm[i].matstiff (d,ipp,ido);
04165 break;
04166 }
04167 case druckerprager:{
04168 i=ip[ipp].idm[im];
04169 drprm[i].matstiff (d,ipp,ido);
04170 break;
04171 }
04172 case druckerprager2:{
04173 i=ip[ipp].idm[im];
04174 drprm2[i].matstiff (d,ipp,ido);
04175 break;
04176 }
04177 case chenplast:{
04178 i=ip[ipp].idm[im];
04179 chplast[i].matstiff (d,ipp,ido);
04180 break;
04181 }
04182 case modcamclaymat:{
04183 i=ip[ipp].idm[im];
04184 cclay[i].matstiff (d,ipp,im,ido);
04185 break;
04186 }
04187 case modcamclaycoupmat:{
04188 i=ip[ipp].idm[im];
04189 cclayc[i].matstiff (d,ipp,ido);
04190 break;
04191 }
04192 case hypoplastmat:{
04193 i=ip[ipp].idm[im];
04194 hypopl[i].matstiff (d,ipp,ido);
04195 break;
04196 }
04197 case shefpl:{
04198 i=ip[ipp].idm[im];
04199 shpl[i].matstiff (d,ipp,ido);
04200 break;
04201 }
04202 case hissplasmat:{
04203 i=ip[ipp].idm[im];
04204 hisspl[i].matstiff (d,ipp,ido);
04205 break;
04206 }
04207 case glasgowmechmat:{
04208 i=ip[ipp].idm[im];
04209 glasgmat[i].matstiff (d,ipp,ido);
04210 break;
04211 }
04212 case glasgowdamage:{
04213 i=ip[ipp].idm[im];
04214 glasgdam[i].matstiff (d,ipp,ido);
04215 break;
04216 }
04217 case creep_damage:{
04218 i=ip[ipp].idm[im];
04219 crdam[i].matstiff (d,ipp,im,ido);
04220 break;
04221 }
04222 case time_switchmat:{
04223 i=ip[ipp].idm[im];
04224 tswmat[i].matstiff (d,ipp,im,ido);
04225 break;
04226 }
04227 case effective_stress:{
04228 matstiff (d,ipp,im+1,ido);
04229 break;
04230 }
04231
04232
04233
04234
04235
04236
04237
04238 case scaldamage:{
04239 i=ip[ipp].idm[im];
04240 scdam[i].matstiff (d,ipp,ido);
04241 break;
04242 }
04243 case scaldamagecc:{
04244 i=ip[ipp].idm[im];
04245 scdamcc[i].matstiff (d,ipp,ido);
04246 break;
04247 }
04248 case ortodamage:{
04249 i=ip[ipp].idm[im];
04250 ortdam[i].matstiff (d,ipp,ido);
04251 break;
04252 }
04253 case ortodamagerot:{
04254 i=ip[ipp].idm[im];
04255 ortdamrot[i].matstiff (d,ipp,ido);
04256 break;
04257 }
04258 case anisodamage:{
04259 i=ip[ipp].idm[im];
04260 anidam[i].matstiff (d,ipp,ido);
04261 break;
04262 }
04263 case anisodamagerot:{
04264 i=ip[ipp].idm[im];
04265 anidamrot[i].matstiff (d,ipp,ido);
04266 break;
04267 }
04268 case graphm:{
04269 i=ip[ipp].idm[im];
04270 grmat[i].matstiff (d,ipp);
04271 break;
04272 }
04273 case varelisomat:{
04274 i=ip[ipp].idm[im];
04275 veliso[i].matstiff (d,ipp);
04276 break;
04277 }
04278 case geoelast:{
04279 i=ip[ipp].idm[im];
04280 geoel[i].matstiff (d,ipp,ido);
04281 break;
04282 }
04283 case creeprs:
04284 case creepb3:{
04285 creep_matstiff (d,ipp,im,ido);
04286 break;
04287 }
04288 case creepdpl:{
04289 creep_matstiff (d,ipp,im,ido);
04290 break;
04291 }
04292 case creepbaz:{
04293 i=ip[ipp].idm[im];
04294 crbaz[i].matstiff (d,ipp);
04295 break;
04296 }
04297 case creepeffym:{
04298 i=ip[ipp].idm[im];
04299 creffym[i].matstiff (d,ipp,im,ido);
04300 break;
04301 }
04302 case consolidation:{
04303 i=ip[ipp].idm[im];
04304 csol[i].matstiff (d,ipp);
04305 break;
04306 }
04307 case winklerpasternak:{
04308 i=ip[ipp].idm[im];
04309 wpast[i].matstiff (d,ipp);
04310 break;
04311 }
04312
04313 case nonlocdamgmat:{
04314 matstiff (d,ipp,im+1,ido);
04315 break;
04316 }
04317 case nonlocplastmat:{
04318 matstiff (d,ipp,im+1,ido);
04319 break;
04320 }
04321 case nonlocalmicroM4:{
04322 matstiff (d,ipp,im+1,ido);
04323 break;
04324 }
04325
04326 case damage_plasticity:{
04327 i=ip[ipp].idm[im];
04328 dampl[i].matstiff (d,ipp);
04329 break;
04330 }
04331
04332 case viscoplasticity:{
04333 i=ip[ipp].idm[im];
04334 visplas[i].matstiff (d,ipp,im,ido);
04335 break;
04336 }
04337
04338 case contmat:{
04339 i=ip[ipp].idm[im];
04340 conmat[i].matstiff (d,ipp);
04341 break;
04342 }
04343
04344 case cebfipcontmat:{
04345 i=ip[ipp].idm[im];
04346 cebfipconmat[i].matstiff (d,ipp);
04347 break;
04348 }
04349 case shrinkagemat:{
04350 i=ip[ipp].idm[im];
04351 shmat[i].matstiff (d,ipp,im,ido);
04352 break;
04353 }
04354
04355 default:
04356 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
04357
04358 }
04359 }
04360
04361
04362
04363
04364
04365
04366
04367
04368
04369
04370
04371
04372
04373
04374
04375 void mechmat::elmatstiff (matrix &d,long ipp)
04376
04377 {
04378 long i,idem;
04379
04380 idem = ip[ipp].gemid();
04381 switch (ip[ipp].tm[idem]){
04382 case elisomat:{
04383 i=ip[ipp].idm[idem];
04384 eliso[i].elmatstiff (d,ip[ipp].ssst);
04385 break;
04386 }
04387 case elgmat3d:{
04388 i=ip[ipp].idm[idem];
04389 elgm3d[i].elmatstiff (d);
04390 break;
04391 }
04392 case varelisomat:{
04393 i=ip[ipp].idm[idem];
04394 veliso[i].elmatstiff (d,ipp);
04395 break;
04396 }
04397 default:
04398 print_err("Unknown elastic material type is required",__FILE__,__LINE__,__func__);
04399 }
04400 }
04401
04402
04403
04404
04405
04406
04407
04408
04409
04410
04411
04412
04413
04414
04415 void mechmat::elmatcompl (matrix &c,long ipp)
04416
04417 {
04418 long i, idem;
04419
04420 idem = ip[ipp].gemid();
04421 switch (ip[ipp].tm[idem]){
04422 case elisomat:{
04423 i=ip[ipp].idm[idem];
04424 eliso[i].matcompl (c,ip[ipp].ssst);
04425 break;
04426 }
04427 case varelisomat:{
04428 i=ip[ipp].idm[idem];
04429 veliso[i].matcompl (c,ipp);
04430 break;
04431 }
04432 default:
04433 print_err("unknown elastic material type is required",__FILE__,__LINE__,__func__);
04434
04435 }
04436 }
04437
04438
04439
04440
04441
04442
04443
04444
04445
04446
04447
04448
04449
04450
04451 void mechmat::eigenstresses (vector &sig,long ipp,long id)
04452 {
04453 long i,n;
04454 n=sig.n;
04455 vector eps;
04456 allocv (n,eps);
04457
04458 switch (ip[ipp].tm[id]){
04459 case relaxationeuro:{
04460
04461 i=ip[ipp].idm[id];
04462
04463 relaxec[i].stress (sig,eps,ip[ipp].ssst);
04464
04465 storeeigstrain (ipp,eps);
04466
04467 storeeigstress (ipp,sig);
04468 break;
04469 }
04470 default:
04471 print_err("Unknown model for eigenstresses computation is required",__FILE__,__LINE__,__func__);
04472 }
04473 destrv (eps);
04474 }
04475
04476
04477
04478
04479
04480
04481
04482
04483
04484
04485
04486
04487
04488
04489
04490
04491 void mechmat::computenlstresses (long ipp,long im,long ido)
04492 {
04493 switch (ip[ipp].tm[im]){
04494 case elisomat:{
04495 eliso[ip[ipp].idm[im]].nlstresses (ipp);
04496 break;
04497 }
04498 case elgmat3d:{
04499 elgm3d[ip[ipp].idm[im]].nlstresses (ipp);
04500 break;
04501 }
04502 case simplas1d:{
04503 spl1d[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04504 break;
04505 }
04506 case jflow:{
04507 j2f[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04508 break;
04509 }
04510 case jflow2:{
04511 j2f2[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04512 break;
04513 }
04514 case microplaneM4:{
04515 mpM4[ip[ipp].idm[im]].nlstresses (ipp,ido);
04516 break;
04517 }
04518 case microsimp:{
04519 mpSIM[ip[ipp].idm[im]].nlstresses (ipp,ido);
04520 break;
04521 }
04522 case microfibro:{
04523 mpfib[ip[ipp].idm[im]].nlstresses (ipp,ido);
04524 break;
04525 }
04526 case mohrcoul:{
04527 mcoul[ip[ipp].idm[im]].nlstresses (ipp,ido);
04528 break;
04529 }
04530 case mohrcoulparab:{
04531 mcpar[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04532 break;
04533 }
04534 case boermaterial:{
04535 boerm[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04536 break;
04537 }
04538 case druckerprager:{
04539 drprm[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04540 break;
04541 }
04542 case druckerprager2:{
04543 drprm2[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04544 break;
04545 }
04546 case chenplast:{
04547 chplast[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04548 break;
04549 }
04550 case modcamclaymat:{
04551 cclay[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04552 break;
04553 }
04554 case modcamclaycoupmat:{
04555 cclayc[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04556 break;
04557 }
04558 case hypoplastmat:{
04559 hypopl[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04560 break;
04561 }
04562 case hissplasmat:{
04563 hisspl[ip[ipp].idm[im]].nlstresses (ipp,ido);
04564 break;
04565 }
04566 case shefpl:{
04567 shpl[ip[ipp].idm[im]].nlstresses (ipp,ido);
04568 break;
04569 }
04570 case glasgowmechmat:{
04571 glasgmat[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04572 break;
04573 }
04574 case glasgowdamage:{
04575 glasgdam[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04576 break;
04577 }
04578 case creep_damage:{
04579 crdam[ip[ipp].idm[im]].nlstresses (ipp, im, ido);
04580 break;
04581 }
04582 case shrinkagemat:{
04583 shmat[ip[ipp].idm[im]].nlstresses (ipp, im, ido);
04584 break;
04585 }
04586 case time_switchmat:{
04587 tswmat[ip[ipp].idm[im]].nlstresses (ipp, im, ido);
04588 break;
04589 }
04590 case effective_stress:{
04591 effstr[ip[ipp].idm[im]].nlstresses (ipp, im, ido);
04592 break;
04593 }
04594
04595
04596
04597
04598
04599
04600
04601
04602
04603
04604
04605
04606 case scaldamage:{
04607 scdam[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04608 break;
04609 }
04610 case scaldamagecc:{
04611 scdamcc[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04612 break;
04613 }
04614 case ortodamage:{
04615 ortdam[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04616 break;
04617 }
04618 case ortodamagerot:{
04619 ortdamrot[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04620 break;
04621 }
04622 case anisodamage:{
04623 anidam[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04624 break;
04625 }
04626 case anisodamagerot:{
04627 anidamrot[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04628 break;
04629 }
04630 case graphm:{
04631 grmat[ip[ipp].idm[im]].nlstresses (ipp);
04632 break;
04633 }
04634 case varelisomat:{
04635 veliso[ip[ipp].idm[im]].nlstresses (ipp);
04636 break;
04637 }
04638 case geoelast:{
04639 geoel[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04640 break;
04641 }
04642 case creepb3:
04643 case creeprs:
04644 case creepdpl:{
04645 creep_nlstresses(ipp,im,ido);
04646 break;
04647 }
04648 case creepbaz:{
04649 crbaz[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04650 break;
04651 }
04652 case creepeffym:{
04653 creffym[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04654 break;
04655 }
04656 case consolidation:{
04657 csol[ip[ipp].idm[im]].nlstresses (ipp);
04658 break;
04659 }
04660 case winklerpasternak:{
04661 wpast[ip[ipp].idm[im]].nlstresses (ipp);
04662 break;
04663 }
04664
04665 case nonlocdamgmat:{
04666 computenlstresses (ipp,im+1,ido);
04667 break;
04668 }
04669 case nonlocplastmat:{
04670 computenlstresses (ipp,im+1,ido);
04671 break;
04672 }
04673 case nonlocalmicroM4:{
04674 computenlstresses (ipp,im+1,ido);
04675 break;
04676 }
04677 case damage_plasticity:{
04678 dampl[ip[ipp].idm[im]].nlstresses (ipp);
04679 break;
04680 }
04681 case viscoplasticity:{
04682 visplas[ip[ipp].idm[im]].nlstresses (ipp,ido);
04683 break;
04684 }
04685
04686
04687 case contmat:{
04688 conmat[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04689 break;
04690 }
04691
04692 case cebfipcontmat:{
04693 cebfipconmat[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04694 break;
04695 }
04696
04697 default:{
04698 print_err("unknown material type is required", __FILE__,__LINE__,__func__);
04699 }
04700 }
04701
04702 }
04703
04704
04705
04706
04707
04708
04709
04710
04711
04712
04713
04714
04715
04716
04717
04718 void mechmat::computenlstressesincr (long ipp,long im,long ido)
04719 {
04720 switch (ip[ipp].tm[im]){
04721
04722 case creep_damage:{
04723 crdam[ip[ipp].idm[im]].nlstressesincr (ipp, im, ido);
04724 break;
04725 }
04726 case shrinkagemat:{
04727 shmat[ip[ipp].idm[im]].nlstressesincr (ipp, im, ido);
04728 break;
04729 }
04730
04731
04732
04733
04734
04735
04736
04737
04738
04739
04740
04741
04742 case creepb3:
04743 case creeprs:
04744 case creepdpl:{
04745 creep_nlstressesincr (ipp,im,ido);
04746 break;
04747 }
04748
04749
04750
04751
04752
04753
04754
04755
04756
04757
04758
04759 case consolidation:{
04760 csol[ip[ipp].idm[im]].nlstressesincr (ipp);
04761 break;
04762 }
04763
04764
04765 case viscoplasticity:{
04766 visplas[ip[ipp].idm[im]].nlstressesincr (ipp,im,ido);
04767 break;
04768 }
04769
04770 case elisomat:
04771 case simplas1d:
04772 case jflow:
04773 case jflow2:
04774 case microplaneM4:
04775 case microsimp:
04776 case microfibro:
04777 case mohrcoul:
04778 case mohrcoulparab:
04779 case boermaterial:
04780 case druckerprager:
04781 case druckerprager2:
04782 case chenplast:
04783 case modcamclaymat:
04784 case modcamclaycoupmat:
04785 case hypoplastmat:
04786 case shefpl:
04787 case effective_stress:
04788 case glasgowmechmat:
04789 case glasgowdamage:
04790 case scaldamage:
04791 case scaldamagecc:
04792 case ortodamage:
04793 case ortodamagerot:
04794 case anisodamage:
04795 case anisodamagerot:
04796 case graphm:
04797 case varelisomat:
04798 case geoelast:
04799 case creepeffym:
04800 case winklerpasternak:
04801 case nonlocdamgmat:
04802 case nonlocplastmat:
04803 case nonlocalmicroM4:
04804 case damage_plasticity:
04805 case cebfipcontmat:
04806 case contmat:
04807 if (im == 0)
04808 ip[ipp].clean_stresses(Mb->nlc);
04809 break;
04810 default:{
04811 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
04812 }
04813 }
04814
04815 }
04816
04817
04818
04819
04820
04821
04822
04823
04824
04825
04826
04827
04828
04829
04830
04831 void mechmat::compnonloc_nlstresses (long ipp,long im,long ido)
04832 {
04833 switch (ip[ipp].tm[im]){
04834 case elisomat:{
04835 eliso[ip[ipp].idm[im]].nlstresses (ipp);
04836 break;
04837 }
04838 case simplas1d:{
04839 if (ip[ipp].hmt & 2)
04840 spl1d[ip[ipp].idm[im]].nonloc_nlstresses (ipp,im,ido);
04841 else
04842 spl1d[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04843 break;
04844 }
04845 case jflow:{
04846 if (ip[ipp].hmt & 2)
04847 j2f[ip[ipp].idm[im]].nonloc_nlstresses (ipp,im,ido);
04848 else
04849 j2f[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04850 break;
04851 }
04852 case microplaneM4:{
04853 if (ip[ipp].hmt & 2)
04854 mpM4[ip[ipp].idm[im]].nonloc_nlstresses (ipp,ido);
04855 else
04856 mpM4[ip[ipp].idm[im]].nlstresses (ipp,ido);
04857 break;
04858 }
04859 case microsimp:{
04860
04861 break;
04862 }
04863 case microfibro:{
04864
04865 break;
04866 }
04867 case mohrcoul:{
04868 if (ip[ipp].hmt & 2)
04869 mcoul[ip[ipp].idm[im]].nonloc_nlstresses (ipp,ido);
04870 else
04871 mcoul[ip[ipp].idm[im]].nlstresses (ipp,ido);
04872 break;
04873 }
04874 case mohrcoulparab:{
04875 if (ip[ipp].hmt & 2)
04876 mcpar[ip[ipp].idm[im]].nonloc_nlstresses (ipp,im,ido);
04877 else
04878 mcpar[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04879 break;
04880 }
04881 case boermaterial:{
04882 if (ip[ipp].hmt & 2)
04883 boerm[ip[ipp].idm[im]].nonloc_nlstresses (ipp,im,ido);
04884 else
04885 boerm[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04886 break;
04887 }
04888 case druckerprager:{
04889 if (ip[ipp].hmt & 2)
04890 drprm[ip[ipp].idm[im]].nonloc_nlstresses (ipp,im,ido);
04891 else
04892 drprm[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04893 break;
04894 }
04895
04896 case druckerprager2:{
04897
04898
04899
04900
04901 break;
04902 }
04903 case modcamclaymat:{
04904
04905 break;
04906 }
04907 case hissplasmat:{
04908
04909 break;
04910 }
04911 case effective_stress:{
04912 compnonloc_nlstresses(ipp, im+1, ido);
04913 break;
04914 }
04915 case simvisplas:{
04916
04917 break;
04918 }
04919 case lemaitr:{
04920
04921 break;
04922 }
04923 case scaldamage:{
04924 scdam[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04925 break;
04926 }
04927 case scaldamagecc:{
04928 scdamcc[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04929 break;
04930 }
04931 case ortodamage:{
04932 ortdam[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04933 break;
04934 }
04935 case ortodamagerot:{
04936 ortdamrot[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04937 break;
04938 }
04939 case anisodamage:{
04940 anidam[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04941 break;
04942 }
04943 case anisodamagerot:{
04944 anidamrot[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04945 break;
04946 }
04947 case graphm:{
04948 grmat[ip[ipp].idm[im]].nlstresses (ipp);
04949 break;
04950 }
04951 case varelisomat:{
04952 veliso[ip[ipp].idm[im]].nlstresses (ipp);
04953 break;
04954 }
04955 case creepbaz:{
04956
04957 break;
04958 }
04959 case consolidation:{
04960 csol[ip[ipp].idm[im]].nlstresses (ipp);
04961 break;
04962 }
04963 case winklerpasternak:{
04964 wpast[ip[ipp].idm[im]].nlstresses (ipp);
04965 break;
04966 }
04967
04968 case creep_damage:{
04969 crdam[ip[ipp].idm[im]].nlstresses (ipp, im, ido);
04970 break;
04971 }
04972 case creepb3:
04973 case creeprs:
04974 case creepdpl:{
04975 creep_nlstresses(ipp,im,ido);
04976 break;
04977 }
04978 case creepeffym:{
04979 creffym[ip[ipp].idm[im]].nlstresses (ipp,im,ido);
04980 break;
04981 }
04982 case nonlocdamgmat:{
04983 compnonloc_nlstresses (ipp,im+1,ido);
04984 break;
04985 }
04986 case nonlocplastmat:{
04987 compnonloc_nlstresses (ipp,im+1,ido);
04988 break;
04989 }
04990 case nonlocalmicroM4:{
04991 compnonloc_nlstresses (ipp,im+1,ido);
04992 break;
04993 }
04994
04995 default:{
04996 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
04997 }
04998 }
04999
05000 }
05001
05002
05003
05004
05005
05006
05007
05008
05009
05010
05011
05012
05013 void mechmat::compute_averstrains ()
05014 {
05015 long i,ne;
05016 double volume;
05017 elemtype te;
05018 vector averstra (6);
05019
05020 volume=0.0;
05021 fillv (0.0,averstra);
05022
05023 ne=Mt->ne;
05024 for (i=0;i<ne;i++){
05025 if (Gtm->leso[i]==1){
05026
05027 te = Mt->give_elem_type (i);
05028
05029 switch (te){
05030 case linearhex:{
05031 Lhex->aver_strains (0,i,0,0,averstra,volume);
05032 break;
05033 }
05034 case quadrhex:{
05035 Qhex->aver_strains (0,i,0,0,averstra,volume);
05036 break;
05037 }
05038
05039 default:{
05040
05041 }
05042 }
05043 }
05044 }
05045
05046 for (i=0;i<6;i++){
05047 averstra[i]/=volume;
05048 }
05049
05050 printf ("\n\n prumerne deformace");
05051 for (i=0;i<6;i++){
05052 printf ("\n slozka %ld %f",i+1,averstra[i]);
05053 }
05054
05055 }
05056
05057
05058
05059
05060
05061
05062
05063
05064
05065
05066
05067
05068 void mechmat::temprstrains (long lcid)
05069 {
05070 long i,nm,ii;
05071
05072 for (i=0;i<tnip;i++){
05073
05074 ii = ip[i].hmt & 1;
05075 if (ii > 0){
05076 nm=ip[i].nm-1;
05077 switch (ip[i].tm[nm]){
05078 case therisodilat:{
05079 tidilat[ip[i].idm[nm]].temprstrains (lcid,i);
05080 break;
05081 }
05082 default:
05083 print_err("unknown thermal dilatation material is required",__FILE__,__LINE__,__func__);
05084 }
05085 }
05086
05087 }
05088 }
05089
05090
05091
05092
05093
05094
05095
05096
05097
05098
05099 void mechmat::eigstrmod ()
05100 {
05101 long i,j,k,nc,ipp,nip;
05102
05103 if (eigstrains != NULL){
05104 for (i=0;i<Mt->ne;i++){
05105 if (Gtm->leso[i]==1){
05106
05107 ipp=Mt->elements[i].ipp[0][0];
05108
05109 nip = Mt->give_tnip (i);
05110 for (j=0;j<nip;j++){
05111
05112 nc=ip[ipp].ncompstr;
05113 for (k=0;k<nc;k++)
05114 ip[ipp].strain[k]-=eigstrains[ipp][k];
05115 ipp++;
05116 }
05117 }
05118 }
05119 }
05120 if (tempstrains != NULL){
05121 for (i=0;i<Mt->ne;i++){
05122 if (Gtm->leso[i]==1){
05123
05124 ipp=Mt->elements[i].ipp[0][0];
05125
05126 nip = Mt->give_tnip (i);
05127 for (j=0;j<nip;j++){
05128
05129 nc=ip[ipp].ncompstr;
05130 for (k=0;k<nc;k++)
05131 ip[ipp].strain[k]-=tempstrains[ipp][k];
05132 ipp++;
05133 }
05134 }
05135 }
05136 }
05137 }
05138
05139
05140
05141
05142
05143
05144
05145
05146
05147
05148
05149 void mechmat::add_macro_strains (vector ¯ostrains)
05150 {
05151 long i,j,ncompstr;
05152
05153 for (i=0;i<tnip;i++){
05154 ncompstr = ip[i].ncompstr;
05155 if (ncompstr != macrostrains.n){
05156 print_err("the number of strain components (%ld) in integration point %ld\n"
05157 "is not equal to the number of macro-strain component (%ld) in homogenization",
05158 __FILE__, __LINE__, __func__,ncompstr,i,macrostrains.n);
05159 abort ();
05160 }
05161 for (j=0;j<ncompstr;j++){
05162 ip[i].strain[j] += macrostrains[j];
05163 }
05164 }
05165 }
05166
05167
05168
05169
05170
05171
05172
05173
05174 void mechmat::totstrains ()
05175 {
05176 long i,j,k,nc,ipp,nip;
05177
05178 if (eigstrains != NULL){
05179 for (i=0;i<Mt->ne;i++){
05180 if (Gtm->leso[i]==1){
05181
05182 ipp=Mt->elements[i].ipp[0][0];
05183
05184 nip = Mt->give_tnip (i);
05185 for (j=0;j<nip;j++){
05186
05187 nc=ip[ipp].ncompstr;
05188 for (k=0;k<nc;k++)
05189 ip[ipp].strain[k]+=eigstrains[ipp][k];
05190 ipp++;
05191 }
05192 }
05193 }
05194 }
05195 if (tempstrains != NULL){
05196 for (i=0;i<Mt->ne;i++){
05197 if (Gtm->leso[i]==1){
05198
05199 ipp=Mt->elements[i].ipp[0][0];
05200
05201 nip = Mt->give_tnip (i);
05202 for (j=0;j<nip;j++){
05203
05204 nc=ip[ipp].ncompstr;
05205 for (k=0;k<nc;k++)
05206 ip[ipp].strain[k]+=tempstrains[ipp][k];
05207 ipp++;
05208 }
05209 }
05210 }
05211 }
05212 }
05213
05214
05215
05216
05217
05218
05219
05220
05221
05222
05223
05224
05225
05226
05227 void mechmat::matdilat (matrix &d,long ipp)
05228 {
05229 long i, iddm = ip[ipp].nm-1;
05230
05231 switch (ip[ipp].tm[iddm]){
05232 case therisodilat:{
05233 i=ip[ipp].idm[iddm];
05234 tidilat[i].matdilat (d,ip[ipp].ssst);
05235 break;
05236 }
05237 default:{
05238 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
05239 }
05240 }
05241
05242 }
05243
05244
05245
05246
05247
05248
05249
05250
05251
05252
05253
05254
05255
05256 void mechmat::storeipvol (long ipp,double vol)
05257 {
05258 ipv[ipp]=vol;
05259 }
05260
05261
05262
05263
05264
05265
05266
05267
05268
05269
05270
05271
05272
05273
05274 void mechmat::nonlocaverage (long ipp,long im,long ido)
05275 {
05276 if (ip[ipp].hmt & 2)
05277 {
05278 switch (ip[ipp].tm[im]){
05279 case nonlocplastmat :
05280 nlplast[ip[ipp].idm[im]].average(ipp, ido);
05281 break;
05282 case nonlocdamgmat :
05283 nldamg[ip[ipp].idm[im]].average(ipp,im,ido);
05284 break;
05285 case nonlocalmicroM4 :
05286 nmpM4[ip[ipp].idm[im]].average(ipp);
05287 break;
05288 case creep_damage:
05289 long tncompo;
05290 long nldamncompo;
05291 tncompo = givencompeqother(ipp,im);
05292 nldamncompo = givencompeqother(ipp,im+2);
05293 nonlocaverage(ipp, im+2, ido+tncompo-nldamncompo);
05294 break;
05295 default:{
05296 print_err("unknown type of averaging model is required", __FILE__, __LINE__, __func__);
05297 }
05298 }
05299 }
05300 }
05301
05302
05303
05304
05305
05306
05307
05308
05309
05310
05311
05312
05313
05314
05315 long mechmat::give_num_averq(long ipp,long im)
05316 {
05317 long ret = 0;
05318 if (ip[ipp].hmt & 2)
05319 {
05320 switch (ip[ipp].tm[im]){
05321 case nonlocalmicroM4 :
05322 ret = nmpM4[ip[ipp].idm[im]].give_num_averq(ipp);
05323 break;
05324 case nonlocplastmat :
05325 ret = nlplast[ip[ipp].idm[im]].give_num_averq(ipp);
05326 break;
05327 case nonlocdamgmat :
05328 ret = nldamg[ip[ipp].idm[im]].give_num_averq(ipp, im);
05329 break;
05330 case creep_damage:
05331 ret = give_num_averq(ipp, im+2);
05332 break;
05333 default:{
05334 print_err("unknown type of averaging model is required", __FILE__, __LINE__, __func__);
05335 }
05336 }
05337 }
05338 return ret;
05339 }
05340
05341
05342
05343
05344
05345
05346
05347
05348
05349
05350
05351
05352
05353
05354 void mechmat::give_aver_quantv(long ipp,long im,long ido, vector &qv)
05355 {
05356 if (ip[ipp].hmt & 2)
05357 {
05358 switch (ip[ipp].tm[im]){
05359 case nonlocplastmat :
05360 nlplast[ip[ipp].idm[im]].give_aver_quantv(ipp, im, ido, qv);
05361 break;
05362 case nonlocdamgmat :
05363 nldamg[ip[ipp].idm[im]].give_aver_quantv(ipp, im, ido, qv);
05364 break;
05365 case creep_damage:
05366 long tncompo;
05367 long nldamncompo;
05368 tncompo = givencompeqother(ipp,im);
05369 nldamncompo = givencompeqother(ipp,im+2);
05370 give_aver_quantv(ipp, im+2, ido+tncompo-nldamncompo, qv);
05371 break;
05372 default:{
05373 print_err("unknown type of averaging model is required", __FILE__, __LINE__, __func__);
05374 }
05375 }
05376 }
05377 return;
05378 }
05379
05380
05381
05382
05383
05384
05385
05386
05387
05388
05389
05390
05391
05392 double mechmat::nonlocradius (long ipp,long im)
05393 {
05394 double r;
05395
05396 switch (ip[ipp].tm[im]){
05397 case nonlocplastmat :
05398 {
05399 r = nlplast[ip[ipp].idm[im]].r;
05400 break;
05401 }
05402 case nonlocdamgmat :
05403 {
05404 r = nldamg[ip[ipp].idm[im]].r;
05405 break;
05406 }
05407 case nonlocalmicroM4 :
05408 {
05409 r = nmpM4[ip[ipp].idm[im]].r;
05410 break;
05411 }
05412 default:{
05413 print_err("wrong material type is required",__FILE__,__LINE__,__func__);
05414 }
05415 }
05416
05417 return r;
05418 }
05419
05420
05421
05422
05423
05424
05425
05426
05427
05428
05429
05430
05431
05432 long mechmat::givenonlocid(long ipp)
05433 {
05434 long i, ret=-1;
05435 for(i=0; i<ip[ipp].nm; i++)
05436 {
05437 switch(ip[ipp].tm[i])
05438 {
05439 case nonlocplastmat:
05440 case nonlocdamgmat:
05441 case nonlocalmicroM4:
05442 ret = i;
05443 break;
05444 default:
05445 break;
05446 }
05447 if (ret >= 0)
05448 break;
05449 }
05450 return ret;
05451 }
05452
05453
05454
05455
05456
05457
05458
05459
05460
05461
05462
05463
05464 void mechmat::updateipval (void)
05465 {
05466 long i,j,ipp,nip;
05467
05468
05469 for (i=0;i<Mt->ne;i++){
05470 if (Gtm->leso[i]==1){
05471 ipp=Mt->elements[i].ipp[0][0];
05472 nip=Mt->give_tnip (i);
05473 for (j=0;j<nip;j++){
05474 updateipvalmat (ipp,0,0);
05475 ipp++;
05476 }
05477 }
05478 }
05479 }
05480
05481
05482
05483
05484
05485
05486
05487
05488
05489
05490
05491
05492
05493
05494 void mechmat::updateipvalmat (long ipp,long im,long ido)
05495 {
05496 switch (ip[ipp].tm[im]){
05497 case elisomat :
05498 break;
05499 case simplas1d:{
05500 spl1d[ip[ipp].idm[im]].updateval (ipp,ido);
05501 break;
05502 }
05503 case jflow:{
05504 j2f[ip[ipp].idm[im]].updateval (ipp,im,ido);
05505 break;
05506 }
05507 case jflow2:{
05508 j2f2[ip[ipp].idm[im]].updateval (ipp,im,ido);
05509 break;
05510 }
05511 case mohrcoul:{
05512 mcoul[ip[ipp].idm[im]].updateval (ipp,im,ido);
05513 break;
05514 }
05515 case mohrcoulparab:{
05516 mcpar[ip[ipp].idm[im]].updateval (ipp,im,ido);
05517 break;
05518 }
05519 case boermaterial:{
05520 boerm[ip[ipp].idm[im]].updateval (ipp,im,ido);
05521 break;
05522 }
05523 case druckerprager:{
05524 drprm[ip[ipp].idm[im]].updateval (ipp,im,ido);
05525 break;
05526 }
05527 case druckerprager2:{
05528 drprm2[ip[ipp].idm[im]].updateval (ipp,im,ido);
05529 break;
05530 }
05531 case chenplast:{
05532 chplast[ip[ipp].idm[im]].updateval (ipp,im,ido);
05533 break;
05534 }
05535 case modcamclaymat:{
05536 cclay[ip[ipp].idm[im]].updateval (ipp,ido);
05537 break;
05538 }
05539 case modcamclaycoupmat:{
05540 cclayc[ip[ipp].idm[im]].updateval (ipp,ido);
05541 break;
05542 }
05543 case hypoplastmat:{
05544 hypopl[ip[ipp].idm[im]].updateval (ipp,im,ido);
05545 break;
05546 }
05547 case hissplasmat:{
05548 hisspl[ip[ipp].idm[im]].updateval (ipp,ido);
05549 break;
05550 }
05551 case effective_stress:
05552 updateipvalmat(ipp, im+1, ido);
05553 break;
05554 case microplaneM4:{
05555 mpM4[ip[ipp].idm[im]].updateval (ipp,im,ido);
05556 }
05557 case microsimp:{
05558 mpSIM[ip[ipp].idm[im]].updateval (ipp,im,ido);
05559 }
05560 case microfibro:{
05561 mpfib[ip[ipp].idm[im]].updateval (ipp,im,ido);
05562 }
05563 case shefpl:{
05564 shpl[ip[ipp].idm[im]].updateval (ipp,ido);
05565 break;
05566 }
05567 case glasgowdamage:{
05568 glasgdam[ip[ipp].idm[im]].updateval (ipp,im,ido);
05569 break;
05570 }
05571 case creep_damage:{
05572 crdam[ip[ipp].idm[im]].updateval (ipp, im, ido);
05573 break;
05574 }
05575 case shrinkagemat:{
05576 shmat[ip[ipp].idm[im]].updateval (ipp, im, ido);
05577 break;
05578 }
05579 case time_switchmat:{
05580 tswmat[ip[ipp].idm[im]].updateval (ipp, im, ido);
05581 break;
05582 }
05583 case graphm:{
05584 break;
05585 }
05586 case varelisomat :
05587 break;
05588 case geoelast:{
05589 geoel[ip[ipp].idm[im]].updateval (ipp,im,ido);
05590 break;
05591 }
05592 case scaldamage:{
05593 scdam[ip[ipp].idm[im]].updateval (ipp,im,ido);
05594 break;
05595 }
05596 case scaldamagecc:{
05597 scdamcc[ip[ipp].idm[im]].updateval (ipp,im,ido);
05598 break;
05599 }
05600 case ortodamage:{
05601 ortdam[ip[ipp].idm[im]].updateval (ipp,im,ido);
05602 break;
05603 }
05604 case ortodamagerot:{
05605 ortdamrot[ip[ipp].idm[im]].updateval (ipp,im,ido);
05606 break;
05607 }
05608 case anisodamage:{
05609 anidam[ip[ipp].idm[im]].updateval (ipp,im,ido);
05610 break;
05611 }
05612 case anisodamagerot:{
05613 anidamrot[ip[ipp].idm[im]].updateval (ipp,im,ido);
05614 break;
05615 }
05616
05617 case nonlocdamgmat:{
05618 updateipvalmat (ipp,im+1,ido);
05619 break;
05620 }
05621 case nonlocplastmat:{
05622 updateipvalmat (ipp,im+1,ido);
05623 break;
05624 }
05625 case nonlocalmicroM4:{
05626 updateipvalmat (ipp,im+1,ido);
05627 break;
05628 }
05629
05630 case damage_plasticity:{
05631 dampl[ip[ipp].idm[im]].updateval (ipp);
05632 break;
05633 }
05634 case simvisc:
05635
05636 case lemaitr:
05637 case viscoplasticity:
05638 visplas[ip[ipp].idm[im]].updateval (ipp,im,ido);
05639 break;
05640 case creeprs:
05641 case creepb3:
05642 case creepdpl:{
05643 creep_updateval (ipp,im,ido);
05644 break;
05645 }
05646 case creepbaz:
05647 crbaz[ip[ipp].idm[im]].updateval ();
05648 break;
05649 case creepeffym:
05650 creffym[ip[ipp].idm[im]].updateval (ipp,im,ido);
05651 break;
05652 case lenjonespot:
05653 break;
05654
05655 case contmat:
05656 conmat[ip[ipp].idm[im]].updateval (ipp,im,ido);
05657 break;
05658
05659 case cebfipcontmat:
05660 cebfipconmat[ip[ipp].idm[im]].updateval (ipp,im,ido);
05661 break;
05662
05663
05664 case consolidation:{
05665 csol[ip[ipp].idm[im]].updateval (ipp,im,ido);
05666 break;
05667 }
05668 default:{
05669 print_err ("unknown material type is required",__FILE__,__LINE__,__func__);
05670 break;
05671 }
05672 }
05673 }
05674
05675
05676
05677
05678
05679
05680
05681
05682
05683
05684 void mechmat::clean_ip ()
05685 {
05686 long i;
05687 for (i=0;i<tnip;i++){
05688 ip[i].clean (Mb->nlc);
05689 }
05690 }
05691
05692
05693
05694
05695
05696
05697
05698
05699
05700
05701
05702
05703
05704
05705
05706
05707 double mechmat::dstep_red_mat(long im, long ido)
05708 {
05709 long i;
05710 double ret = 1.0;
05711 double dstep_r;
05712
05713 for (i=0; i<tnip; i++)
05714 {
05715 dstep_r = dstep_red_ip(i, 0, 0);
05716 if (dstep_r < ret)
05717 ret = dstep_r;
05718 }
05719 return ret;
05720 }
05721
05722
05723
05724
05725
05726
05727
05728
05729
05730
05731
05732
05733
05734
05735
05736
05737 double mechmat::dstep_red_ip(long ipp, long im, long ido)
05738 {
05739 double ret = 1.0;
05740
05741 switch (ip[ipp].tm[im])
05742 {
05743 case elisomat :
05744 case simplas1d:
05745 case jflow:
05746 case jflow2:
05747 case mohrcoul:
05748 case mohrcoulparab:
05749 case boermaterial:
05750 case druckerprager:
05751 case druckerprager2:
05752 case chenplast:
05753 case modcamclaymat:
05754 case modcamclaycoupmat:
05755 case hissplasmat:
05756 case microplaneM4:
05757 case microsimp:
05758 case microfibro:
05759 case shefpl:
05760 case glasgowdamage:
05761 case creep_damage:
05762 case shrinkagemat:
05763 case time_switchmat:
05764 case graphm:
05765 case varelisomat :
05766 case geoelast:
05767 case scaldamage:
05768 case scaldamagecc:
05769 case ortodamage:
05770 case ortodamagerot:
05771 case anisodamage:
05772 case anisodamagerot:
05773 case damage_plasticity:
05774 case simvisc:
05775 case lemaitr:
05776 case viscoplasticity:
05777 case creeprs:
05778 case creepb3:
05779 case creepdpl:
05780 case creepbaz:
05781 case creepeffym:
05782 case lenjonespot:
05783 case contmat:
05784 case cebfipcontmat:
05785 case consolidation:
05786 ret = 1.0;
05787 break;
05788
05789 case effective_stress:
05790 case nonlocdamgmat:
05791 case nonlocplastmat:
05792 case nonlocalmicroM4:
05793 ret = dstep_red_ip(ipp,im+1,ido);
05794 break;
05795
05796 case hypoplastmat:
05797 ret = hypopl[ip[ipp].idm[im]].dstep_red(ipp,im,ido);
05798 break;
05799 default:{
05800 print_err ("unknown material type is required",__FILE__,__LINE__,__func__);
05801 break;
05802 }
05803 }
05804 return ret;
05805 }
05806
05807
05808
05809
05810
05811
05812
05813
05814
05815
05816
05817
05818
05819
05820
05821
05822
05823
05824
05825
05826
05827
05828
05829
05830
05831
05832 double mechmat::yieldfunction (long ipp, long idpm, matrix &sig,vector &q)
05833 {
05834 double f;
05835
05836 switch (ip[ipp].tm[idpm]){
05837 case simplas1d:{
05838 f = spl1d[ip[ipp].idm[idpm]].yieldfunction (sig,q);
05839 break;
05840 }
05841 case jflow:{
05842 f = j2f[ip[ipp].idm[idpm]].yieldfunction (sig,q);
05843 break;
05844 }
05845 case jflow2:{
05846 f = j2f2[ip[ipp].idm[idpm]].yieldfunction (sig,q);
05847 break;
05848 }
05849 case mohrcoulparab:{
05850 f = mcpar[ip[ipp].idm[idpm]].yieldfunction (sig);
05851 break;
05852 }
05853 case boermaterial:{
05854 f = boerm[ip[ipp].idm[idpm]].yieldfunction (sig);
05855 break;
05856 }
05857 case druckerprager:{
05858 f = drprm[ip[ipp].idm[idpm]].yieldfunction (sig,q);
05859 break;
05860 }
05861
05862
05863
05864
05865 case chenplast:{
05866 f = chplast[ip[ipp].idm[idpm]].yieldfunction (sig,q);
05867 break;
05868 }
05869 case modcamclaymat:{
05870 f = cclay[ip[ipp].idm[idpm]].yieldfunction (sig,q);
05871 break;
05872 }
05873 case modcamclaycoupmat:{
05874 f = cclayc[ip[ipp].idm[idpm]].yieldfunction (sig,q);
05875 break;
05876 }
05877 case hissplasmat:{
05878 f = hisspl[ip[ipp].idm[idpm]].yieldfunction (sig,q);
05879 break;
05880 }
05881 default:
05882 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
05883
05884 }
05885 return f;
05886 }
05887
05888
05889
05890
05891
05892
05893
05894
05895
05896
05897
05898
05899
05900
05901
05902
05903
05904
05905
05906 void mechmat::dfdsigma (long ipp, long idpm, matrix &sig, vector &q, vector &dfds)
05907 {
05908 matrix dfdst(3,3);
05909
05910 switch (ip[ipp].tm[idpm]){
05911 case simplas1d:{
05912 spl1d[ip[ipp].idm[idpm]].dfdsigma (sig,dfdst);
05913 break;
05914 }
05915 case jflow:{
05916 j2f[ip[ipp].idm[idpm]].deryieldfsigma (sig,dfdst);
05917 break;
05918 }
05919 case jflow2:{
05920 j2f2[ip[ipp].idm[idpm]].deryieldfsigma (sig,dfdst);
05921 break;
05922 }
05923 case mohrcoulparab:{
05924 mcpar[ip[ipp].idm[idpm]].deryieldfsigma (sig,dfdst);
05925 break;
05926 }
05927 case boermaterial:{
05928 boerm[ip[ipp].idm[idpm]].deryieldfsigma (sig,dfdst);
05929 break;
05930 }
05931 case druckerprager:{
05932 drprm[ip[ipp].idm[idpm]].deryieldfsigma (sig,dfdst,q);
05933 break;
05934 }
05935 case druckerprager2:{
05936
05937 break;
05938 }
05939 case chenplast:{
05940 chplast[ip[ipp].idm[idpm]].deryieldfdsigma (sig,q,dfdst);
05941 break;
05942 }
05943 case modcamclaymat:{
05944 cclay[ip[ipp].idm[idpm]].deryieldfsigma (sig, q, dfdst);
05945 break;
05946 }
05947 case modcamclaycoupmat:{
05948 cclayc[ip[ipp].idm[idpm]].deryieldfsigma (sig, q, dfdst);
05949 break;
05950 }
05951 default:
05952 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
05953
05954 }
05955
05956
05957 tensor_vector (dfds,dfdst,ip[ipp].ssst,stress);
05958
05959 }
05960
05961
05962
05963
05964
05965
05966
05967
05968
05969
05970
05971
05972
05973
05974
05975
05976
05977
05978
05979 void mechmat::dgdsigma (long ipp, long idpm, matrix &sig, vector &q, vector &dgds)
05980 {
05981 matrix dgdst(3,3);
05982
05983 switch (ip[ipp].tm[idpm]){
05984 case simplas1d:{
05985 spl1d[ip[ipp].idm[idpm]].dfdsigma (sig,dgdst);
05986 break;
05987 }
05988 case jflow:{
05989 j2f[ip[ipp].idm[idpm]].deryieldfsigma (sig,dgdst);
05990 break;
05991 }
05992 case jflow2:{
05993 j2f2[ip[ipp].idm[idpm]].deryieldfsigma (sig,dgdst);
05994 break;
05995 }
05996 case mohrcoulparab:{
05997 mcpar[ip[ipp].idm[idpm]].derplaspotsigma (sig,dgdst);
05998 break;
05999 }
06000 case boermaterial:{
06001 boerm[ip[ipp].idm[idpm]].deryieldfsigma (sig,dgdst);
06002 break;
06003 }
06004 case druckerprager:{
06005 drprm[ip[ipp].idm[idpm]].deryieldfsigma (sig,dgdst,q);
06006 break;
06007 }
06008 case druckerprager2:{
06009
06010 break;
06011 }
06012 case chenplast:{
06013 chplast[ip[ipp].idm[idpm]].deryieldfdsigma (sig,q,dgdst);
06014 break;
06015 }
06016 case modcamclaymat:{
06017 cclay[ip[ipp].idm[idpm]].deryieldfsigma (sig, q, dgdst);
06018 break;
06019 }
06020 case modcamclaycoupmat:{
06021 cclayc[ip[ipp].idm[idpm]].deryieldfsigma (sig, q, dgdst);
06022 break;
06023 }
06024 default:
06025 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06026
06027 }
06028
06029
06030 tensor_vector (dgds,dgdst,ip[ipp].ssst,stress);
06031 }
06032
06033
06034
06035
06036
06037
06038
06039
06040
06041
06042
06043
06044
06045
06046
06047
06048
06049 void mechmat::dfdsigmadsigma (long ipp, long idpm, matrix &sig, matrix &dfdsds)
06050 {
06051 matrix dfdsdst(6,6);
06052
06053 switch (ip[ipp].tm[idpm]){
06054 case simplas1d:{
06055 spl1d[ip[ipp].idm[idpm]].dfdsigmadsigma (dfdsdst);
06056 break;
06057 }
06058 case jflow2:{
06059 j2f2[ip[ipp].idm[idpm]].deryielddsigmadsigma (sig,dfdsdst);
06060 break;
06061 }
06062 case modcamclaymat:{
06063 cclay[ip[ipp].idm[idpm]].dderyieldfsigma (dfdsdst);
06064 break;
06065 }
06066 case modcamclaycoupmat:{
06067 cclayc[ip[ipp].idm[idpm]].dderyieldfsigma (dfdsdst);
06068 break;
06069 }
06070 case chenplast:{
06071 chplast[ip[ipp].idm[idpm]].deryieldfdsigmadsigma (sig,dfdsdst);
06072 break;
06073 }
06074 default:
06075 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06076
06077 }
06078
06079
06080 tensor4_matrix (dfdsds,dfdsdst,ip[ipp].ssst);
06081
06082 }
06083
06084
06085
06086
06087
06088
06089
06090
06091
06092
06093
06094
06095
06096
06097
06098
06099
06100 void mechmat::dgdsigmadsigma (long ipp, long idpm, matrix &sig, matrix &dgdsds)
06101 {
06102 matrix dgdsdst(6,6);
06103
06104 switch (ip[ipp].tm[idpm]){
06105 case simplas1d:{
06106 spl1d[ip[ipp].idm[idpm]].dfdsigmadsigma (dgdsdst);
06107 break;
06108 }
06109 case jflow2:{
06110 j2f2[ip[ipp].idm[idpm]].deryielddsigmadsigma (sig,dgdsdst);
06111 break;
06112 }
06113 case modcamclaymat:{
06114 cclay[ip[ipp].idm[idpm]].dderyieldfsigma (dgdsdst);
06115 break;
06116 }
06117 case modcamclaycoupmat:{
06118 cclayc[ip[ipp].idm[idpm]].dderyieldfsigma (dgdsdst);
06119 break;
06120 }
06121 case chenplast:{
06122 chplast[ip[ipp].idm[idpm]].deryieldfdsigmadsigma (sig,dgdsdst);
06123 break;
06124 }
06125 default:
06126 print_err("unknown plasticity model is required", __FILE__, __LINE__, __func__);
06127
06128 }
06129
06130
06131 tensor4_matrix (dgdsds,dgdsdst,ip[ipp].ssst);
06132
06133 }
06134
06135
06136
06137
06138
06139
06140
06141
06142
06143
06144
06145
06146
06147
06148
06149
06150
06151
06152 void mechmat::dfdqpar (long ipp, long idpm,matrix &sig,vector &q,vector &dq)
06153 {
06154 switch (ip[ipp].tm[idpm]){
06155 case simplas1d:{
06156 spl1d[ip[ipp].idm[idpm]].dfdqpar (dq);
06157 break;
06158 }
06159 case jflow:{
06160 j2f[ip[ipp].idm[idpm]].deryieldfq (dq);
06161 break;
06162 }
06163 case jflow2:{
06164 j2f2[ip[ipp].idm[idpm]].deryieldfq (dq);
06165 break;
06166 }
06167 case mohrcoulparab:{
06168 break;
06169 }
06170 case boermaterial:{
06171 break;
06172 }
06173 case druckerprager:{
06174 break;
06175 }
06176 case druckerprager2:{
06177 break;
06178 }
06179 case modcamclaymat:{
06180 cclay[ip[ipp].idm[idpm]].deryieldfq (sig, dq);
06181 break;
06182 }
06183 case modcamclaycoupmat:{
06184 cclayc[ip[ipp].idm[idpm]].deryieldfq (sig, q, dq);
06185 break;
06186 }
06187 case chenplast:{
06188 chplast[ip[ipp].idm[idpm]].deryieldfdq (sig,q,dq);
06189 break;
06190 }
06191 default:
06192 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06193
06194 }
06195 }
06196
06197
06198
06199
06200
06201
06202
06203
06204
06205
06206
06207
06208
06209
06210
06211
06212
06213
06214 void mechmat::dgdqpar (long ipp, long idpm,matrix &sig,vector &q,vector &dq)
06215 {
06216 switch (ip[ipp].tm[idpm]){
06217 case simplas1d:{
06218 spl1d[ip[ipp].idm[idpm]].dfdqpar (dq);
06219 break;
06220 }
06221 case jflow:{
06222 j2f[ip[ipp].idm[idpm]].deryieldfq (dq);
06223 break;
06224 }
06225 case jflow2:{
06226 j2f2[ip[ipp].idm[idpm]].deryieldfq (dq);
06227 break;
06228 }
06229 case mohrcoulparab:{
06230 break;
06231 }
06232 case boermaterial:{
06233 break;
06234 }
06235 case druckerprager:{
06236 break;
06237 }
06238 case druckerprager2:{
06239 break;
06240 }
06241 case modcamclaymat:{
06242 cclay[ip[ipp].idm[idpm]].deryieldfq (sig, dq);
06243 break;
06244 }
06245 case modcamclaycoupmat:{
06246 cclayc[ip[ipp].idm[idpm]].deryieldfq (sig, q, dq);
06247 break;
06248 }
06249 case chenplast:{
06250 chplast[ip[ipp].idm[idpm]].deryieldfdq (sig,q,dq);
06251 break;
06252 }
06253 default:
06254 print_err ("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06255
06256 }
06257 }
06258
06259
06260
06261
06262
06263
06264
06265
06266
06267
06268
06269
06270
06271
06272
06273
06274
06275
06276
06277
06278 void mechmat::dfdsigmadq (long ipp, long idpm,matrix &sig,vector &q,matrix &dfdsdq)
06279 {
06280
06281 long nh=dfdsdq.n;
06282
06283 matrix dfdsdqt(6,nh);
06284
06285 switch (ip[ipp].tm[idpm]){
06286 case simplas1d:{
06287 spl1d[ip[ipp].idm[idpm]].dfdsigmadq (dfdsdqt);
06288 break;
06289 }
06290 case jflow2:{
06291 j2f2[ip[ipp].idm[idpm]].deryieldfdsigmadq (dfdsdqt);
06292 break;
06293 }
06294 case modcamclaymat:{
06295 cclay[ip[ipp].idm[idpm]].deryieldfdsdq (dfdsdqt);
06296 break;
06297 }
06298 case modcamclaycoupmat:{
06299 cclayc[ip[ipp].idm[idpm]].deryieldfdsdq (dfdsdqt);
06300 break;
06301 }
06302 case chenplast:{
06303 chplast[ip[ipp].idm[idpm]].deryieldfdsigmadq (sig,q,dfdsdqt);
06304 break;
06305 }
06306 default:
06307 print_err("unknown plasticity model is required", __FILE__,__LINE__,__func__);
06308
06309 }
06310
06311
06312 gentensor_matrix (dfdsdq,dfdsdqt,ip[ipp].ssst);
06313
06314 }
06315
06316
06317
06318
06319
06320
06321
06322
06323
06324
06325
06326
06327
06328
06329
06330
06331
06332
06333
06334
06335 void mechmat::dgdsigmadq (long ipp, long idpm,matrix &sig,vector &q,matrix &dgdsdq)
06336 {
06337
06338 long nh=dgdsdq.n;
06339
06340 matrix dgdsdqt(6,nh);
06341
06342 switch (ip[ipp].tm[idpm]){
06343 case simplas1d:{
06344 spl1d[ip[ipp].idm[idpm]].dfdsigmadq (dgdsdqt);
06345 break;
06346 }
06347 case jflow2:{
06348 j2f2[ip[ipp].idm[idpm]].deryieldfdsigmadq (dgdsdqt);
06349 break;
06350 }
06351 case modcamclaymat:{
06352 cclay[ip[ipp].idm[idpm]].deryieldfdsdq (dgdsdqt);
06353 break;
06354 }
06355 case modcamclaycoupmat:{
06356 cclayc[ip[ipp].idm[idpm]].deryieldfdsdq (dgdsdqt);
06357 break;
06358 }
06359 case chenplast:{
06360 chplast[ip[ipp].idm[idpm]].deryieldfdsigmadq (sig,q,dgdsdqt);
06361 break;
06362 }
06363 default:
06364 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06365
06366 }
06367
06368
06369 gentensor_matrix (dgdsdq,dgdsdqt,ip[ipp].ssst);
06370 }
06371
06372
06373
06374
06375
06376
06377
06378
06379
06380
06381
06382
06383
06384
06385
06386
06387 void mechmat::dhdsigma (long ipp,long idpm,matrix &sigt,vector &q,matrix &dhds)
06388 {
06389
06390 long nh=dhds.m;
06391
06392 matrix dhdst(6,nh);
06393
06394 switch (ip[ipp].tm[idpm]){
06395 case simplas1d:{
06396
06397 break;
06398 }
06399 case chenplast:{
06400 chplast[ip[ipp].idm[idpm]].dhdsigma (sigt,q,dhdst);
06401 break;
06402 }
06403 default:
06404 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06405
06406 }
06407
06408
06409 gentensor_matrix (dhds,dhdst,ip[ipp].ssst);
06410 }
06411
06412
06413
06414
06415
06416
06417
06418
06419
06420
06421
06422
06423
06424
06425 void mechmat::dhdgamma (long ipp,long idpm,vector &epsp,matrix &sig,vector &q,vector &dhdc)
06426 {
06427 switch (ip[ipp].tm[idpm]){
06428 case simplas1d:{
06429 spl1d[ip[ipp].idm[idpm]].dhdgamma (ipp,epsp,sig,dhdc);
06430 break;
06431 }
06432 case chenplast:{
06433 chplast[ip[ipp].idm[idpm]].dhdgamma (dhdc);
06434 break;
06435 }
06436 default:
06437 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06438
06439 }
06440 }
06441
06442
06443
06444
06445
06446
06447
06448
06449
06450
06451
06452
06453
06454
06455
06456
06457
06458
06459 void mechmat::dhdqpar (long ipp,long idpm,long ido,matrix &sigt,double dgamma, vector &q,matrix &dhdq)
06460 {
06461 switch (ip[ipp].tm[idpm]){
06462 case simplas1d:{
06463
06464 break;
06465 }
06466 case chenplast:{
06467 chplast[ip[ipp].idm[idpm]].dhdqpar (sigt,q,dhdq);
06468 break;
06469 }
06470 case modcamclaymat:{
06471 cclay[ip[ipp].idm[idpm]].dhardfdq (ipp, ido, dgamma, q, dhdq);
06472 break;
06473 }
06474 default:
06475 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06476
06477 }
06478 }
06479
06480
06481
06482
06483
06484
06485
06486
06487
06488
06489
06490
06491
06492
06493
06494
06495 void mechmat::hardsoftfunction (long ipp,long idpm,matrix &sigt,vector &q,vector &h)
06496 {
06497 switch (ip[ipp].tm[idpm]){
06498 case chenplast:{
06499 chplast[ip[ipp].idm[idpm]].hvalues (sigt,q,h);
06500 break;
06501 }
06502 default:
06503 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06504
06505 }
06506 }
06507
06508
06509
06510
06511
06512
06513
06514
06515
06516
06517
06518
06519
06520
06521
06522
06523
06524 void mechmat::plasticmoduli (long ipp,long idpm,vector &epsp,matrix &sig,matrix &h)
06525 {
06526 switch (ip[ipp].tm[idpm]){
06527 case simplas1d:{
06528 spl1d[ip[ipp].idm[idpm]].plasmod (h);
06529 break;
06530 }
06531 case jflow:{
06532 j2f[ip[ipp].idm[idpm]].plasmod (h);
06533 break;
06534 }
06535 case jflow2:{
06536 j2f2[ip[ipp].idm[idpm]].plasmod (h);
06537 break;
06538 }
06539 case chenplast:{
06540
06541 break;
06542 }
06543 default:
06544 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06545
06546 }
06547
06548 }
06549
06550
06551
06552
06553
06554
06555
06556
06557
06558
06559
06560
06561
06562
06563
06564
06565
06566
06567
06568
06569
06570
06571
06572
06573
06574
06575
06576
06577
06578
06579
06580
06581
06582
06583
06584
06585
06586
06587
06588
06589
06590
06591
06592
06593
06594 void mechmat::plasmodscalar(double &denom,long ipp, long idpm, long ido, matrix &sig, vector &eps, vector &epsp, vector &qtr,double gamma)
06595 {
06596 double ret = 0.0;
06597 switch (ip[ipp].tm[idpm]){
06598 case simplas1d:{
06599
06600 ret = spl1d[ip[ipp].idm[idpm]].plasmodscalar (sig,epsp,qtr,gamma);
06601 break;
06602 }
06603 case jflow:{
06604 ret = j2f[ip[ipp].idm[idpm]].plasmodscalar (qtr);
06605 break;
06606 }
06607 case jflow2:{
06608 ret = j2f2[ip[ipp].idm[idpm]].plasmodscalar (qtr);
06609 break;
06610 }
06611 case mohrcoulparab:{
06612 break;
06613 }
06614 case boermaterial:{
06615 break;
06616 }
06617 case druckerprager:{
06618 ret = drprm[ip[ipp].idm[idpm]].plasmodscalar(qtr);
06619 break;
06620 }
06621 case druckerprager2:{
06622
06623 break;
06624 }
06625 case chenplast:{
06626
06627 if (ret<1.0e-6)
06628 denom*=2.0;
06629
06630
06631
06632
06633
06634 return;
06635
06636 }
06637 case modcamclaymat:{
06638 ret = cclay[ip[ipp].idm[idpm]].plasmodscalar(ipp, ido, sig, qtr);
06639 break;
06640 }
06641 case modcamclaycoupmat:{
06642 ret = cclayc[ip[ipp].idm[idpm]].plasmodscalar(ipp, ido, sig, qtr);
06643 break;
06644 }
06645 default:
06646 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06647
06648 }
06649 denom = ret;
06650
06651 }
06652
06653
06654
06655
06656
06657
06658
06659
06660
06661
06662
06663
06664
06665 long mechmat::give_num_interparam (long ipp,long idpm)
06666 {
06667 long np;
06668
06669 switch (ip[ipp].tm[idpm]){
06670 case simplas1d:{
06671 np = spl1d[ip[ipp].idm[idpm]].give_num_interparam ();
06672 break;
06673 }
06674 case jflow:{
06675 np = j2f[ip[ipp].idm[idpm]].give_num_interparam ();
06676 break;
06677 }
06678 case mohrcoul:{
06679 break;
06680 }
06681 case mohrcoulparab:{
06682 break;
06683 }
06684 case boermaterial:{
06685 break;
06686 }
06687 case druckerprager:{
06688 break;
06689 }
06690 case druckerprager2:{
06691 break;
06692 }
06693 case modcamclaymat:{
06694 break;
06695 }
06696 case modcamclaycoupmat:{
06697 break;
06698 }
06699 default:
06700 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06701
06702 }
06703
06704 return np;
06705 }
06706
06707
06708
06709
06710
06711
06712
06713
06714
06715
06716
06717
06718
06719
06720
06721 void mechmat::give_interparam (long ipp,long im,long ido,vector &q)
06722 {
06723 switch (ip[ipp].tm[im]){
06724 case simplas1d:{
06725 spl1d[ip[ipp].idm[im]].give_interparam (ipp,ido,q);
06726 break;
06727 }
06728 case jflow:{
06729 j2f[ip[ipp].idm[im]].give_interparam (ipp,ido,q);
06730 break;
06731 }
06732 case mohrcoul:{
06733 break;
06734 }
06735 case mohrcoulparab:{
06736 break;
06737 }
06738 case boermaterial:{
06739 break;
06740 }
06741 case druckerprager:{
06742 break;
06743 }
06744 case druckerprager2:{
06745 break;
06746 }
06747 case modcamclaymat:{
06748 break;
06749 }
06750 case modcamclaycoupmat:{
06751 break;
06752 }
06753
06754
06755
06756
06757
06758 default:
06759 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
06760
06761 }
06762
06763 }
06764
06765
06766
06767
06768
06769
06770
06771
06772
06773
06774
06775
06776
06777
06778
06779
06780
06781
06782
06783
06784 void mechmat::updateq (long ipp, long idpm, long ido,double dgamma, vector &eps, vector &epsp, matrix &sig, vector &q)
06785 {
06786 switch (ip[ipp].tm[idpm]){
06787 case simplas1d:{
06788
06789 spl1d[ip[ipp].idm[idpm]].updateq (ipp,dgamma, epsp,q);
06790 break;
06791 }
06792 case jflow:{
06793 j2f[ip[ipp].idm[idpm]].updateq (dgamma, q);
06794 break;
06795 }
06796 case jflow2:{
06797 j2f2[ip[ipp].idm[idpm]].updateq (dgamma, q);
06798 break;
06799 }
06800 case mohrcoul:{
06801 break;
06802 }
06803 case mohrcoulparab:{
06804 break;
06805 }
06806 case boermaterial:{
06807 break;
06808 }
06809 case druckerprager:{
06810 drprm[ip[ipp].idm[idpm]].updateq(ipp, epsp, q);
06811 break;
06812 }
06813
06814
06815
06816
06817 case chenplast:{
06818
06819 break;
06820 }
06821 case modcamclaymat:{
06822 cclay[ip[ipp].idm[idpm]].updateq(ipp, ido, eps, epsp, sig, q);
06823 break;
06824 }
06825 case modcamclaycoupmat:{
06826 cclayc[ip[ipp].idm[idpm]].updateq(ipp, ido, eps, epsp, sig, q);
06827 break;
06828 }
06829 default:
06830 print_err("unknown plasticity model is required",__FILE__,__LINE__,__func__);
06831
06832 }
06833 }
06834
06835
06836
06837
06838
06839
06840
06841
06842
06843
06844
06845
06846
06847
06848
06849
06850
06851
06852
06853
06854
06855
06856
06857
06858
06859
06860 void mechmat::cutting_plane (long ipp,long im,long ido,double &gamma,vector &epsn,vector &epsp,vector &q,long ni,double err)
06861 {
06862 long i,j,ncomp=epsn.n,idem;
06863 double f,denom,denomh,dgamma,nu;
06864 vector epsa(ncomp),sig(ncomp),dfds(ncomp),dgds(ncomp);
06865 matrix d(ncomp,ncomp),sigt(3,3);
06866
06867
06868 dgamma=0.0;
06869
06870
06871 elmatstiff (d,ipp);
06872 if (ip[ipp].ssst == planestrain)
06873 {
06874 d[0][3] = d[0][1]; d[1][3] = d[1][0];
06875 d[3][0] = d[1][0]; d[3][1] = d[1][0]; d[3][3] = d[1][1];
06876 }
06877 idem = ip[ipp].gemid();
06878 if (ip[ipp].ssst == planestress)
06879 {
06880 nu = Mm->give_actual_nu(ipp);
06881 epsn[3] = -nu / (1.0 - nu) * (epsn[0]-epsp[0]+epsn[1]-epsp[1]);
06882 }
06883
06884
06885
06886
06887
06888 for (i=0;i<ni;i++){
06889
06890
06891
06892
06893 subv (epsn,epsp,epsa);
06894
06895 mxv (d,epsa,sig);
06896
06897 vector_tensor (sig,sigt,ip[ipp].ssst,stress);
06898
06899 updateq(ipp, im, ido, dgamma, epsn, epsp, sigt, q);
06900
06901 f = yieldfunction (ipp,im,sigt,q);
06902
06903 if (f<err) break;
06904 if (i==ni-1 && f>err){
06905 print_err("yield surface was not reached",__FILE__,__LINE__,__func__);
06906 fprintf (Out,"\n yield surface was not reached");
06907 }
06908
06909 dfdsigma (ipp, im, sigt, q, dfds);
06910 dgdsigma (ipp, im, sigt, q, dgds);
06911
06912 mxv (d,dgds,epsa);
06913 scprd (dfds,epsa,denom);
06914
06915 denomh=0;
06916 plasmodscalar(denomh,ipp, im, ido, sigt, epsn, epsp, q,gamma);
06917 denom += denomh;
06918
06919
06920 dgamma = f/denom;
06921
06922 gamma+=dgamma;
06923
06924 for (j=0;j<ncomp;j++){
06925 epsp[j]+=dgamma*dgds[j];
06926 }
06927 if (ip[ipp].ssst == planestress)
06928 epsp[3] = epsp[0]+epsp[1];
06929 }
06930
06931
06932 subv (epsn,epsp,epsa);
06933
06934 subv (epsn,epsp,epsa);
06935
06936 mxv (d,epsa,sig);
06937
06938 storestress (0,ipp,sig);
06939 }
06940
06941
06942
06943
06944
06945
06946
06947
06948
06949
06950
06951
06952
06953
06954
06955
06956
06957
06958
06959
06960
06961
06962
06963
06964
06965 void mechmat::mult_surf_cutting_plane (long ipp,vector &epsn,vector &epsp,vector &gamma,vector &q)
06966 {
06967 long ni,i,j,k,nas,ncomp=epsn.n,nhp=q.n,nyf=gamma.n,*stat;
06968 double err;
06969 vector epsptr(ncomp),epsa(ncomp),sig(ncomp),f(nyf),qtr(nhp);
06970 vector dq,dgamma,ff;
06971 matrix d(ncomp,ncomp),h(nhp,nhp);
06972 matrix fsig,fq,cpm,amq,am,acpm;
06973
06974
06975
06976 ni=0;
06977 err=10.0;
06978
06979
06980 copyv (epsp,epsptr);
06981 copyv (q,qtr);
06982
06983 stat = new long [nyf];
06984
06985
06986 for (i=0;i<ni;i++){
06987
06988 subv (epsn,epsptr,epsa);
06989
06990
06991 elmatstiff (d,ipp);
06992
06993
06994 mxv (d,epsa,sig);
06995
06996
06997
06998
06999
07000
07001
07002
07003
07004
07005
07006
07007
07008
07009
07010
07011 nas=0;
07012 for (j=0;j<nyf;j++){
07013 if (f[j]<err){ stat[j]=0; }
07014 else{ nas++; stat[j]=1; }
07015 }
07016 if (nas==0) break;
07017
07018
07019 allocm (nas,6,fsig); allocm (nas,6,am); allocm (nas,nas,cpm);
07020 allocm (nas,nhp,fq); allocm (nhp,nas,amq); allocm (nas,nas,acpm);
07021 allocv (nas,ff); allocv (nas,dgamma); allocv (nas,dq);
07022
07023
07024
07025
07026
07027
07028
07029
07030
07031
07032
07033
07034
07035
07036
07037
07038
07039
07040
07041
07042
07043
07044 mxm (fsig,d,am);
07045 mxmt (am,fsig,cpm);
07046
07047
07048 mxmt (h,fq,amq);
07049 mxm (fq,amq,acpm);
07050
07051 addm (cpm,acpm,cpm);
07052
07053
07054 gemp (cpm.a,dgamma.a,ff.a,nas,1,Mp->zero,1);
07055
07056
07057 k=0;
07058 for (j=0;j<nyf;j++){
07059 if (stat[j]==1){
07060 gamma[j]+=dgamma[k];
07061 k++;
07062 }
07063 }
07064
07065 mtxv (fsig,dgamma,sig);
07066 mxv (amq,dgamma,dq);
07067
07068
07069 for (j=0;j<ncomp;j++){
07070 epsptr[j]+=sig[j];
07071 }
07072
07073 for (j=0;j<nhp;j++){
07074 qtr[j]-=dq[j];
07075 }
07076
07077
07078 destrm (fsig); destrm (am); destrm (cpm);
07079 destrm (fq); destrm (amq); destrm (acpm);
07080 destrv (ff); destrv (dgamma); destrv (dq);
07081
07082 }
07083
07084 copyv (epsptr,epsp);
07085 copyv (qtr,q);
07086
07087
07088
07089 subv (epsn,epsp,epsa);
07090
07091
07092 elmatstiff (d,ipp);
07093
07094
07095 mxv (d,epsa,sig);
07096
07097 storestress (0,ipp,sig);
07098
07099 }
07100
07101
07102
07103
07104
07105
07106
07107
07108
07109
07110
07111
07112
07113
07114
07115
07116
07117
07118
07119
07120
07121
07122
07123
07124
07125 void mechmat::newton_stress_return (long ipp,long im,long ido,double &gamma,vector &epsn,vector &epsp,vector &q,long ni,double err)
07126 {
07127 long i,j,k,n,ncompstr,ncomphard;
07128 double zero,f,dgamma,norstr,norstrincr,norconspar,norhardpar,norhardparincr;
07129
07130
07131 zero=Mp->zero;
07132
07133
07134 ncompstr = epsn.n;
07135
07136 ncomphard = q.n;
07137
07138 n = ncompstr + 1 + ncomphard;
07139
07140
07141
07142
07143
07144
07145
07146
07147 double *a;
07148
07149 double *x;
07150
07151 double *y;
07152
07153 matrix d(ncompstr,ncompstr);
07154
07155 matrix c(ncompstr,ncompstr);
07156
07157 vector epse(ncompstr);
07158
07159 vector sig(ncompstr);
07160
07161 matrix sigt(3,3);
07162
07163 vector sigtrial(ncompstr);
07164
07165
07166
07167 matrix dgdsds(ncompstr,ncompstr);
07168
07169 vector dgds(ncompstr);
07170
07171 vector dfds(ncompstr);
07172
07173 matrix dgdsdq(ncompstr,ncomphard);
07174
07175 vector dfdq(ncomphard);
07176
07177 matrix dhds (ncomphard,ncompstr);
07178
07179 vector dhdc (ncomphard);
07180
07181 matrix dhdq(ncomphard,ncomphard);
07182
07183 vector hs(ncomphard);
07184
07185
07186 matrix auxm(ncompstr,ncompstr),auxm2(ncompstr,ncomphard);
07187
07188 vector auxv(ncompstr),auxv2(ncompstr);
07189
07190
07191
07192 a = new double [n*n];
07193 x = new double [n];
07194 y = new double [n];
07195
07196
07197 elmatstiff (d,ipp);
07198
07199 elmatcompl (c,ipp);
07200
07201
07202 dgamma=0.0;
07203
07204
07205 subv (epsn,epsp,epse);
07206
07207 mxv (d,epse,sig);
07208
07209
07210 if (ic){
07211 for (j=0; j < ncompstr; j++)
07212 sig[j] += ic[ipp][3+j];
07213 }
07214
07215
07216 vector_tensor (sig,sigt,ip[ipp].ssst,stress);
07217
07218
07219
07220
07221 f = yieldfunction (ipp,im,sigt,q);
07222
07223
07224 hardsoftfunction (ipp,im,sigt,q,hs);
07225
07226 if (f>err){
07227
07228
07229
07230
07231 copyv (sig,sigtrial);
07232
07233
07234 for (i=0;i<ni;i++){
07235
07236
07237
07238 nullv (a,n*n);
07239 nullv (x,n);
07240 nullv (y,n);
07241
07242
07243
07244
07245
07246 dgdsigmadsigma (ipp, im, sigt, dgdsds);
07247
07248 mxm (d,dgdsds,auxm);
07249 cmulm (dgamma,auxm);
07250
07251
07252
07253 for (j=0;j<ncompstr;j++){
07254 for (k=0;k<ncompstr;k++){
07255 a[j*n+k]=auxm[j][k];
07256 }
07257 a[j*n+j]+=1.0;
07258 }
07259
07260
07261
07262
07263 dgdsigma (ipp, im, sigt, q, dgds);
07264
07265 mxv (d,dgds,auxv);
07266
07267
07268
07269 for (j=0;j<ncompstr;j++){
07270 a[j*n+ncompstr]=auxv[j];
07271 }
07272
07273
07274
07275
07276 dfdsigma (ipp, im, sigt, q, dfds);
07277
07278
07279
07280 for (j=0;j<ncompstr;j++){
07281 a[ncompstr*n+j]=dfds[j];
07282 }
07283
07284
07285
07286
07287
07288
07289 if (ncomphard>0){
07290
07291
07292
07293
07294
07295
07296
07297
07298 for (j=0;j<ncompstr;j++){
07299 for (k=0;k<ncomphard;k++){
07300 a[j*n+ncompstr+1+k]=auxm2[j][k]*dgamma;
07301 }
07302 }
07303
07304
07305
07306
07307
07308
07309
07310
07311
07312
07313 for (j=0;j<ncomphard;j++){
07314 for (k=0;k<ncompstr;k++){
07315 a[(ncompstr+1+j)*n+k]=dhds[j][k];
07316 }
07317 }
07318
07319
07320
07321
07322
07323
07324
07325
07326 for (j=0;j<ncomphard;j++){
07327 a[ncompstr*n+ncompstr+1+j]=dfdq[j];
07328 }
07329
07330
07331
07332
07333
07334
07335
07336
07337 for (j=0;j<ncomphard;j++){
07338 a[(ncompstr+1+j)*n+ncompstr]=0.0-dhdc[j];
07339 }
07340
07341
07342
07343
07344
07345
07346
07347
07348
07349 for (j=0;j<ncomphard;j++){
07350 for (k=0;k<ncomphard;k++){
07351 a[(ncompstr+1+j)*n+ncompstr+1+k]=dhdq[j][k];
07352 }
07353 a[(ncompstr+1+j)*n+ncompstr+1+j]-=1.0;
07354 }
07355 }
07356
07357
07358
07359
07360
07361 for (j=0;j<ncompstr;j++){
07362 y[j]=sigtrial[j]-sig[j]-dgamma*auxv[j];
07363 }
07364
07365
07366 y[ncompstr]=0.0-f;
07367
07368
07369 for (j=0;j<ncomphard;j++){
07370
07371
07372 }
07373
07374
07375
07376
07377
07378
07379
07380
07381
07382
07383
07384
07385
07386
07387
07388 gemp (a,x,y,n,1,1.0e-10,1);
07389
07390
07391
07392
07393
07394
07395
07396 norstr=0.0;
07397
07398 norstrincr=0.0;
07399
07400 for (j=0;j<ncompstr;j++){
07401 norstrincr+=x[j]*x[j];
07402 norstr+=sig[j]*sig[j];
07403 }
07404
07405 norconspar=x[ncompstr]*x[ncompstr];
07406
07407
07408 norhardpar=0.0;
07409
07410 norhardparincr=0.0;
07411
07412 for (j=0;j<ncomphard;j++){
07413 norhardpar+=q[j]*q[j];
07414 }
07415 for (j=ncompstr+1;j<ncompstr+1+ncomphard;j++){
07416 norhardparincr+=x[j]*x[j];
07417 }
07418
07419
07420
07421 for (j=0;j<ncompstr;j++){
07422 sig[j]+=x[j];
07423 auxv[j]=x[j];
07424 }
07425
07426
07427 dgamma+=x[ncompstr];
07428
07429
07430 for (j=0;j<ncomphard;j++){
07431 q[j]+=x[ncompstr+1+j];
07432 }
07433
07434
07435
07436
07437
07438
07439
07440
07441
07442
07443
07444
07445 mxv (c,auxv,auxv2);
07446
07447 for (j=0;j<ncompstr;j++){
07448 epsp[j]-=auxv2[j];
07449 }
07450
07451
07452
07453 vector_tensor (sig,sigt,ip[ipp].ssst,stress);
07454
07455
07456
07457
07458 f = yieldfunction (ipp,im,sigt,q);
07459
07460
07461 hardsoftfunction (ipp,im,sigt,q,hs);
07462
07463
07464
07465
07466
07467
07468
07469
07470
07471 j=0; k=0;
07472 if (norstr>zero){
07473 j++;
07474 if (norstrincr/norstr<err) k++;
07475 }
07476 if (norhardpar>zero){
07477 j++;
07478 if (norhardparincr/norhardpar<err) k++;
07479 }
07480 if (dgamma>zero){
07481 j++;
07482 if (norconspar/dgamma<err) k++;
07483 }
07484
07485
07486 if (j==k){
07487
07488
07489
07490
07491 break;
07492 }
07493 }
07494 }
07495
07496
07497 gamma+=dgamma;
07498
07499
07500 storestress (0,ipp,sig);
07501
07502 delete [] a;
07503 delete [] x;
07504 delete [] y;
07505 }
07506
07507
07508
07509
07510
07511
07512
07513
07514
07515
07516
07517
07518
07519
07520
07521
07522
07523
07524
07525
07526
07527
07528
07529
07530
07531 void mechmat::newton_stress_return_2 (long ipp,long im,long ido,double &gamma,vector &epsn,vector &epsp,vector &q,long ni,double err)
07532 {
07533 long i,j,k,n,ncompstr,ncomphard;
07534 double zero,f,dgamma,norstr,norstrincr,norconspar,norhardpar,norhardparincr;
07535
07536
07537
07538 zero=Mp->zero;
07539
07540
07541 ncompstr = epsn.n;
07542
07543 ncomphard = q.n;
07544
07545 n = ncompstr + 1 + ncomphard;
07546
07547
07548
07549
07550
07551
07552 vector qold(ncomphard);
07553
07554
07555 double *a;
07556
07557 double *x;
07558
07559 double *y;
07560
07561 matrix d(ncompstr,ncompstr);
07562
07563 matrix c(ncompstr,ncompstr);
07564
07565 vector epse(ncompstr);
07566
07567 vector sig(ncompstr);
07568
07569 matrix sigt(3,3);
07570
07571 vector sigtrial(ncompstr);
07572
07573
07574
07575 matrix dgdsds(ncompstr,ncompstr);
07576
07577 vector dgds(ncompstr);
07578
07579 vector dfds(ncompstr);
07580
07581 matrix dgdsdq(ncompstr,ncomphard);
07582
07583 vector dfdq(ncomphard);
07584
07585 matrix dhds (ncompstr,ncomphard);
07586
07587 vector dhdc (ncomphard);
07588
07589 matrix dhdq(ncomphard,ncomphard);
07590
07591 vector hs(ncomphard);
07592
07593
07594 matrix auxm(ncompstr,ncompstr),auxm2(ncompstr,ncomphard);
07595
07596 vector auxv(ncompstr),auxv2(ncompstr);
07597
07598
07599
07600 a = new double [n*n];
07601 x = new double [n];
07602 y = new double [n];
07603
07604
07605 elmatstiff (d,ipp);
07606
07607 elmatcompl (c,ipp);
07608
07609
07610 dgamma=0.0;
07611
07612
07613 subv (epsn,epsp,epse);
07614
07615 mxv (d,epse,sig);
07616
07617
07618 if (ic){
07619 for (j=0; j < ncompstr; j++)
07620 sig[j] += ic[ipp][3+j];
07621 }
07622
07623
07624 vector_tensor (sig,sigt,ip[ipp].ssst,stress);
07625
07626 copyv (q,qold);
07627
07628
07629 f = yieldfunction (ipp,im,sigt,q);
07630
07631
07632 hardsoftfunction (ipp,im,sigt,q,hs);
07633
07634 if (f>err){
07635
07636
07637
07638
07639 copyv (sig,sigtrial);
07640
07641
07642 for (i=0;i<ni;i++){
07643
07644
07645
07646 nullv (a,n*n);
07647 nullv (x,n);
07648 nullv (y,n);
07649
07650
07651
07652
07653
07654 dgdsigmadsigma (ipp, im, sigt, dgdsds);
07655
07656 mxm (d,dgdsds,auxm);
07657 cmulm (dgamma,auxm);
07658
07659
07660
07661 for (j=0;j<ncompstr;j++){
07662 for (k=0;k<ncompstr;k++){
07663 a[j*n+k]=auxm[j][k];
07664 }
07665 a[j*n+j]+=1.0;
07666 }
07667
07668
07669
07670
07671 dgdsigma (ipp, im, sigt, q, dgds);
07672
07673 mxv (d,dgds,auxv);
07674
07675
07676
07677 for (j=0;j<ncompstr;j++){
07678 a[j*n+ncompstr]=auxv[j];
07679 }
07680
07681
07682
07683
07684 dfdsigma (ipp, im, sigt, q, dfds);
07685
07686
07687
07688 for (j=0;j<ncompstr;j++){
07689 a[ncompstr*n+j]=dfds[j];
07690 }
07691
07692
07693
07694
07695
07696
07697 if (ncomphard>0){
07698
07699
07700
07701
07702 dgdsigmadq (ipp,im,sigt,q,dgdsdq);
07703
07704 mxm (d,dgdsdq,auxm2);
07705
07706 for (j=0;j<ncompstr;j++){
07707 for (k=0;k<ncomphard;k++){
07708 a[j*n+ncompstr+1+k]=auxm2[j][k]*dgamma;
07709 }
07710 }
07711
07712
07713
07714
07715
07716
07717 dhdsigma (ipp,im,sigt,q,dhds);
07718
07719
07720
07721
07722 for (j=0;j<ncomphard;j++){
07723 for (k=0;k<ncompstr;k++){
07724 a[(ncompstr+1+j)*n+k]=dhds[k][j]*dgamma;
07725 }
07726 }
07727
07728
07729
07730
07731 dfdqpar (ipp,im,sigt,q,dfdq);
07732
07733
07734
07735 for (j=0;j<ncomphard;j++){
07736 a[ncompstr*n+ncompstr+1+j]=dfdq[j];
07737 }
07738
07739
07740
07741
07742 dhdgamma (ipp,im,epsp,sigt,q,dhdc);
07743
07744
07745
07746 for (j=0;j<ncomphard;j++){
07747 a[(ncompstr+1+j)*n+ncompstr]=dhdc[j];
07748 }
07749
07750
07751
07752
07753 dhdqpar (ipp,im,ido,sigt,dgamma,q,dhdq);
07754
07755
07756
07757 for (j=0;j<ncomphard;j++){
07758 for (k=0;k<ncomphard;k++){
07759 a[(ncompstr+1+j)*n+ncompstr+1+k]=dhdq[j][k];
07760 }
07761 a[(ncompstr+1+j)*n+ncompstr+1+j]-=1.0;
07762 }
07763
07764 }
07765
07766
07767
07768
07769 for (j=0;j<ncompstr;j++){
07770 y[j]=sigtrial[j]-sig[j]-dgamma*auxv[j];
07771 }
07772
07773
07774 y[ncompstr]=0.0-f;
07775
07776
07777 for (j=0;j<ncomphard;j++){
07778
07779 }
07780
07781
07782
07783
07784
07785
07786
07787
07788
07789
07790
07791
07792
07793
07794
07795 gemp (a,x,y,n,1,1.0e-10,1);
07796
07797
07798
07799
07800
07801
07802
07803 norstr=0.0;
07804
07805 norstrincr=0.0;
07806
07807 for (j=0;j<ncompstr;j++){
07808 norstrincr+=x[j]*x[j];
07809 norstr+=sig[j]*sig[j];
07810 }
07811
07812 norconspar=x[ncompstr]*x[ncompstr];
07813
07814
07815 norhardpar=0.0;
07816
07817 norhardparincr=0.0;
07818
07819 for (j=0;j<ncomphard;j++){
07820 norhardpar+=q[j]*q[j];
07821 }
07822 for (j=ncompstr+1;j<ncompstr+1+ncomphard;j++){
07823 norhardparincr+=x[j]*x[j];
07824 }
07825
07826
07827
07828 for (j=0;j<ncompstr;j++){
07829 sig[j]+=x[j];
07830 auxv[j]=x[j];
07831 }
07832
07833
07834 dgamma+=x[ncompstr];
07835
07836
07837 for (j=0;j<ncomphard;j++){
07838 q[j]+=x[ncompstr+1+j];
07839 }
07840
07841
07842
07843
07844
07845
07846
07847
07848
07849
07850
07851
07852 mxv (c,auxv,auxv2);
07853
07854 for (j=0;j<ncompstr;j++){
07855 epsp[j]-=auxv2[j];
07856 }
07857
07858
07859
07860 vector_tensor (sig,sigt,ip[ipp].ssst,stress);
07861
07862
07863 f = yieldfunction (ipp,im,sigt,q);
07864
07865
07866 hardsoftfunction (ipp,im,sigt,q,hs);
07867
07868
07869
07870
07871
07872
07873
07874
07875
07876 j=0; k=0;
07877 if (norstr>zero){
07878 j++;
07879 if (norstrincr/norstr<err) k++;
07880 }
07881 if (norhardpar>zero){
07882 j++;
07883 if (norhardparincr/norhardpar<err) k++;
07884 }
07885 if (dgamma>zero){
07886 j++;
07887 if (norconspar/dgamma<err) k++;
07888 }
07889
07890
07891 if (j==k){
07892
07893
07894
07895
07896 break;
07897 }
07898 }
07899 }
07900
07901 fprintf (stdout,"\n\n parametr zpevneni %lf",q[0]);
07902 fprintf (stdout,"\n parametr konzistence %lf\n",gamma);
07903
07904
07905 gamma+=dgamma;
07906
07907
07908 storestress (0,ipp,sig);
07909
07910 delete [] a;
07911 delete [] x;
07912 delete [] y;
07913 }
07914
07915
07916
07917
07918
07919
07920
07921
07922
07923
07924
07925
07926
07927
07928 void mechmat::refresh_plast (long im, long ido)
07929 {
07930 if (plast == 0)
07931 for (long i=0;i<tnip;i++)
07932 if ( give_consparam(i,im,ido) ){
07933 plast = 1;
07934 return;
07935 }
07936 }
07937
07938
07939
07940
07941
07942
07943
07944
07945
07946
07947
07948
07949
07950
07951
07952
07953
07954
07955
07956
07957
07958
07959
07960
07961
07962
07963 void mechmat::damfuncpar(long ipp, long im, vector &eps, vector &kappa)
07964 {
07965 switch (ip[ipp].tm[im])
07966 {
07967 case scaldamage:
07968 {
07969 scdam[ip[ipp].idm[im]].damfuncpar(ipp, eps, kappa);
07970 break;
07971 }
07972 case scaldamagecc:
07973 {
07974 scdamcc[ip[ipp].idm[im]].damfuncpar(ipp, eps, kappa);
07975 break;
07976 }
07977 default:
07978 print_err("unknown scalar damage model is required",__FILE__,__LINE__,__func__);
07979 }
07980 return;
07981 }
07982
07983
07984
07985
07986
07987
07988
07989
07990
07991
07992
07993
07994
07995
07996
07997
07998
07999
08000
08001 double mechmat::damfunction(long ipp, long im, long ido,vector &kappa, vector &eps, vector &omegao)
08002 {
08003 double ret = 0.0;
08004
08005 switch (ip[ipp].tm[im])
08006 {
08007 case scaldamage:
08008 {
08009 ret = scdam[ip[ipp].idm[im]].damfunction(ipp,kappa,omegao);
08010 break;
08011 }
08012 case scaldamagecc:
08013 {
08014 ret = scdamcc[ip[ipp].idm[im]].damfunction(ipp,kappa,eps);
08015 break;
08016 }
08017 default:
08018 print_err("unknown scalar damage model is required",__FILE__,__LINE__,__func__);
08019 }
08020 return ret;
08021 }
08022
08023
08024
08025
08026
08027
08028
08029
08030
08031
08032
08033
08034
08035
08036 double mechmat::epsefunction(long ipp, long im)
08037 {
08038 double ret = 0.0;
08039
08040 switch (ip[ipp].tm[im])
08041 {
08042 case scaldamage:
08043 {
08044 ret = scdam[ip[ipp].idm[im]].epsefunction(ipp);
08045 break;
08046 }
08047 case scaldamagecc:
08048 {
08049 ret = scdamcc[ip[ipp].idm[im]].epsefunction();
08050 break;
08051 }
08052 default:
08053 print_err("unknown scalar damage model in function",__FILE__,__LINE__,__func__);
08054 }
08055 return ret;
08056 }
08057
08058
08059
08060
08061
08062
08063
08064
08065
08066
08067
08068
08069
08070
08071
08072
08073
08074
08075
08076
08077 double mechmat::scal_dam_sol (long ipp,long im,long ido,vector &eps,vector &kappa, vector &omegao, vector &sigma, matrix &d)
08078 {
08079 long nk = kappa.n;
08080
08081 double epse, omega;
08082 vector updkappa(nk);
08083
08084
08085
08086 if (nk > 1)
08087 {
08088 print_err("unsupported number of componetnts of vector kappa", __FILE__, __LINE__,__func__);
08089 return 0.0;
08090 }
08091
08092 if ((Mp->matmodel == local) || ((Mm->ip[ipp].hmt & 2) == 0))
08093
08094 {
08095 damfuncpar(ipp, im, eps, updkappa);
08096 if (updkappa[0] > kappa[0])
08097 {
08098 kappa[0] = updkappa[0];
08099 epse = epsefunction(ipp,im) ;
08100 if (kappa[0] > epse)
08101 omega = damfunction(ipp, im, ido, kappa, eps, omegao);
08102 else
08103 omega = 0.0;
08104 if (omega < omegao[0])
08105 omega = omegao[0];
08106 }
08107 else
08108 omega = omegao[0];
08109
08110 mxv(d, eps, sigma);
08111 cmulv(1.0-omega, sigma);
08112 return omega;
08113 }
08114
08115 if (Mp->matmodel == nonlocal)
08116
08117 {
08118 if (Mp->nonlocphase == 1)
08119
08120 {
08121 damfuncpar(ipp, im, eps, updkappa);
08122 kappa[0] = updkappa[0];
08123 fillv(0.0, sigma);
08124 return 0.0;
08125 }
08126 if (Mp->nonlocphase == 2)
08127
08128 {
08129 updkappa[0] = ip[ipp].nonloc[0];
08130 if (updkappa[0] > kappa[0])
08131 {
08132 kappa[0] = updkappa[0];
08133 epse = epsefunction(ipp, im);
08134 if (kappa[0] > epse)
08135 omega = damfunction(ipp, im, ido, kappa, eps, omegao);
08136 else
08137 omega = 0.0;
08138 if (omega < omegao[0])
08139 omega = omegao[0];
08140 }
08141 else
08142 omega = omegao[0];
08143
08144 mxv(d, eps, sigma);
08145 cmulv(1.0-omega, sigma);
08146 return omega;
08147 }
08148 }
08149
08150
08151 return 0.0;
08152 }
08153
08154
08155
08156
08157
08158
08159
08160
08161
08162
08163
08164
08165
08166
08167
08168
08169
08170
08171
08172
08173
08174
08175
08176
08177
08178
08179
08180
08181 void mechmat::givestressincr (long lcid,long ipp,long im,long ido,long fi,vector &sig)
08182 {
08183 long i;
08184
08185 switch (ip[ipp].tm[im]){
08186 case effective_stress:{
08187 givestressincr (lcid, ipp, im+1, ido, fi, sig);
08188 break;
08189 }
08190 case shrinkagemat:{
08191 i=ip[ipp].idm[im];
08192 shmat[i].givestressincr (lcid, ipp, im, ido, fi, sig);
08193 break;
08194 }
08195 case creep_damage:{
08196 givestressincr (lcid, ipp, im+1, ido, fi, sig);
08197 break;
08198 }
08199
08200
08201
08202
08203
08204
08205 case viscoplasticity:{
08206 i=ip[ipp].idm[im];
08207 visplas[i].givestressincr (ipp,ido,fi,sig);
08208 break;
08209 }
08210 case creepb3:
08211 case creeprs:
08212 case creepdpl:{
08213 creep_givestressincr (ipp,im,ido,fi,sig);
08214 break;
08215 }
08216 case hypoplastmat:{
08217 i=ip[ipp].idm[im];
08218 hypopl[i].givestressincr(ipp, im, ido, sig);
08219 break;
08220 }
08221
08222
08223
08224
08225
08226
08227
08228
08229
08230
08231
08232
08233
08234
08235
08236
08237
08238
08239
08240
08241
08242
08243
08244
08245
08246
08247
08248
08249
08250
08251
08252
08253
08254 case elisomat:
08255 case simplas1d:
08256 case simvisc:
08257 case jflow:
08258 case jflow2:
08259 case microplaneM4:
08260 case microsimp:
08261 case microfibro:
08262 case mohrcoul:
08263 case mohrcoulparab:
08264 case boermaterial:
08265 case druckerprager:
08266 case druckerprager2:
08267 case chenplast:
08268 case modcamclaymat:
08269 case modcamclaycoupmat:
08270 case shefpl:
08271 case glasgowmechmat:
08272 case glasgowdamage:
08273 case scaldamage:
08274 case scaldamagecc:
08275 case ortodamage:
08276 case ortodamagerot:
08277 case anisodamage:
08278 case anisodamagerot:
08279 case graphm:
08280 case varelisomat:
08281 case geoelast:
08282 case creepeffym:
08283 case winklerpasternak:
08284 case nonlocdamgmat:
08285 case nonlocplastmat:
08286 case nonlocalmicroM4:
08287 case damage_plasticity:
08288 case cebfipcontmat:
08289 case contmat:
08290 if (im == 0)
08291 fillv(0.0, sig);
08292 break;
08293 default:{
08294 print_err("unknown material model is required",__FILE__,__LINE__,__func__);
08295 }
08296 }
08297 }
08298
08299
08300
08301
08302
08303
08304
08305
08306
08307
08308
08309
08310
08311
08312
08313
08314
08315
08316 void mechmat::stiff_deps_vispl (long ipp,long im,long ido,vector &sig,vector &q,double dt)
08317 {
08318 long i,ncomp,ncomph,ncompo;
08319
08320 ncomp=sig.n;
08321
08322 ncomph=q.n;
08323
08324
08325 ncompo = givencompother(ipp,im+1);
08326
08327 double f,g,dlambda,cs;
08328 vector epsp(ncomp),dq(ncomph),dh(ncomph),dfds(ncomp);
08329 matrix d(ncomp,ncomp),sigt(3,3),epspt(3,3), h(ncomph, ncomph);
08330
08331 vector_tensor (sig,sigt,ip[ipp].ssst,stress);
08332
08333 f = yieldfunction (ipp,im+2,sigt,q);
08334
08335 dfdsigma (ipp,im+2,sigt,q,dfds);
08336
08337 dfdqpar (ipp,im+2,sigt,q,dq);
08338
08339 switch (ip[ipp].tm[im+1]){
08340
08341
08342
08343
08344
08345 case simvisc:{
08346
08347 g = svis[ip[ipp].idm[im+1]].gfun (f);
08348 break;
08349 }
08350 case lemaitr:{
08351
08352 cs=ip[ipp].other[ido+3*ncomp];
08353
08354 g = lmtr[ip[ipp].idm[im+1]].gfun (f,cs);
08355 break;
08356 }
08357 default:{
08358 print_err("unknown material type is required",__FILE__,__LINE__,__func__);
08359 }
08360 }
08361
08362 dlambda = 0.0;
08363 if (g>0.0){
08364
08365 dlambda = dt*g;
08366
08367 for (i=0;i<ncomp;i++){
08368 epsp[i]=dlambda*dfds[i];
08369 }
08370
08371
08372 plasticmoduli (ipp, im+2, epsp, sigt, h);
08373 mxv (h,dq,dh);
08374 for (i=0;i<ncomph;i++){
08375 dh[i]=-dlambda*dh[i];
08376 }
08377 }
08378
08379
08380 switch (ip[ipp].tm[im+1]){
08381 case simvisc:{
08382
08383 for (i=0;i<ncomp;i++){
08384 ip[ipp].eqother[ido+ncompo+i]+=epsp[i];
08385 }
08386
08387 ip[ipp].eqother[ido+ncompo+ncomp]+=dlambda;
08388
08389 for (i=0;i<ncomph;i++){
08390 ip[ipp].eqother[ido+ncompo+ncomp+1+i]+=dh[i];
08391 }
08392
08393 break;
08394 }
08395
08396
08397
08398
08399
08400
08401
08402
08403
08404 case lemaitr:{
08405
08406 for (i=0;i<ncomp;i++){
08407 ip[ipp].eqother[ido+ncompo+i]+=epsp[i];
08408 }
08409
08410 vector_tensor (epsp,epspt,ip[ipp].ssst,strain);
08411 cumulstrain (epspt,cs);
08412 ip[ipp].eqother[ido+ncomp*2]=cs;
08413
08414
08415 ip[ipp].eqother[ido+ncompo+ncomp]+=dlambda;
08416
08417
08418 for (i=0;i<ncomph;i++){
08419 ip[ipp].eqother[ido+ncompo+ncomp+1+i]+=dh[i];
08420 }
08421 break;
08422 }
08423 default:{
08424 print_err("unknown material type is requiered",__FILE__,__LINE__,__func__);
08425 }
08426 }
08427
08428
08429 matstiff (d,ipp);
08430 mxv (d,epsp,sig);
08431
08432 }
08433
08434
08435
08436
08437
08438
08439
08440
08441
08442
08443
08444
08445
08446
08447 void mechmat::stiff_eps (long ipp)
08448 {
08449 long i,ncomp=ip[ipp].ncompstr;
08450 double time,
08451 complian,
08452 strength,
08453 shrink;
08454 vector eps(ncomp),sig(ncomp);
08455
08456 time=Mp->time/3600.0/24.0;
08457
08458 switch (ip[ipp].tm[0]){
08459 case aci:{
08460 aci78mod[ip[ipp].idm[0]].compliance (time,complian,strength,shrink);
08461 break;
08462 }
08463 case cebfip:{
08464 cebfip78mod[ip[ipp].idm[0]].compliance (time,complian,strength,shrink);
08465 break;
08466 }
08467 default:
08468 print_err("unknown material type is requiered",__FILE__,__LINE__,__func__);
08469 }
08470
08471
08472
08473 for (i=0;i<ncomp;i++){
08474 ip[ipp].strain[i]=ip[ipp].stress[i]*complian;
08475 }
08476 }
08477
08478
08479
08480
08481
08482
08483
08484
08485
08486
08487
08488
08489
08490
08491
08492
08493
08494
08495
08496
08497
08498
08499
08500 double mechmat::give_first_derivative (long ipp,double r)
08501 {
08502 double d;
08503
08504 switch (ip[ipp].tm[0]){
08505 case lenjonespot:{
08506 d = lenjon[ip[ipp].idm[0]].first_derivative (r);
08507 break;
08508 }
08509 default:{
08510 print_err("unknown material model is required",__FILE__,__LINE__,__func__);
08511 }
08512 }
08513
08514 return d;
08515 }
08516
08517
08518
08519
08520
08521
08522
08523
08524
08525
08526
08527
08528
08529 double mechmat::give_second_derivative (long ipp,double r)
08530 {
08531 double d;
08532
08533 switch (ip[ipp].tm[0]){
08534 case lenjonespot:{
08535 d = lenjon[ip[ipp].idm[0]].second_derivative (r);
08536 break;
08537 }
08538 default:{
08539 print_err("unknown material model is required",__FILE__,__LINE__,__func__);
08540 }
08541 }
08542
08543 return d;
08544 }
08545
08546
08547
08548
08549
08550
08551
08552
08553
08554
08555
08556
08557
08558
08559
08560
08561
08562
08563
08564
08565
08566
08567
08568
08569
08570
08571
08572
08573
08574
08575
08576
08577
08578
08579
08580
08581
08582
08583
08584
08585
08586
08587
08588
08589
08590
08591
08592
08593
08594
08595
08596
08597
08598
08599
08600
08601
08602
08603
08604
08605
08606
08607
08608
08609
08610
08611
08612
08613
08614
08615
08616
08617
08618
08619
08620
08621
08622
08623
08624
08625
08626
08627
08628
08629
08630
08631
08632
08633
08634
08635
08636
08637
08638
08639
08640
08641
08642
08643
08644
08645
08646
08647
08648
08649
08650
08651 void mechmat::save_intpoints_txt (FILE *aux, sel &selelems, sel *selother)
08652 {
08653 long i,j,k,n,nip,ir,ipp;
08654 sel selno;
08655 int prec = (int)Mp->hdbcont.prec;
08656
08657 selno.n=1;
08658 fprintf(aux, "\n");
08659
08660 if (nonmechq)
08661 {
08662 for (i=0; i<Mm->tnknmq; i++)
08663 fprintf(aux, "%ld ", Mm->nmqid[i]);
08664 fprintf(aux, "\n");
08665 }
08666
08667 for (i=0;i<Mt->ne;i++)
08668 {
08669 ipp = Mt->elements[i].ipp[0][0];
08670 nip = Mt->give_totnip(i);
08671 selelems.presence_id(i, ir);
08672 if (ir < 0)
08673 n = 0;
08674 else
08675 n = selother[ir].give_nselcomp(ip[ipp].ncompeqother);
08676 for (j=0; j<nip; j++)
08677 {
08678 fprintf (aux,"%ld %ld %ld %ld\n",ipp, Mb->nlc, ip[ipp].ncompstr, n);
08679 if (ir < 0)
08680 ip[ipp].save_data_txt(aux,Mb->nlc, selno);
08681 else
08682 ip[ipp].save_data_txt(aux,Mb->nlc, selother[ir]);
08683
08684
08685 if (nonmechq)
08686 {
08687 for (k=0; k<nnmq; k++)
08688 fprintf (aux,"%.*le ",prec, Mm->nonmechq[k][ipp]);
08689 fprintf(aux, "\n");
08690 }
08691 ipp++;
08692 }
08693 }
08694 }
08695
08696
08697
08698
08699
08700
08701
08702
08703
08704
08705
08706
08707
08708
08709
08710
08711
08712
08713
08714 void mechmat::save_intpoints_txt (long ni, sel &selelems, sel *selother)
08715 {
08716 long i,j,k,ir,n,nip,ipp;
08717 sel selno;
08718 selno.n=1;
08719 FILE *aux;
08720 char name[FNAMELEN+20];
08721 char emsg[FNAMELEN+100];
08722 int prec = (int)Mp->hdbcont.prec;
08723
08724 sprintf(name, "%s.%ld.strain.bac",Mp->hdbcont.hdbnames, ni);
08725 aux = fopen(name,"wt");
08726 if (aux==NULL)
08727 {
08728 sprintf(emsg, "cannot open backup file %s", name);
08729 print_err(emsg, __FILE__, __LINE__, __func__);
08730 abort();
08731 }
08732 for (i=0;i<tnip;i++)
08733 {
08734 fprintf (aux,"%ld %ld %ld\n",i, Mb->nlc, ip[i].ncompstr);
08735 for (j=0;j<Mb->nlc*ip[i].ncompstr;j++)
08736 fprintf (aux, "%.*le\n",prec , ip[i].strain[j]);
08737 }
08738 fclose(aux);
08739
08740 sprintf(name, "%s.%ld.stress.bac",Mp->hdbcont.hdbnames, ni);
08741 aux = fopen(name,"wt");
08742 if (aux==NULL)
08743 {
08744 sprintf(emsg, "cannot open backup file %s", name);
08745 print_err(emsg, __FILE__, __LINE__, __func__);
08746 abort();
08747 }
08748 for (i=0;i<tnip;i++)
08749 {
08750 fprintf (aux,"%ld %ld %ld\n",i, Mb->nlc, ip[i].ncompstr);
08751 for (j=0;j<Mb->nlc*ip[i].ncompstr;j++)
08752 fprintf (aux, "%.*le\n",prec , ip[i].stress[j]);
08753 }
08754 fclose(aux);
08755
08756 sprintf(name, "%s.%ld.other.bac",Mp->hdbcont.hdbnames, ni);
08757 aux = fopen(name,"wt");
08758 if (aux==NULL)
08759 {
08760 sprintf(emsg, "cannot open backup file %s", name);
08761 print_err(emsg, __FILE__, __LINE__, __func__);
08762 abort();
08763 }
08764 for (i=0;i<Mt->ne;i++)
08765 {
08766 ipp = Mt->elements[i].ipp[0][0];
08767 nip = Mt->give_totnip(i);
08768 selelems.presence_id(i, ir);
08769 if (ir < 0)
08770 n = 0;
08771 else
08772 n = selother[ir].give_nselcomp(ip[ipp].ncompeqother);
08773 for (j=0; j<nip; j++)
08774 {
08775 fprintf (aux,"%ld %ld\n",ipp, n);
08776 if (n==0)
08777 {
08778 ipp++;
08779 continue;
08780 }
08781 for (k=0;k<ip[ipp].ncompeqother;k++)
08782 {
08783 if (selother[ir].presence_id(k))
08784 fprintf (aux, "%.*le\n",prec , ip[ipp].eqother[k]);
08785 }
08786 ipp++;
08787 }
08788 }
08789 fclose(aux);
08790
08791
08792 if (nonmechq)
08793 {
08794 sprintf(name, "%s.%ld.trfquant.bac",Mp->hdbcont.hdbnames, ni);
08795 aux = fopen(name,"wt");
08796 if (aux==NULL)
08797 {
08798 sprintf(emsg, "cannot open backup file %s", name);
08799 print_err(emsg, __FILE__, __LINE__, __func__);
08800 abort();
08801 }
08802
08803 for (i=0; i<Mm->tnknmq; i++)
08804 fprintf(aux, "%ld ", Mm->nmqid[i]);
08805 fprintf(aux, "\n");
08806
08807 for (i=0;i<tnip;i++)
08808 {
08809 fprintf (aux,"%ld ",i);
08810 for (j=0; j<nnmq; j++)
08811 fprintf (aux,"%.*le ",prec , nonmechq[j][i]);
08812
08813 fprintf(aux, "\n");
08814 }
08815 fclose(aux);
08816 }
08817 }
08818
08819
08820
08821
08822
08823
08824
08825
08826
08827
08828
08829
08830
08831
08832
08833
08834
08835 void mechmat::save_intpoints_bin (FILE *aux, sel &selelems, sel *selother)
08836 {
08837 long i,j,k,n,nip,ir,ipp;
08838 sel selno;
08839
08840 selno.n=1;
08841
08842 if (nonmechq)
08843 fwrite (nmqid, sizeof(*nmqid), sizeof(nmqid), aux);
08844
08845 for (i=0;i<Mt->ne;i++)
08846 {
08847 ipp = Mt->elements[i].ipp[0][0];
08848 nip = Mt->give_totnip(i);
08849 selelems.presence_id(i, ir);
08850
08851 if (ir < 0)
08852 n = 0;
08853 else
08854 n = selother[ir].give_nselcomp(ip[ipp].ncompeqother);
08855 for (j=0; j<nip; j++)
08856 {
08857 fwrite(&ipp, sizeof(ipp), 1, aux);
08858 fwrite(&Mb->nlc, sizeof(Mb->nlc), 1, aux);
08859 fwrite(&ip[ipp].ncompstr, sizeof(ip[ipp].ncompstr), 1, aux);
08860 fwrite(&n, sizeof(n), 1, aux);
08861 if (ir < 0)
08862 ip[ipp].save_data_bin(aux,Mb->nlc, selno);
08863 else
08864 ip[ipp].save_data_bin(aux,Mb->nlc, selother[ir]);
08865
08866
08867 if (nonmechq)
08868 {
08869 for (k=0; k<nnmq; k++)
08870 fwrite (&Mm->nonmechq[k][ipp], sizeof(Mm->nonmechq[k][ipp]), 1, aux);
08871 }
08872 ipp++;
08873 }
08874 }
08875 }
08876
08877
08878
08879
08880
08881
08882
08883
08884
08885
08886
08887
08888
08889
08890
08891
08892
08893 void mechmat::save_intpoints_bin (long ni, sel &selelems, sel *selother)
08894 {
08895 long i,j,k,ir,n,nip,ipp;
08896 sel selno;
08897 selno.n=1;
08898 FILE *aux;
08899 char name[FNAMELEN+20];
08900 char emsg[FNAMELEN+100];
08901
08902 sprintf(name, "%s.%ld.strain.bac",Mp->hdbcont.hdbnames, ni);
08903 aux = fopen(name,"wb");
08904 if (aux==NULL)
08905 {
08906 sprintf(emsg, "cannot open backup file %s", name);
08907 print_err(emsg, __FILE__, __LINE__, __func__);
08908 abort();
08909 }
08910 for (i=0;i<tnip;i++){
08911 fwrite(&i, sizeof(i), 1, aux);
08912 fwrite(&Mb->nlc, sizeof(Mb->nlc), 1, aux);
08913 fwrite(&ip[i].ncompstr, sizeof(ip[i].ncompstr), 1, aux);
08914 fwrite (ip[i].strain, sizeof(*ip[i].strain),Mb->nlc*ip[i].ncompstr, aux);
08915 }
08916 fclose(aux);
08917
08918 sprintf(name, "%s.%ld.stress.bac",Mp->hdbcont.hdbnames, ni);
08919 aux = fopen(name,"wb");
08920 if (aux==NULL)
08921 {
08922 sprintf(emsg, "cannot open backup file %s", name);
08923 print_err(emsg, __FILE__, __LINE__, __func__);
08924 abort();
08925 }
08926 for (i=0;i<tnip;i++)
08927 {
08928 fwrite(&i, sizeof(i), 1, aux);
08929 fwrite(&Mb->nlc, sizeof(Mb->nlc), 1, aux);
08930 fwrite(&ip[i].ncompstr, sizeof(ip[i].ncompstr), 1, aux);
08931 fwrite(ip[i].stress, sizeof(*ip[i].stress),Mb->nlc*ip[i].ncompstr, aux);
08932 }
08933 fclose(aux);
08934
08935 sprintf(name, "%s.%ld.other.bac",Mp->hdbcont.hdbnames, ni);
08936 aux = fopen(name,"wb");
08937 if (aux==NULL)
08938 {
08939 sprintf(emsg, "cannot open backup file %s", name);
08940 print_err(emsg, __FILE__, __LINE__, __func__);
08941 abort();
08942 }
08943 for (i=0;i<Mt->ne;i++)
08944 {
08945 ipp = Mt->elements[i].ipp[0][0];
08946 nip = Mt->give_totnip(i);
08947 selelems.presence_id(i, ir);
08948 if (ir < 0)
08949 n = 0;
08950 else
08951 n = selother[ir].give_nselcomp(ip[ipp].ncompeqother);
08952 for (j=0; j<nip; j++)
08953 {
08954 fwrite(&ipp, sizeof(ipp), 1, aux);
08955 fwrite(&n, sizeof(n), 1, aux);
08956 if (n==0)
08957 {
08958 ipp++;
08959 continue;
08960 }
08961 for (k=0;k<ip[ipp].ncompeqother;k++)
08962 {
08963 if (selother[ir].presence_id(k))
08964 fwrite (ip[ipp].eqother+k, sizeof(*ip[ipp].eqother), 1, aux);
08965 }
08966 ipp++;
08967 }
08968 }
08969 fclose(aux);
08970
08971 if (nonmechq)
08972 {
08973
08974 sprintf(name, "%s.%ld.trfquant.bac",Mp->hdbcont.hdbnames, ni);
08975 aux = fopen(name,"wb");
08976 if (aux==NULL)
08977 {
08978 sprintf(emsg, "cannot open backup file %s", name);
08979 print_err(emsg, __FILE__, __LINE__, __func__);
08980 abort();
08981 }
08982
08983 fwrite (nmqid, sizeof(*nmqid), sizeof(nmqid), aux);
08984
08985 for (ipp=0;ipp<tnip;ipp++)
08986 {
08987 fwrite(&ipp, sizeof(ipp), 1, aux);
08988 for (j=0; j<nnmq; j++)
08989 fwrite(&Mm->nonmechq[j][ipp], sizeof(Mm->nonmechq[j][ipp]), 1, aux);
08990 }
08991 fclose(aux);
08992 }
08993 }
08994
08995
08996
08997
08998
08999
09000
09001
09002
09003
09004
09005
09006
09007
09008
09009
09010
09011
09012
09013
09014 void mechmat::restore_intpoints_txt (FILE *aux, sel &selelems, sel *selother, long **selid)
09015 {
09016 long i,j,n,ir,ipp;
09017 long nlc, ncompstr;
09018 sel selno;
09019
09020 selno.n=1;
09021
09022 if (nonmechq)
09023 {
09024 for (i=0; i<Mm->tnknmq; i++)
09025 fscanf(aux, "%ld", &Mm->nmqid[i]);
09026 }
09027
09028 for (i=0;i<tnip;i++)
09029 {
09030 fscanf (aux,"%ld %ld %ld %ld",&ipp, &nlc, &ncompstr, &n);
09031 if ((ipp < 0) || (ipp >= tnip))
09032 {
09033 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09034 abort();
09035 }
09036 if (nlc != Mb->nlc)
09037 {
09038 print_err("invalid number of load cases", __FILE__, __LINE__, __func__);
09039 abort();
09040 }
09041 if (ncompstr != ip[ipp].ncompstr)
09042 {
09043 print_err("incompatible number of stress/strain components", __FILE__, __LINE__, __func__);
09044 abort();
09045 }
09046 selelems.presence_id(elip[ipp], ir);
09047 if (ir < 0)
09048 ip[ipp].restore_data_txt(aux, nlc, n, selno, NULL);
09049 else
09050 ip[ipp].restore_data_txt(aux, nlc, n, selother[ir], selid[ir]);
09051
09052
09053 if (nonmechq)
09054 {
09055 for (j=0; j<nnmq; j++)
09056 fscanf(aux,"%le", &Mm->nonmechq[j][ipp]);
09057 }
09058 }
09059 }
09060
09061
09062
09063
09064
09065
09066
09067
09068
09069
09070
09071
09072
09073
09074
09075
09076
09077
09078
09079
09080 void mechmat::restore_intpoints_txt (sel &selelems, sel *selother, long **selid)
09081 {
09082 long i,j,k,n,ir,ik,is,ipp;
09083 long nlc, ncompstr;
09084 double tmp;
09085 FILE *aux;
09086 char name[FNAMELEN+20];
09087 char emsg[FNAMELEN+100];
09088
09089
09090 sprintf(name, "%s.strain.bac",Mp->hdbcont.hdbnamer);
09091 aux = fopen(name,"rt");
09092 if (aux==NULL)
09093 {
09094 sprintf(emsg, "cannot open backup file %s", name);
09095 print_err(emsg, __FILE__, __LINE__, __func__);
09096 abort();
09097 }
09098 for (i=0;i<tnip;i++){
09099 fscanf (aux,"%ld %ld %ld",&ipp, &nlc, &ncompstr);
09100 if ((ipp < 0) || (ipp >= tnip))
09101 {
09102 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09103 abort();
09104 }
09105 if (nlc != Mb->nlc)
09106 {
09107 print_err("invalid number of load cases", __FILE__, __LINE__, __func__);
09108 abort();
09109 }
09110 if (ncompstr != ip[ipp].ncompstr)
09111 {
09112 print_err("incompatible number of stress/strain components", __FILE__, __LINE__, __func__);
09113 abort();
09114 }
09115 for (j=0; j<Mb->nlc*ip[i].ncompstr; j++)
09116 fscanf (aux, "%le", ip[i].strain+j);
09117 }
09118
09119 sprintf(name, "%s.stress.bac",Mp->hdbcont.hdbnamer);
09120 aux = fopen(name,"rt");
09121 if (aux==NULL)
09122 {
09123 sprintf(emsg, "cannot open backup file %s", name);
09124 print_err(emsg, __FILE__, __LINE__, __func__);
09125 abort();
09126 }
09127 for (i=0;i<tnip;i++){
09128 fscanf (aux,"%ld %ld %ld",&ipp, &nlc, &ncompstr);
09129 if ((ipp < 0) || (ipp >= tnip))
09130 {
09131 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09132 abort();
09133 }
09134 if (nlc != Mb->nlc)
09135 {
09136 print_err("invalid number of load cases", __FILE__, __LINE__, __func__);
09137 abort();
09138 }
09139 if (ncompstr != ip[ipp].ncompstr)
09140 {
09141 print_err("incompatible number of stress/strain components", __FILE__, __LINE__, __func__);
09142 abort();
09143 }
09144 for (j=0;j<nlc*ncompstr;j++)
09145 fscanf (aux, "%le", ip[i].stress+j);
09146 }
09147
09148 sprintf(name, "%s.other.bac",Mp->hdbcont.hdbnamer);
09149 aux = fopen(name,"rt");
09150 if (aux==NULL)
09151 {
09152 sprintf(emsg, "cannot open backup file %s", name);
09153 print_err(emsg, __FILE__, __LINE__, __func__);
09154 abort();
09155 }
09156 for (i=0;i<tnip;i++)
09157 {
09158 fscanf (aux,"%ld %ld",&ipp, &n);
09159 if ((ipp < 0) || (ipp >= tnip))
09160 {
09161 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09162 abort();
09163 }
09164 selelems.presence_id(elip[ipp], ir);
09165 for (k=0;k<n;k++)
09166 {
09167 fscanf (aux, "%le", &tmp);
09168 if (ir < 0)
09169 continue;
09170 if (selother[ir].presence_id(k,ik))
09171 {
09172 is = selid[ir][ik]+k-selother[ir].id1[ik];
09173 if (is >= ip[ipp].ncompeqother)
09174 print_err("invalid index for eqother restoring is required", __FILE__, __LINE__, __func__);
09175 else
09176 ip[ipp].eqother[is] = tmp;
09177 }
09178 }
09179 ipp++;
09180 }
09181 fclose(aux);
09182
09183
09184 if (nonmechq)
09185 {
09186 sprintf(name, "%s.trfquant.bac",Mp->hdbcont.hdbnamer);
09187 aux = fopen(name,"rt");
09188 if (aux==NULL)
09189 {
09190 sprintf(emsg, "cannot open backup file %s", name);
09191 print_err(emsg, __FILE__, __LINE__, __func__);
09192 abort();
09193 }
09194
09195
09196 for (i=0; i<Mm->tnknmq; i++)
09197 fscanf(aux, "%ld ", &Mm->nmqid[i]);
09198
09199 for (i=0;i<tnip;i++)
09200 {
09201 fscanf (aux, "%ld", &ipp);
09202 if ((ipp < 0) || (ipp >= tnip))
09203 {
09204 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09205 abort();
09206 }
09207 for (j=0; j<nnmq; j++)
09208 fscanf (aux,"%le", &nonmechq[j][ipp]);
09209 }
09210 fclose(aux);
09211 }
09212 }
09213
09214
09215
09216
09217
09218
09219
09220
09221
09222
09223
09224
09225
09226
09227
09228
09229
09230
09231
09232
09233 void mechmat::restore_intpoints_bin (FILE *aux, sel &selelems, sel *selother, long **selid)
09234 {
09235 long i, j, ipp, nlc, ncompstr, n;
09236 long ir;
09237 sel selno;
09238
09239 if (nonmechq)
09240 fread (nmqid, sizeof(*nmqid), sizeof(nmqid), aux);
09241
09242 for (i=0;i<tnip;i++){
09243 fread(&ipp, sizeof(ipp), 1, aux);
09244 fread(&nlc, sizeof(nlc), 1, aux);
09245 fread(&ncompstr, sizeof(ncompstr), 1, aux);
09246 fread(&n, sizeof(n), 1, aux);
09247 if ((ipp < 0) || (ipp >= tnip))
09248 {
09249 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09250 abort();
09251 }
09252 if (nlc != Mb->nlc)
09253 {
09254 print_err("invalid number of load cases", __FILE__, __LINE__, __func__);
09255 abort();
09256 }
09257 if (ncompstr != ip[ipp].ncompstr)
09258 {
09259 print_err("incompatible number of stress/strain components", __FILE__, __LINE__, __func__);
09260 abort();
09261 }
09262 selelems.presence_id(elip[ipp], ir);
09263 if (ir < 0)
09264 ip[ipp].restore_data_bin(aux, Mb->nlc, n, selno, NULL);
09265 else
09266 ip[ipp].restore_data_bin(aux, Mb->nlc, n, selother[ir], selid[ir]);
09267
09268
09269 if (nonmechq)
09270 {
09271 for (j=0; j<nnmq; j++)
09272 fread(&Mm->nonmechq[j][ipp], sizeof(Mm->nonmechq[j][ipp]), 1, aux);
09273 }
09274 }
09275 }
09276
09277
09278
09279
09280
09281
09282
09283
09284
09285
09286
09287
09288
09289
09290
09291
09292
09293
09294
09295 void mechmat::restore_intpoints_bin (sel &selelems, sel *selother, long **selid)
09296 {
09297 long i,k,n,ir,ik,is,ipp;
09298 long nlc, ncompstr;
09299 double tmp;
09300 FILE *aux;
09301 char name[FNAMELEN+20];
09302 char emsg[FNAMELEN+100];
09303
09304
09305 sprintf(name, "%s.strain.bac",Mp->hdbcont.hdbnamer);
09306 aux = fopen(name,"rb");
09307 if (aux==NULL)
09308 {
09309 sprintf(emsg, "cannot open backup file %s", name);
09310 print_err(emsg, __FILE__, __LINE__, __func__);
09311 abort();
09312 }
09313 for (i=0;i<tnip;i++){
09314 fread(&ipp, sizeof(ipp), 1, aux);
09315 fread(&nlc, sizeof(nlc), 1, aux);
09316 fread(&ncompstr, sizeof(ncompstr), 1, aux);
09317 if ((ipp < 0) || (ipp >= tnip))
09318 {
09319 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09320 abort();
09321 }
09322 if (nlc != Mb->nlc)
09323 {
09324 print_err("invalid number of load cases", __FILE__, __LINE__, __func__);
09325 abort();
09326 }
09327 if (ncompstr != ip[ipp].ncompstr)
09328 {
09329 print_err("incompatible number of stress/strain components", __FILE__, __LINE__, __func__);
09330 abort();
09331 }
09332 fread (ip[ipp].strain, sizeof(*ip[ipp].strain),Mb->nlc*ip[ipp].ncompstr, aux);
09333 }
09334
09335 sprintf(name, "%s.stress.bac",Mp->hdbcont.hdbnamer);
09336 aux = fopen(name,"rb");
09337 if (aux==NULL)
09338 {
09339 sprintf(emsg, "cannot open backup file %s", name);
09340 print_err(emsg, __FILE__, __LINE__, __func__);
09341 abort();
09342 }
09343 for (i=0;i<tnip;i++){
09344 fread(&ipp, sizeof(ipp), 1, aux);
09345 fread(&nlc, sizeof(nlc), 1, aux);
09346 fread(&ncompstr, sizeof(ncompstr), 1, aux);
09347 if ((ipp < 0) || (ipp >= tnip))
09348 {
09349 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09350 abort();
09351 }
09352 if (nlc != Mb->nlc)
09353 {
09354 print_err("invalid number of load cases", __FILE__, __LINE__, __func__);
09355 abort();
09356 }
09357 if (ncompstr != ip[ipp].ncompstr)
09358 {
09359 print_err("incompatible number of stress/strain components", __FILE__, __LINE__, __func__);
09360 abort();
09361 }
09362 fread (ip[ipp].stress, sizeof(*ip[ipp].stress),Mb->nlc*ip[ipp].ncompstr, aux);
09363 }
09364
09365 sprintf(name, "%s.other.bac",Mp->hdbcont.hdbnamer);
09366 aux = fopen(name,"rb");
09367 if (aux==NULL)
09368 {
09369 sprintf(emsg, "cannot open backup file %s", name);
09370 print_err(emsg, __FILE__, __LINE__, __func__);
09371 abort();
09372 }
09373 for (i=0;i<tnip;i++){
09374 fread (&ipp, sizeof(ipp), 1, aux);
09375 fread (&n, sizeof(n), 1, aux);
09376 if ((ipp < 0) || (ipp >= tnip))
09377 {
09378 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09379 abort();
09380 }
09381 selelems.presence_id(ipp, ir);
09382 for (k=0;k<n;k++)
09383 {
09384 fread (&tmp, sizeof(tmp), 1, aux);
09385 if (ir < 0)
09386 continue;
09387 if (selother[ir].presence_id(k,ik))
09388 {
09389 is = selid[ir][ik]+k-selother[ir].id1[ik];
09390 if (is >= ip[ipp].ncompeqother)
09391 print_err("invalid index for eqother restoring is required", __FILE__, __LINE__, __func__);
09392 else
09393 ip[ipp].eqother[is] = tmp;
09394 }
09395 }
09396 }
09397 fclose(aux);
09398
09399
09400 if (nonmechq)
09401 {
09402 sprintf(name, "%s.trfquant.bac",Mp->hdbcont.hdbnamer);
09403 aux = fopen(name,"rb");
09404 if (aux==NULL)
09405 {
09406 sprintf(emsg, "cannot open backup file %s", name);
09407 print_err(emsg, __FILE__, __LINE__, __func__);
09408 abort();
09409 }
09410
09411
09412 fread (nmqid, sizeof(*nmqid), sizeof(nmqid), aux);
09413
09414 for (i=0;i<tnip;i++)
09415 {
09416 fread (&ipp, sizeof(ipp), 1, aux);
09417 if ((ipp < 0) || (ipp >= tnip))
09418 {
09419 print_err("invalid integration point number", __FILE__, __LINE__, __func__);
09420 abort();
09421 }
09422 fread(&Mm->nonmechq[k][ipp], sizeof(Mm->nonmechq[k][ipp]), 1, aux);
09423 }
09424 fclose(aux);
09425 }
09426 }