00001 #include <math.h>
00002 #include "barel2d.h"
00003 #include "global.h"
00004 #include "globmat.h"
00005 #include "genfile.h"
00006 #include "element.h"
00007 #include "node.h"
00008 #include "intpoints.h"
00009 #include "loadcase.h"
00010
00011
00012
00013 barel2d::barel2d (void)
00014 {
00015 long i;
00016
00017
00018 nne=2;
00019
00020 ndofe=4;
00021
00022 tncomp=1;
00023
00024 napfun=2;
00025
00026 ssst = bar;
00027
00028
00029 nb=1;
00030
00031
00032 ncomp = new long [nb];
00033 ncomp[0]=1;
00034
00035
00036 cncomp = new long [nb];
00037 cncomp[0]=0;
00038
00039
00040
00041 nip = new long* [nb];
00042 intordsm = new long* [nb];
00043 for (i=0;i<nb;i++){
00044 nip[i] = new long [nb];
00045 intordsm[i] = new long [nb];
00046 }
00047
00048 nip[0][0]=1;
00049
00050
00051 tnip=0;
00052 for (i=0;i<nb;i++){
00053 tnip+=nip[i][i];
00054 }
00055
00056 intordsm[0][0]=1;
00057 }
00058
00059
00060
00061 barel2d::~barel2d (void)
00062 {
00063 long i;
00064
00065 for (i=0;i<nb;i++){
00066 delete [] nip[i];
00067 delete [] intordsm[i];
00068 }
00069 delete [] nip;
00070 delete [] intordsm;
00071
00072 delete [] cncomp;
00073 delete [] ncomp;
00074 }
00075
00076
00077
00078
00079
00080
00081
00082
00083
00084
00085
00086
00087
00088
00089
00090 void barel2d::geom_matrix (matrix &gm,double s,double c,double l)
00091 {
00092 gm[0][0]=c/l*(-1.0);
00093 gm[0][1]=s/l*(-1.0);
00094 gm[0][2]=c/l;
00095 gm[0][3]=s/l;
00096 }
00097
00098
00099
00100
00101
00102
00103
00104
00105
00106
00107
00108
00109
00110
00111
00112
00113 void barel2d::geom_matrix (matrix &gm,vector &x,double &jac)
00114 {
00115 double xi=0.0,n[2];
00116 vector dx(2);
00117
00118 dx_bf_lin_1d (n);
00119 derivatives_1d (dx,jac,x,xi);
00120
00121 gm[0][0]= dx[0];
00122 gm[0][2]= dx[1];
00123 }
00124
00125
00126
00127
00128
00129
00130
00131
00132
00133
00134
00135
00136
00137 double barel2d::approx (double xi,vector &nodval)
00138 {
00139 double f;
00140 vector bf(nne);
00141
00142 bf_lin_1d (bf.a,xi);
00143
00144 scprd (bf,nodval,f);
00145
00146 return f;
00147 }
00148
00149
00150
00151
00152
00153
00154
00155
00156
00157
00158
00159
00160 void barel2d::transf_matrix (ivector &nodes,matrix &tmat)
00161 {
00162 long i;
00163 fillm (0.0,tmat);
00164
00165 for (i=0;i<tmat.m;i++){
00166 tmat[i][i]=1.0;
00167 }
00168
00169 for (i=0;i<nodes.n;i++){
00170 if (Mt->nodes[nodes[i]].transf>0){
00171 tmat[i*2+0][i*2]=Mt->nodes[nodes[i]].e1[0]; tmat[i*2+0][i*2+1]=Mt->nodes[nodes[i]].e2[0];
00172 tmat[i*2+1][i*2]=Mt->nodes[nodes[i]].e1[1]; tmat[i*2+1][i*2+1]=Mt->nodes[nodes[i]].e2[1];
00173 }
00174 }
00175
00176 }
00177
00178
00179
00180
00181
00182
00183
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194 void barel2d::tran_mat (vector &x, matrix &tran,vector &gx,vector &gy)
00195 {
00196 double c,l;
00197
00198 fillm(0.0, tran);
00199
00200 tran[0][0]=gx[1]-gx[0];
00201 tran[1][0]=gy[1]-gy[0];
00202
00203
00204 l=sqrt(tran[0][0]*tran[0][0]+tran[1][0]*tran[1][0]);
00205
00206
00207 tran[0][0]/=l;
00208 tran[1][0]/=l;
00209
00210
00211 tran[0][1]=0.0-tran[1][0];
00212 tran[1][1]=tran[0][0];
00213
00214
00215 x[0]=tran[0][0]*gx[0]+tran[1][0]*gy[0];
00216 c =tran[0][1]*gx[0]+tran[1][1]*gy[0];
00217
00218
00219 x[1]=tran[0][0]*gx[1]+tran[1][0]*gy[1];
00220 c =tran[0][1]*gx[1]+tran[1][1]*gy[1];
00221 }
00222
00223
00224
00225
00226
00227
00228
00229
00230
00231
00232
00233
00234 void barel2d::stiffness_matrix (long eid,long ri,long ci,matrix &sm)
00235 {
00236 long ipid;
00237 double c,s,l,a,jac;
00238 ivector nodes (nne);
00239 vector x(nne),y(nne);
00240 matrix gm(tncomp,ndofe),d(tncomp,tncomp),tmat (ndofe,ndofe);
00241
00242 fillm (0.0,sm);
00243
00244
00245 Mt->give_node_coord2d (x,y,eid);
00246
00247
00248 l=sqrt ((x[1]-x[0])*(x[1]-x[0])+(y[1]-y[0])*(y[1]-y[0]));
00249 if (l<Mp->zero)
00250 print_err("zero length of the truss element",__FILE__,__LINE__,__func__);
00251
00252 s=(y[1]-y[0])/l;
00253 c=(x[1]-x[0])/l;
00254
00255
00256 geom_matrix (gm,s,c,l);
00257
00258 ipid=Mt->elements[eid].ipp[ri][ci];
00259
00260 Mm->matstiff (d,ipid);
00261
00262 Mc->give_area (eid,a);
00263
00264 jac=a*l;
00265
00266 bdbj (sm.a,gm.a,d.a,jac,gm.m,gm.n);
00267 }
00268
00269
00270
00271
00272
00273
00274
00275
00276
00277
00278
00279
00280
00281
00282 void barel2d::res_stiffness_matrix (long eid,matrix &sm)
00283 {
00284 long transf;
00285 ivector nodes(nne);
00286
00287 stiffness_matrix (eid,0,0,sm);
00288
00289
00290
00291 Mt->give_elemnodes (eid,nodes);
00292 transf = Mt->locsystems (nodes);
00293 if (transf>0){
00294 matrix tmat (ndofe,ndofe);
00295 transf_matrix (nodes,tmat);
00296 glmatrixtransf (sm,tmat);
00297 }
00298 }
00299
00300
00301
00302
00303
00304
00305
00306
00307
00308
00309
00310
00311
00312 void barel2d::mass_matrix (long eid,matrix &mm)
00313 {
00314 double l,a,rho;
00315 vector x(nne),y(nne);
00316
00317 fillm (0.0,mm);
00318
00319 Mt->give_node_coord2d (x,y,eid);
00320
00321
00322 l=sqrt ((x[1]-x[0])*(x[1]-x[0])+(y[1]-y[0])*(y[1]-y[0]));
00323 if (l<Mp->zero){
00324 fprintf (stderr,"\n\n zero length of the %ld truss element",eid);
00325 fprintf (stderr,"\n in function barel2d::mass_matrix (file %s, line %d).\n",__FILE__,__LINE__);
00326 }
00327
00328
00329 Mc->give_area (eid,a);
00330
00331 Mc->give_densitye (eid,rho);
00332
00333 mm[0][0]=a*rho*l*1.0/3.0;
00334 mm[0][1]=0.0;
00335 mm[0][2]=a*rho*l*1.0/6.0;
00336 mm[0][3]=0.0;
00337
00338 mm[1][0]=0.0;
00339 mm[1][1]=a*rho*l*1.0/3.0;
00340 mm[1][2]=0.0;
00341 mm[1][3]=a*rho*l*1.0/6.0;
00342
00343 mm[2][0]=a*rho*l*1.0/6.0;
00344 mm[2][1]=0.0;
00345 mm[2][2]=a*rho*l*1.0/3.0;
00346 mm[2][3]=0.0;
00347
00348 mm[3][0]=0.0;
00349 mm[3][1]=a*rho*l*1.0/6.0;
00350 mm[3][2]=0.0;
00351 mm[3][3]=a*rho*l*1.0/3.0;
00352 }
00353
00354
00355
00356
00357
00358
00359
00360
00361
00362
00363
00364
00365 void barel2d::res_mass_matrix (long eid,matrix &mm)
00366 {
00367 vector gx(nne),gy(nne),x(nne);
00368 long transf;
00369
00370 Mt->give_node_coord2d (gx,gy,eid);
00371 giveloccoord (gx,gy,x);
00372
00373 mass_matrix (eid,mm);
00374
00375
00376
00377 ivector nodes (nne);
00378 Mt->give_elemnodes (eid,nodes);
00379 transf = Mt->locsystems (nodes);
00380 if (transf>0){
00381 matrix tmat (ndofe,ndofe);
00382 transf_matrix (nodes,tmat);
00383 glmatrixtransf (mm,tmat);
00384 }
00385 }
00386
00387
00388
00389
00390
00391
00392
00393
00394
00395
00396
00397
00398
00399
00400 void barel2d::load_matrix (long eid,matrix &lm,vector &x)
00401 {
00402 double l,a;
00403
00404 fillm (0.0,lm);
00405
00406
00407 l=sqrt ((x[1]-x[0])*(x[1]-x[0]));
00408 if (l<Mp->zero){
00409 fprintf (stderr,"\n\n zero length of the %ld truss element",eid);
00410 fprintf (stderr,"\n in function barel2d::load_matrix (file %s, line %d).\n",__FILE__,__LINE__);
00411 }
00412
00413
00414 Mc->give_area (eid,a);
00415
00416 lm[0][0]=a*l*1.0/3.0;
00417 lm[0][1]=0.0;
00418 lm[0][2]=a*l*1.0/6.0;
00419 lm[0][3]=0.0;
00420
00421 lm[1][0]=0.0;
00422 lm[1][1]=a*l*1.0/3.0;
00423 lm[1][2]=0.0;
00424 lm[1][3]=a*l*1.0/6.0;
00425
00426 lm[2][0]=a*l*1.0/6.0;
00427 lm[2][1]=0.0;
00428 lm[2][2]=a*l*1.0/3.0;
00429 lm[2][3]=0.0;
00430
00431 lm[3][0]=0.0;
00432 lm[3][1]=a*l*1.0/6.0;
00433 lm[3][2]=0.0;
00434 lm[3][3]=a*l*1.0/3.0;
00435 }
00436
00437
00438
00439
00440
00441
00442
00443
00444
00445
00446
00447
00448
00449
00450 void barel2d::res_load_matrix (long eid,matrix &lm)
00451 {
00452 long transf;
00453 ivector nodes(nne);
00454 vector x(nne),gx(nne),gy(nne);
00455 matrix tran(2,2);
00456
00457
00458 Mt->give_node_coord2d (gx,gy,eid);
00459
00460 tran_mat (x,tran,gx,gy);
00461
00462
00463 load_matrix (eid,lm,x);
00464
00465
00466 lgmatrixtransfblock(lm, tran);
00467
00468
00469
00470 Mt->give_elemnodes (eid,nodes);
00471 transf = Mt->locsystems (nodes);
00472 if (transf>0){
00473 matrix tmat (ndofe,ndofe);
00474 transf_matrix (nodes,tmat);
00475 glmatrixtransf (lm,tmat);
00476 }
00477
00478 }
00479
00480
00481
00482
00483
00484
00485
00486
00487
00488
00489
00490
00491 void barel2d::giveloccoord(vector &x, vector &y,vector &lx)
00492 {
00493 double l;
00494
00495 lx[0] = 0.0;
00496 l = sqr(x[1] - x[0]) + sqr(y[1] - y[0]);
00497 l = sqrt(l);
00498 lx[1] = l;
00499 }
00500
00501
00502
00503
00504
00505
00506
00507
00508
00509
00510
00511 void barel2d::res_ip_strains (long lcid,long eid)
00512 {
00513 long transf;
00514 vector x(nne),r(ndofe),aux(ndofe);
00515 ivector nodes(nne);
00516 matrix tmat(ndofe,ndofe);
00517
00518
00519 eldispl (lcid,eid,r.a);
00520
00521
00522
00523 Mt->give_elemnodes (eid,nodes);
00524 transf = Mt->locsystems (nodes);
00525 if (transf>0){
00526 transf_matrix (nodes,tmat);
00527 lgvectortransf (aux,r,tmat);
00528 copyv (aux,r);
00529 }
00530
00531 ip_strains (lcid,eid,0,0,r);
00532
00533 }
00534
00535
00536
00537
00538
00539
00540
00541
00542
00543
00544
00545
00546
00547
00548 void barel2d::ip_strains (long lcid,long eid,long ri,long ci,vector &r)
00549 {
00550 long ipid;
00551 double s,c,l;
00552 vector x(nne),y(nne),eps(tncomp);
00553 matrix gm(tncomp,ndofe);
00554
00555
00556 Mt->give_node_coord2d (x,y,eid);
00557
00558
00559 l=sqrt ((x[1]-x[0])*(x[1]-x[0])+(y[1]-y[0])*(y[1]-y[0]));
00560 if (l<Mp->zero){
00561 fprintf (stderr,"\n\n zero length of the %ld truss element",eid);
00562 fprintf (stderr,"\n in function barel2d::ip_strains (barel2d.cpp).\n");
00563 }
00564
00565 s=(y[1]-y[0])/l;
00566 c=(x[1]-x[0])/l;
00567
00568
00569 geom_matrix (gm,s,c,l);
00570
00571 mxv (gm,r,eps);
00572
00573 ipid=Mt->elements[eid].ipp[ri][ci];
00574 Mm->storestrain (lcid,ipid,eps);
00575 }
00576
00577
00578
00579
00580
00581
00582
00583
00584
00585
00586
00587
00588 void barel2d::nod_strains_ip (long lcid,long eid,long ri,long ci)
00589 {
00590 long i,j,ipid;
00591 ivector ipnum(nne),nod(nne);
00592 vector eps(tncomp);
00593
00594
00595
00596 ipid=Mt->elements[eid].ipp[ri][ci];
00597 nodip_bar (ipid,intordsm[0][0],ipnum);
00598
00599
00600 Mt->give_elemnodes (eid,nod);
00601
00602 for (i=0;i<nne;i++){
00603
00604 Mm->givestrain (lcid,ipnum[i],eps);
00605
00606
00607 j=nod[i];
00608 Mt->nodes[j].storestrain (lcid,0,eps);
00609 }
00610 }
00611
00612
00613
00614
00615
00616
00617
00618
00619
00620
00621
00622 void barel2d::nod_strains_comp(long lcid,long eid)
00623 {
00624 long i;
00625 double s,c,l;
00626 ivector enod(nne);
00627 vector x(nne),y(nne),eps(tncomp), r(ndofe);
00628 matrix gm(tncomp,ndofe);
00629
00630
00631 eldispl (lcid,eid,r.a);
00632
00633
00634 Mt->give_node_coord2d (x,y,eid);
00635
00636 Mt->give_elemnodes (eid, enod);
00637
00638
00639 l=sqrt ((x[1]-x[0])*(x[1]-x[0])+(y[1]-y[0])*(y[1]-y[0]));
00640 if (l<Mp->zero){
00641 fprintf (stderr,"\n\n zero length of the %ld truss element",eid);
00642 fprintf (stderr,"\n in function barel2d::ip_strains (barel2d.cpp).\n");
00643 }
00644
00645 s=(y[1]-y[0])/l;
00646 c=(x[1]-x[0])/l;
00647
00648
00649 geom_matrix (gm,s,c,l);
00650
00651 mxv (gm,r,eps);
00652
00653 for (i=0;i<nne;i++)
00654 Mt->nodes[enod[i]].storestrain (lcid,0,eps);
00655 }
00656
00657
00658
00659
00660
00661
00662
00663
00664
00665
00666
00667
00668 void barel2d::strains (long lcid,long eid,long ri,long ci)
00669 {
00670 long i,naep,ncp,sid;
00671 vector coord,eps;
00672
00673 switch (Mm->stra.tape[eid]){
00674 case nowhere:{
00675 break;
00676 }
00677 case intpts:{
00678
00679 break;
00680 }
00681 case enodes:{
00682 nod_strains_ip (lcid,eid,ri,ci);
00683 break;
00684 }
00685 case userdefined:{
00686
00687 naep = Mm->stra.give_naep (eid);
00688 ncp = Mm->stra.give_ncomp (eid);
00689 sid = Mm->stra.give_sid (eid);
00690 allocv (ncp,eps);
00691 allocv (1,coord);
00692 for (i=0;i<naep;i++){
00693 Mm->stra.give_aepcoord (sid,i,coord);
00694
00695 if (Mp->strainaver==0)
00696
00697 if (Mp->strainaver==1)
00698
00699
00700 Mm->stra.storevalues(lcid,eid,i,eps);
00701 }
00702 destrv (eps);
00703 destrv (coord);
00704 break;
00705 }
00706 default:
00707 print_err("unknown strain point is required",__FILE__,__LINE__,__func__);
00708 }
00709
00710 }
00711
00712
00713
00714
00715
00716
00717
00718
00719
00720
00721
00722 void barel2d::res_ip_stresses (long lcid,long eid)
00723 {
00724 ip_stresses (lcid,eid,0,0);
00725 }
00726
00727
00728
00729
00730
00731
00732
00733
00734
00735
00736
00737
00738
00739 void barel2d::ip_stresses (long lcid,long eid,long ri,long ci)
00740 {
00741 long ipid;
00742
00743 ipid=Mt->elements[eid].ipp[ri][ci];
00744
00745
00746 if (Mp->strcomp==1)
00747 Mm->computenlstresses (ipid);
00748 }
00749
00750
00751
00752
00753
00754
00755
00756
00757
00758
00759
00760
00761
00762 void barel2d::ip_elast_stresses (long lcid,long eid,long ri,long ci)
00763 {
00764 long i,ipid;
00765 vector eps(tncomp),sig(tncomp);
00766 matrix d(tncomp,tncomp);
00767
00768 ipid=Mt->elements[eid].ipp[ri][ci];
00769
00770 for (i=0;i<intordsm[0][0];i++){
00771
00772 Mm->matstiff (d,ipid);
00773
00774 Mm->givestrain (lcid,ipid,eps);
00775
00776 mxv (d,eps,sig);
00777
00778 Mm->storestress (lcid,ipid,sig);
00779
00780 ipid++;
00781 }
00782 }
00783
00784
00785
00786
00787
00788
00789
00790
00791
00792
00793
00794
00795 void barel2d::nod_stresses_ip (long lcid,long eid,long ri,long ci)
00796 {
00797 long i,j,ipid;
00798 ivector ipnum(nne),nod(nne);
00799 vector sig(tncomp);
00800
00801
00802
00803 ipid=Mt->elements[eid].ipp[ri][ci];
00804 nodip_bar (ipid,intordsm[0][0],ipnum);
00805
00806
00807
00808 Mt->give_elemnodes (eid,nod);
00809
00810 for (i=0;i<nne;i++){
00811
00812 Mm->givestress (lcid,ipnum[i],sig);
00813
00814
00815 j=nod[i];
00816 Mt->nodes[j].storestress (lcid,0,sig);
00817 }
00818 }
00819
00820
00821
00822
00823
00824
00825
00826
00827
00828
00829
00830 void barel2d::stresses (long lcid,long eid,long ri,long ci)
00831 {
00832 long i,naep,ncp,sid;
00833 vector coord,sig;
00834
00835 switch (Mm->stre.tape[eid]){
00836 case nowhere:{
00837 break;
00838 }
00839 case intpts:{
00840
00841 break;
00842 }
00843 case enodes:{
00844 nod_stresses_ip (lcid,eid,ri,ci);
00845 break;
00846 }
00847 case userdefined:{
00848
00849 naep = Mm->stre.give_naep (eid);
00850 ncp = Mm->stre.give_ncomp (eid);
00851 sid = Mm->stre.give_sid (eid);
00852 allocv (ncp,sig);
00853 allocv (1,coord);
00854 for (i=0;i<naep;i++){
00855 Mm->stre.give_aepcoord (sid,i,coord);
00856
00857 if (Mp->stressaver==0)
00858
00859 if (Mp->stressaver==1)
00860
00861
00862 Mm->stre.storevalues(lcid,eid,i,sig);
00863 }
00864 destrv (sig);
00865 destrv (coord);
00866 break;
00867 }
00868 default:
00869 print_err("unknown stress point is required",__FILE__,__LINE__,__func__);
00870 }
00871 }
00872
00873
00874
00875
00876
00877
00878
00879
00880
00881
00882 void barel2d::nod_eqother_ip (long eid,long ri,long ci)
00883 {
00884 long i,j,ipid,ncompo;
00885 ivector ipnum(nne),nod(nne);
00886 vector eqother;
00887
00888
00889
00890 ipid=Mt->elements[eid].ipp[ri][ci];
00891 nodip_bar (ipid,intordsm[0][0],ipnum);
00892
00893
00894
00895 Mt->give_elemnodes (eid,nod);
00896
00897 for (i=0;i<nne;i++){
00898 ncompo = Mm->givencompeqother (ipnum[i],0);
00899 allocv (ncompo,eqother);
00900 Mm->giveeqother (ipnum[i],0,ncompo,eqother.a);
00901
00902
00903 j=nod[i];
00904 Mt->nodes[j].storeother (0,ncompo,eqother);
00905
00906 destrv (eqother);
00907 }
00908 }
00909
00910
00911
00912
00913
00914
00915
00916
00917
00918
00919
00920
00921
00922
00923
00924 void barel2d::internal_forces (long lcid,long eid,long ri,long ci,vector &ifor)
00925 {
00926 integratedquant iq;
00927 iq=locstress;
00928 vector x(nne),y(nne);
00929
00930
00931 compute_nlstress (lcid,eid,ri,ci);
00932
00933 Mt->give_node_coord2d (x,y,eid);
00934
00935
00936 elem_integration (iq,lcid,eid,ri,ci,ifor,x,y);
00937
00938 }
00939
00940
00941
00942
00943
00944
00945
00946
00947
00948
00949
00950
00951
00952
00953
00954
00955 void barel2d::nonloc_internal_forces (long lcid,long eid,long ri,long ci,vector &ifor)
00956 {
00957 integratedquant iq=nonlocstress;
00958 vector x(nne),y(nne);
00959 Mt->give_node_coord2d (x,y,eid);
00960
00961
00962 compute_nonloc_nlstress (lcid,eid,ri,ci);
00963
00964
00965 elem_integration (iq,lcid,eid,ri,ci,ifor,x,y);
00966 }
00967
00968
00969
00970
00971
00972
00973
00974
00975
00976
00977
00978
00979
00980
00981
00982 void barel2d::incr_internal_forces (long lcid,long eid,long ri,long ci,vector &ifor)
00983 {
00984 integratedquant iq;
00985 iq=stressincr;
00986 vector x(nne),y(nne);
00987
00988
00989 compute_nlstressincr (lcid,eid,ri,ci);
00990
00991 Mt->give_node_coord2d (x,y,eid);
00992
00993
00994 elem_integration (iq,lcid,eid,ri,ci,ifor,x,y);
00995
00996 }
00997
00998
00999
01000
01001
01002
01003
01004
01005
01006
01007
01008
01009
01010
01011
01012 void barel2d::eigstrain_forces (long lcid,long eid,long ri,long ci,vector &nfor)
01013 {
01014 integratedquant iq;
01015 iq=eigstress;
01016 vector x(nne),y(nne);
01017
01018 Mt->give_node_coord2d (x,y,eid);
01019
01020
01021 if (Mp->eigstrcomp)
01022 compute_eigstress (lcid,eid,ri,ci);
01023
01024
01025 elem_integration (iq,lcid,eid,ri,ci,nfor,x,y);
01026 }
01027
01028
01029
01030
01031
01032
01033
01034
01035
01036
01037
01038
01039
01040
01041
01042 void barel2d::res_internal_forces (long lcid,long eid,vector &ifor)
01043 {
01044 long transf;
01045 ivector nodes(nne);
01046
01047 internal_forces (lcid,eid,0,0,ifor);
01048
01049
01050 Mt->give_elemnodes (eid,nodes);
01051
01052
01053 transf = Mt->locsystems (nodes);
01054 if (transf>0){
01055 matrix tmat (ndofe,ndofe);
01056 vector v(ndofe);
01057 transf_matrix (nodes,tmat);
01058 glvectortransf (ifor,v,tmat);
01059 copyv (v,ifor);
01060 }
01061 }
01062
01063
01064
01065
01066
01067
01068
01069
01070
01071
01072
01073
01074
01075
01076 void barel2d::res_nonloc_internal_forces (long lcid,long eid,vector &ifor)
01077 {
01078 long transf;
01079 ivector nodes(nne);
01080
01081 nonloc_internal_forces (lcid,eid,0,0,ifor);
01082
01083
01084 Mt->give_elemnodes (eid,nodes);
01085
01086
01087 transf = Mt->locsystems (nodes);
01088 if (transf>0){
01089 matrix tmat (ndofe,ndofe);
01090 vector v(ndofe);
01091 transf_matrix (nodes,tmat);
01092 glvectortransf (ifor,v,tmat);
01093 copyv (v,ifor);
01094 }
01095 }
01096
01097
01098
01099
01100
01101
01102
01103
01104
01105
01106
01107
01108
01109
01110 void barel2d::res_incr_internal_forces (long lcid,long eid,vector &ifor)
01111 {
01112 long transf;
01113 ivector nodes(nne);
01114
01115 incr_internal_forces (lcid,eid,0,0,ifor);
01116
01117
01118 Mt->give_elemnodes (eid,nodes);
01119
01120
01121 transf = Mt->locsystems (nodes);
01122 if (transf>0){
01123 matrix tmat (ndofe,ndofe);
01124 vector v(ndofe);
01125 transf_matrix (nodes,tmat);
01126 glvectortransf (ifor,v,tmat);
01127 copyv (v,ifor);
01128 }
01129 }
01130
01131
01132
01133
01134
01135
01136
01137
01138
01139
01140
01141
01142
01143
01144 void barel2d::res_eigstrain_forces (long lcid,long eid,vector &nfor)
01145 {
01146 long transf;
01147 ivector nodes(nne);
01148
01149 eigstrain_forces (lcid,eid,0,0,nfor);
01150
01151
01152 Mt->give_elemnodes (eid,nodes);
01153
01154
01155 transf = Mt->locsystems (nodes);
01156 if (transf>0){
01157 matrix tmat (ndofe,ndofe);
01158 vector v(ndofe);
01159 transf_matrix (nodes,tmat);
01160 glvectortransf (nfor,v,tmat);
01161 copyv (v,nfor);
01162 }
01163 }
01164
01165
01166
01167
01168
01169
01170
01171
01172
01173
01174
01175
01176
01177 void barel2d::intpointval (long eid,vector &nodval,vector &ipval)
01178 {
01179 long ii,jj,i,k,ipp;
01180 double xi;
01181 vector w,gp;
01182
01183 k=0;
01184 for (ii = 0; ii < nb; ii++)
01185 {
01186 for (jj = 0; jj < nb; jj++)
01187 {
01188 ipp=Mt->elements[eid].ipp[ii][jj];
01189 if (intordsm[ii][jj] == 0)
01190 continue;
01191 allocv (intordsm[ii][jj],gp);
01192 allocv (intordsm[ii][jj],w);
01193 gauss_points (gp.a,w.a,intordsm[ii][jj]);
01194 for (i = 0; i < intordsm[ii][jj]; i++)
01195 {
01196 xi=gp[i];
01197
01198 ipval[k] = approx (xi,nodval);
01199 k++;
01200 ipp++;
01201 }
01202 destrv (gp); destrv (w);
01203 }
01204 }
01205 }
01206
01207
01208
01209
01210
01211
01212
01213
01214
01215
01216
01217
01218
01219
01220
01221
01222 void barel2d::ipcoord (long eid,long ipp,long ri,long ci,vector &coord)
01223 {
01224 long i,j,ii;
01225 double xi;
01226 vector x(nne),y(nne),w(intordsm[ri][ci]),gp(intordsm[ri][ci]);
01227
01228 gauss_points (gp.a,w.a,intordsm[ri][ci]);
01229 Mt->give_node_coord2d (x,y,eid);
01230 ii=Mt->elements[eid].ipp[ri][ci];
01231
01232 for (i=0;i<intordsm[ri][ci];i++){
01233 for (j=0;j<intordsm[ri][ci];j++){
01234 xi=gp[i];
01235 if (ii==ipp){
01236 coord[0]=approx (xi,x);
01237 coord[1]=approx (xi,y);
01238 coord[2]=0.0;
01239 }
01240 ii++;
01241 }
01242 }
01243 }
01244
01245
01246
01247
01248
01249
01250
01251
01252
01253
01254
01255 void barel2d::ipvolume (long eid,long ri,long ci)
01256 {
01257 long ipp;
01258 double a,l,vol;
01259 vector x(nne),y(nne);
01260
01261 Mt->give_node_coord2d (x,y,eid);
01262
01263 l=sqrt ((x[1]-x[0])*(x[1]-x[0])+(y[1]-y[0])*(y[1]-y[0]));
01264 if (l<Mp->zero){
01265 char emsg[100];
01266 sprintf(emsg, "zero length of the %ld truss element", eid);
01267 print_err(emsg, __FILE__, __LINE__, __func__);
01268 }
01269
01270 Mc->give_area (eid,a);
01271
01272 vol=a*l;
01273
01274 ipp=Mt->elements[eid].ipp[ri][ci];
01275 Mm->storeipvol (ipp,vol);
01276 }
01277
01278
01279
01280
01281
01282
01283
01284
01285
01286
01287
01288
01289
01290
01291
01292
01293
01294
01295
01296
01297
01298
01299
01300
01301
01302 void barel2d::inicipval(long eid, long ri, long ci, matrix &nodval, inictype *ictn)
01303 {
01304 long i, j, k, ipp;
01305 long ii, jj, nv = nodval.n;
01306 double xi, ipval;
01307 vector w, gp, anv(nne);
01308 long nstra, nstre, ncompstr, ncompeqother;
01309 long idstra, idstre, idoth, idic;
01310 inictype ict;
01311
01312 nstra = idstra = nstre = idstre = idoth = idic = 0;
01313
01314 ict = ictn[0];
01315 for (i=0; i<nne; i++)
01316 {
01317 if (ictn[i] != ict)
01318 {
01319 print_err("Incompatible types of initial conditions on element %ld\n"
01320 " at %ld. and %ld. nodes", __FILE__, __LINE__, __func__, eid+1, 1, i+1);
01321 abort();
01322 }
01323 }
01324 for (j = 0; j < nv; j++)
01325 {
01326 for(i = 0; i < nne; i++)
01327 anv[i] = nodval[i][j];
01328 for (ii = 0; ii < nb; ii++)
01329 {
01330 for (jj = 0; jj < nb; jj++)
01331 {
01332 ipp=Mt->elements[eid].ipp[ri+ii][ci+jj];
01333 if (intordsm[ii][jj] == 0)
01334 continue;
01335 allocv (intordsm[ii][jj],gp);
01336 allocv (intordsm[ii][jj],w);
01337 gauss_points (gp.a,w.a,intordsm[ii][jj]);
01338 for (k = 0; k < intordsm[ii][jj]; k++)
01339 {
01340 xi=gp[k];
01341
01342 ipval = approx (xi,anv);
01343 ncompstr = Mm->ip[ipp].ncompstr;
01344 ncompeqother = Mm->ip[ipp].ncompeqother;
01345 if ((ictn[0] & inistrain) && (j < ncompstr))
01346 {
01347 Mm->ip[ipp].strain[idstra] += ipval;
01348 ipp++;
01349 continue;
01350 }
01351 if ((ictn[0] & inistress) && (j < nstra + ncompstr))
01352 {
01353 Mm->ip[ipp].stress[idstre] += ipval;
01354 ipp++;
01355 continue;
01356 }
01357 if ((ictn[0] & iniother) && (j < nstra+nstre+ncompeqother))
01358 {
01359 Mm->ip[ipp].eqother[idoth] += ipval;
01360 ipp++;
01361 continue;
01362 }
01363 if ((ictn[0] & inicond) && (j < nv))
01364 {
01365 if (Mm->ic[ipp] == NULL)
01366 {
01367 Mm->ic[ipp] = new double[nv-j];
01368 memset(Mm->ic[ipp], 0, sizeof(*Mm->ic[ipp])*(nv-j));
01369 }
01370 Mm->ic[ipp][idic] += ipval;
01371 ipp++;
01372 continue;
01373 }
01374 ipp++;
01375 }
01376 destrv (gp); destrv (w);
01377 }
01378 }
01379 ipp=Mt->elements[eid].ipp[ri][ci];
01380 ncompstr = Mm->ip[ipp].ncompstr;
01381 ncompeqother = Mm->ip[ipp].ncompeqother;
01382 if ((ictn[0] & inistrain) && (j < ncompstr))
01383 {
01384 nstra++;
01385 idstra++;
01386 continue;
01387 }
01388 if ((ictn[0] & inistress) && (j < nstra + ncompstr))
01389 {
01390 nstre++;
01391 idstre++;
01392 continue;
01393 }
01394 if ((ictn[0] & iniother) && (j < nstra + nstre + ncompeqother))
01395 {
01396 idoth++;
01397 continue;
01398 }
01399 if ((ictn[0] & inicond) && (j < nv))
01400 {
01401 idic++;
01402 continue;
01403 }
01404 }
01405 }
01406
01407
01408
01409
01410
01411
01412
01413
01414
01415
01416
01417
01418 void barel2d::compute_nlstress (long lcid,long eid,long ri,long ci)
01419 {
01420 long ipid;
01421
01422 ipid=Mt->elements[eid].ipp[ri][ci];
01423
01424
01425 if (Mp->strcomp==1)
01426 Mm->computenlstresses (ipid);
01427
01428 }
01429
01430
01431
01432
01433
01434
01435
01436
01437
01438
01439
01440
01441 void barel2d::compute_nlstressincr (long lcid,long eid,long ri,long ci)
01442 {
01443 long ipid;
01444
01445 ipid=Mt->elements[eid].ipp[ri][ci];
01446
01447
01448 if (Mp->strcomp==1)
01449 Mm->computenlstressesincr (ipid);
01450
01451 }
01452
01453
01454
01455
01456
01457
01458
01459
01460
01461
01462
01463
01464
01465 void barel2d::local_values (long lcid,long eid,long ri,long ci)
01466 {
01467 long ipp;
01468
01469 ipp=Mt->elements[eid].ipp[ri][ci];
01470
01471
01472
01473 if (Mp->strcomp==1)
01474 Mm->computenlstresses (ipp);
01475
01476 }
01477
01478
01479
01480
01481
01482
01483
01484
01485
01486
01487
01488
01489 void barel2d::compute_nonloc_nlstress (long lcid,long eid,long ri,long ci)
01490 {
01491 long ipid;
01492
01493 ipid=Mt->elements[eid].ipp[ri][ci];
01494
01495
01496 if (Mp->strcomp==1)
01497 Mm->compnonloc_nlstresses (ipid);
01498
01499 }
01500
01501
01502
01503
01504
01505
01506
01507
01508
01509
01510
01511
01512 void barel2d::compute_eigstress (long lcid,long eid,long ri,long ci)
01513 {
01514 long ipid,rmid;
01515 vector eigstr(tncomp),sig(tncomp);
01516 matrix d(tncomp,tncomp);
01517
01518
01519 ipid=Mt->elements[eid].ipp[ri][ci];
01520
01521 rmid = Mm->ip[ipid].grmid ();
01522
01523 if (rmid>-1){
01524
01525 Mm->eigenstresses (sig,ipid,rmid);
01526 }else{
01527
01528
01529 Mm->giveeigstrain (ipid,eigstr);
01530
01531
01532 Mm->matstiff (d,ipid);
01533
01534 mxv (d,eigstr,sig);
01535 }
01536
01537 Mm->storeeigstress (ipid,sig);
01538
01539 }
01540
01541
01542
01543
01544
01545
01546
01547
01548
01549
01550
01551
01552
01553
01554
01555
01556
01557 void barel2d::elem_integration (integratedquant iq,long lcid,long eid,long ri,long ci,vector &nv,vector &x,vector &y)
01558 {
01559 long ipid;
01560 double s,c,l,a;
01561 ivector nodes(nne);
01562 vector ipv(tncomp),contr(ndofe);
01563 matrix gm(tncomp,ndofe);
01564
01565
01566 l=sqrt ((x[1]-x[0])*(x[1]-x[0])+(y[1]-y[0])*(y[1]-y[0]));
01567 if (l<Mp->zero)
01568 print_err("zero length of the truss element",__FILE__,__LINE__,__func__);
01569
01570 s=(y[1]-y[0])/l;
01571 c=(x[1]-x[0])/l;
01572
01573 ipid=Mt->elements[eid].ipp[ri][ci];
01574
01575
01576 Mm->givequantity (iq,lcid,ipid,0,ipv);
01577
01578
01579 geom_matrix (gm,s,c,l);
01580
01581
01582 mtxv (gm,ipv,contr);
01583
01584
01585 Mc->give_area (eid,a);
01586
01587 cmulv (a*l,contr);
01588
01589
01590 addv(contr,nv,nv);
01591
01592 ipid++;
01593
01594 }
01595
01596
01597
01598
01599
01600
01601
01602
01603
01604
01605
01606
01607
01608
01609
01610 void barel2d::mechq_nodval (long eid,vector &nodval,nontransquant qt)
01611 {
01612 long i,ipid;
01613 ivector ipnum(nne);
01614
01615
01616
01617 ipid=Mt->elements[eid].ipp[0][0];
01618 nodip_bar (ipid,intordsm[0][0],ipnum);
01619
01620 for (i=0;i<nne;i++){
01621
01622
01623 nodval[i] = Mm->givemechq(qt, ipnum[i]);
01624
01625 }
01626 }
01627
01628
01629
01630
01631
01632
01633
01634
01635
01636
01637
01638
01639
01640
01641
01642
01643
01644
01645 void barel2d::mechq_nodval_comp (long eid, vector &nodval, long ncnv, long nq, nontransquant *qt)
01646 {
01647 long i, j, ncompstr, ncompo, ncompnl, ipid;
01648 ivector ipnum(nne), enod(nne);
01649 intpoints ipb;
01650
01651
01652 ipid=Mt->elements[eid].ipp[0][0];
01653
01654
01655 Mt->give_elemnodes(eid, enod);
01656
01657
01658
01659 nodip_bar (ipid,intordsm[0][0],ipnum);
01660
01661
01662 nod_strains_comp(0, eid);
01663
01664
01665 ncompnl = Mm->give_num_averq(ipid, Mm->givenonlocid(ipid));
01666
01667
01668
01669 ipb.copy(Mm->ip[ipid], Mb->nlc, ncompnl);
01670
01671
01672
01673 for (i=0;i<ncnv;i++)
01674 {
01675
01676 ncompstr = Mm->ip[ipnum[i]].ncompstr;
01677
01678 Mm->storestrain(0, ipid, 0, ncompstr, Mt->nodes[enod[i]].strain);
01679
01680 ncompo = Mm->ip[ipnum[i]].ncompeqother;
01681 if (ncompo){
01682
01683 Mm->storeeqother(ipid, 0, ncompo, Mm->ip[ipnum[i]].eqother);
01684 }
01685
01686 if (ncompnl){
01687
01688 Mm->storenonloc(ipid, 0, ncompnl, Mm->ip[ipnum[i]].nonloc);
01689 }
01690
01691
01692 switch(Mp->matmodel)
01693 {
01694 case local:
01695 Mm->computenlstresses(ipid);
01696 break;
01697 case nonlocal:
01698 Mm->compnonloc_nlstresses (ipid);
01699 break;
01700 default:
01701 print_err("unknown approach in material model is required (Mp->matmodel=%ld)", __FILE__, __LINE__, __func__, Mp->matmodel);
01702 }
01703
01704
01705
01706 Mm->updateipvalmat(ipid, 0, 0);
01707
01708
01709 for (j=0; j<nq; j++)
01710 nodval[j*ncnv+i] = Mm->givemechq(qt[j], ipid);
01711 }
01712
01713
01714 Mm->storestrain(0, ipid, 0, Mb->nlc*ipb.ncompstr, ipb.strain);
01715 Mm->storestress(0, ipid, 0, Mb->nlc*ipb.ncompstr, ipb.stress);
01716 Mm->storeother(ipid, 0, ipb.ncompother, ipb.other);
01717 Mm->storeeqother(ipid, 0, ipb.ncompeqother, ipb.eqother);
01718 Mm->storenonloc(ipid, 0, ncompnl, ipb.nonloc);
01719 }