00001 #include "axiquadc.h"
00002 #include "quadlineart.h"
00003 #include "axisymlq.h"
00004 #include "axisymqq.h"
00005 #include "global.h"
00006 #include "globalt.h"
00007 #include "globalc.h"
00008 #include "element.h"
00009 #include "intpoints.h"
00010 #include "globmatt.h"
00011
00012
00013 axiquadc::axiquadc (void)
00014 {
00015 long i,j;
00016
00017 if (Cp->lbb==lin_lin){
00018
00019 mnb=Asymlq->nb;
00020
00021 mndofe=Asymlq->ndofe;
00022
00023 tnmcomp=Asymlq->tncomp;
00024
00025 mncomp=Asymlq->ncomp;
00026
00027 nnemp = Asymlq->nne;
00028
00029 nnetp = Lqat->nne;
00030
00031 tndofe = Lqat->ndofe;
00032 }
00033
00034 if (Cp->lbb==quad_lin){
00035
00036 mnb=Asymqq->nb;
00037
00038 mndofe=Asymqq->ndofe;
00039
00040 tnmcomp=Asymqq->tncomp;
00041
00042 mncomp=Asymqq->ncomp;
00043
00044 nnemp = Asymqq->nne;
00045
00046 nnetp = Lqat->nne;
00047
00048 tndofe = Lqat->ndofe;
00049 }
00050
00051 if (Cp->lbb==quad_quad){
00052
00053 mnb=Asymqq->nb;
00054
00055 mndofe=Asymqq->ndofe;
00056
00057 tnmcomp=Asymqq->tncomp;
00058
00059 mncomp=Asymqq->ncomp;
00060
00061 nnemp = Asymqq->nne;
00062
00063 nnetp = Qqat->nne;
00064
00065 tndofe = Qqat->ndofe;
00066 }
00067
00068
00069 ntm=Tp->ntm;
00070
00071
00072
00073 intordvum=NULL; intordvlm=NULL;
00074 nipu=NULL; nipl=NULL;
00075
00076
00077 dofe = new long [ntm];
00078
00079 nipu = new long* [mnb];
00080
00081 nipl = new long* [ntm];
00082
00083 intordvum = new long* [mnb];
00084
00085 intordvlm = new long* [ntm];
00086
00087 for (i=0;i<mnb;i++){
00088 nipu[i] = new long [ntm];
00089 intordvum[i] = new long [ntm];
00090 }
00091 for (i=0;i<ntm;i++){
00092 nipl[i] = new long [mnb];
00093 intordvlm[i] = new long [mnb];
00094 }
00095
00096
00097 mordering = new long [mndofe];
00098
00099 for (i=0;i<mndofe;i++){
00100 mordering[i]=i+1;
00101 }
00102
00103
00104 if (Cp->lbb==lin_lin){
00105 switch (Tp->tmatt){
00106 case onemedium:{
00107 intordvum[0][0]=2; intordvlm[0][0]=2;
00108 nipu[0][0]=4; nipl[0][0]=4;
00109 dofe[0]=nnetp;
00110 break;
00111 }
00112 case twomediacoup:{
00113 intordvum[0][0]=2; intordvum[0][1]=2;
00114 intordvlm[0][0]=2; intordvlm[1][0]=2;
00115
00116 if (Cp->savemode==0){
00117 nipu[0][0]=4; nipu[0][1]=4;
00118 nipl[0][0]=4; nipl[1][0]=4;
00119 }
00120 if (Cp->savemode==1){
00121 nipu[0][0]=4; nipu[0][1]=0;
00122 nipl[0][0]=4; nipl[1][0]=0;
00123 }
00124
00125 dofe[0]=nnetp; dofe[1]=nnetp;
00126 break;
00127 }
00128 case threemediacoup:{
00129 intordvum[0][0]=2; intordvum[0][1]=2; intordvum[0][2]=2;
00130 intordvlm[0][0]=2; intordvlm[1][0]=2; intordvlm[2][0]=2;
00131
00132 if (Cp->savemode==0){
00133 nipu[0][0]=4; nipu[0][1]=4; nipu[0][2]=4;
00134 nipl[0][0]=4; nipl[1][0]=4; nipl[2][0]=4;
00135 }
00136 if (Cp->savemode==1){
00137 nipu[0][0]=4; nipu[0][1]=0; nipu[0][2]=0;
00138 nipl[0][0]=4; nipl[1][0]=0; nipl[2][0]=0;
00139 }
00140
00141 dofe[0]=nnetp; dofe[1]=nnetp; dofe[2]=nnetp;
00142 break;
00143 }
00144 default:{
00145 fprintf (stderr,"\n\n unknown type of transported matter is required");
00146 fprintf (stderr,"\n in function axiquadc::axiquadc (file %s, line %d).\n",__FILE__,__LINE__);
00147 }
00148 }
00149 }
00150
00151 if (Cp->lbb==quad_lin){
00152 switch (Tp->tmatt){
00153 case onemedium:{
00154 intordvum[0][0]=2; intordvlm[0][0]=2;
00155 nipu[0][0]=4; nipl[0][0]=4;
00156 dofe[0]=nnetp;
00157 break;
00158 }
00159 case twomediacoup:{
00160 intordvum[0][0]=2; intordvum[0][1]=2;
00161 intordvlm[0][0]=2; intordvlm[1][0]=2;
00162
00163 if (Cp->savemode==0){
00164 nipu[0][0]=4; nipu[0][1]=4;
00165 nipl[0][0]=4; nipl[1][0]=4;
00166 }
00167 if (Cp->savemode==1){
00168 nipu[0][0]=4; nipu[0][1]=0;
00169 nipl[0][0]=4; nipl[1][0]=0;
00170 }
00171
00172 dofe[0]=nnetp; dofe[1]=nnetp;
00173 break;
00174 }
00175 case threemediacoup:{
00176 intordvum[0][0]=2; intordvum[0][1]=2; intordvum[0][2]=2;
00177 intordvlm[0][0]=2; intordvlm[1][0]=2; intordvlm[2][0]=2;
00178
00179 if (Cp->savemode==0){
00180 nipu[0][0]=4; nipu[0][1]=4; nipu[0][2]=4;
00181 nipl[0][0]=4; nipl[1][0]=4; nipl[2][0]=4;
00182 }
00183 if (Cp->savemode==1){
00184 nipu[0][0]=4; nipu[0][1]=0; nipu[0][2]=0;
00185 nipl[0][0]=4; nipl[1][0]=0; nipl[2][0]=0;
00186 }
00187
00188 dofe[0]=nnetp; dofe[1]=nnetp; dofe[2]=nnetp;
00189 break;
00190 }
00191 default:{
00192 fprintf (stderr,"\n\n unknown type of transported matter is required");
00193 fprintf (stderr,"\n in function axiquadc::axiquadc (file %s, line %d).\n",__FILE__,__LINE__);
00194 }
00195 }
00196 }
00197
00198 if (Cp->lbb==quad_quad){
00199 switch (Tp->tmatt){
00200 case onemedium:{
00201 intordvum[0][0]=3; intordvlm[0][0]=3;
00202 nipu[0][0]=9; nipl[0][0]=9;
00203 dofe[0]=nnetp;
00204 break;
00205 }
00206 case twomediacoup:{
00207 intordvum[0][0]=3; intordvum[0][1]=3;
00208 intordvlm[0][0]=3; intordvlm[1][0]=3;
00209
00210 if (Cp->savemode==0){
00211 nipu[0][0]=9; nipu[0][1]=9;
00212 nipl[0][0]=9; nipl[1][0]=9;
00213 }
00214 if (Cp->savemode==1){
00215 nipu[0][0]=9; nipu[0][1]=0;
00216 nipl[0][0]=9; nipl[1][0]=0;
00217 }
00218
00219 dofe[0]=nnetp; dofe[1]=nnetp;
00220 break;
00221 }
00222 case threemediacoup:{
00223 intordvum[0][0]=3; intordvum[0][1]=3; intordvum[0][2]=3;
00224 intordvlm[0][0]=3; intordvlm[1][0]=3; intordvlm[2][0]=3;
00225
00226 if (Cp->savemode==0){
00227 nipu[0][0]=9; nipu[0][1]=9; nipu[0][2]=9;
00228 nipl[0][0]=9; nipl[1][0]=9; nipl[2][0]=9;
00229 }
00230 if (Cp->savemode==1){
00231 nipu[0][0]=9; nipu[0][1]=0; nipu[0][2]=0;
00232 nipl[0][0]=9; nipl[1][0]=0; nipl[2][0]=0;
00233 }
00234
00235 dofe[0]=nnetp; dofe[1]=nnetp; dofe[2]=nnetp;
00236 break;
00237 }
00238 default:{
00239 fprintf (stderr,"\n\n unknown type of transported matter is required");
00240 fprintf (stderr,"\n in function axiquadc::axiquadc (file %s, line %d).\n",__FILE__,__LINE__);
00241 }
00242 }
00243 }
00244
00245
00246 tnipu=0;
00247 for (i=0;i<mnb;i++){
00248 for (j=0;j<ntm;j++)
00249 tnipu+=nipu[i][j];
00250 }
00251 tnipl=0;
00252 for (i=0;i<mnb;i++){
00253 for (j=0;j<ntm;j++)
00254 tnipl+=nipl[i][j];
00255 }
00256
00257 }
00258
00259
00260
00261
00262 axiquadc::~axiquadc (void)
00263 {
00264 long i;
00265
00266 for (i=0;i<ntm;i++){
00267 delete [] nipl[i];
00268 delete [] intordvlm[i];
00269 }
00270 for (i=0;i<mnb;i++){
00271 delete [] nipu[i];
00272 delete [] intordvum[i];
00273 }
00274 delete [] dofe;
00275 delete [] nipu;
00276 delete [] nipl;
00277 delete [] intordvum;
00278 delete [] intordvlm;
00279 }
00280
00281
00282
00283
00284
00285
00286
00287
00288 void axiquadc::eleminit (long eid)
00289 {
00290 long ii,jj;
00291
00292 Ct->elements[eid].intordvum = new long* [mnb];
00293 Ct->elements[eid].intordvlm = new long* [ntm];
00294 Ct->elements[eid].nipu = new long* [mnb];
00295 Ct->elements[eid].nipl = new long* [ntm];
00296
00297 for (ii=0;ii<mnb;ii++){
00298 Ct->elements[eid].intordvum[ii] = new long [ntm];
00299 Ct->elements[eid].nipu[ii] = new long [ntm];
00300 for (jj=0;jj<ntm;jj++){
00301 Ct->elements[eid].intordvum[ii][jj]=intordvum[ii][jj];
00302 Ct->elements[eid].nipu[ii][jj]=nipu[ii][jj];
00303 }
00304 }
00305 for (ii=0;ii<ntm;ii++){
00306 Ct->elements[eid].intordvlm[ii] = new long [mnb];
00307 Ct->elements[eid].nipl[ii] = new long [mnb];
00308 for (jj=0;jj<mnb;jj++){
00309 Ct->elements[eid].intordvlm[ii][jj]=intordvlm[ii][jj];
00310 Ct->elements[eid].nipl[ii][jj]=nipl[ii][jj];
00311 }
00312 }
00313 }
00314
00315
00316
00317
00318
00319
00320
00321
00322
00323
00324 void axiquadc::dmatblockcol (long ri,long ci,matrix &d, matrix &dd)
00325 {
00326 fillm (0.0,dd);
00327
00328 if (mnb==1){
00329 dd[0][0]=d[0][0];
00330 dd[1][0]=d[1][0];
00331 dd[2][0]=d[2][0];
00332 dd[3][0]=d[3][0];
00333 }
00334 if (mnb==2){
00335 if (ri==0 && ci==0){
00336 dd[0][0]=d[0][0];
00337 dd[1][0]=d[1][0];
00338 }
00339 if (ri==1 && ci==0){
00340 dd[0][0]=d[2][0];
00341 }
00342 if (ri==2 && ci==0){
00343 dd[0][0]=d[3][0];
00344 }
00345 }
00346 }
00347
00348
00349
00350
00351
00352
00353
00354
00355
00356 void axiquadc::dmatblockrow (long ri,long ci,matrix &d, matrix &dd)
00357 {
00358 fillm (0.0,dd);
00359
00360 if (mnb==1){
00361 dd[0][0]=d[0][0];
00362 dd[0][1]=d[0][1];
00363 dd[0][2]=d[0][2];
00364 dd[0][3]=d[0][3];
00365 }
00366 if (mnb==2){
00367 if (ri==0 && ci==0){
00368 dd[0][0]=d[0][0];
00369 dd[0][1]=d[0][1];
00370 }
00371 if (ri==0 && ci==1){
00372 dd[0][0]=d[0][2];
00373 }
00374 if (ri==0 && ci==2){
00375 dd[0][0]=d[0][3];
00376 }
00377 }
00378 }
00379
00380
00381
00382
00383
00384
00385
00386
00387
00388
00389
00390 void axiquadc::upper_cond_coup_matrix (long eid,long ri,long ci,matrix &vm)
00391 {
00392 long i,j,ipp;
00393 double xi,eta,jac,r;
00394 vector x(nnemp),y(nnemp),w(intordvum[ri][ci]),gp(intordvum[ri][ci]);
00395 matrix gm(mncomp[ri],mndofe),d(tnmcomp,1),dd(mncomp[ri],1),n(1,dofe[ci]);
00396
00397
00398 Mt->give_node_coord2d (x,y,eid);
00399
00400
00401 gauss_points (gp.a,w.a,intordvum[ri][ci]);
00402
00403 if (Cp->savemode==0)
00404 ipp=Ct->elements[eid].ippu[ri][ci];
00405 if (Cp->savemode==1)
00406 ipp=Ct->elements[eid].ippu[0][0];
00407
00408 for (i=0;i<intordvum[ri][ci];i++){
00409 xi=gp[i];
00410 for (j=0;j<intordvum[ri][ci];j++){
00411 eta=gp[j];
00412
00413 if (Cp->lbb==lin_lin){
00414
00415 Asymlq->geom_matrix_block (gm,ri,x,y,xi,eta,jac);
00416
00417 Lqat->bf_matrix (n,xi,eta);
00418
00419 r = Asymlq->approx (xi,eta,x);
00420 }
00421 if (Cp->lbb==quad_lin){
00422
00423 Asymqq->geom_matrix_block (gm,ri,x,y,xi,eta,jac);
00424
00425 Lqat->bf_matrix (n,xi,eta);
00426
00427 r = Asymqq->approx (xi,eta,x);
00428 }
00429 if (Cp->lbb==quad_quad){
00430
00431 Asymqq->geom_matrix_block (gm,ri,x,y,xi,eta,jac);
00432
00433 Qqat->bf_matrix (n,xi,eta);
00434
00435 r = Asymqq->approx (xi,eta,x);
00436 }
00437
00438
00439
00440 Cmu->matcond (d,ipp,ri,ci);
00441
00442 jac*=w[i]*w[j]*r;
00443
00444 bdbjac (vm,gm,d,n,jac);
00445
00446 ipp++;
00447 }
00448 }
00449 }
00450
00451
00452
00453
00454
00455
00456
00457
00458
00459
00460
00461 void axiquadc::lower_cond_coup_matrix (long eid,long ri,long ci,matrix &vm)
00462 {
00463 long i,j,ipp;
00464 double xi,eta,jac,r;
00465 vector x(nnemp),y(nnemp),w(intordvum[ri][ci]),gp(intordvum[ri][ci]);
00466 matrix gm(tnmcomp,mndofe),d(1,tnmcomp),n(1,nnetp);
00467
00468
00469 Mt->give_node_coord2d (x,y,eid);
00470
00471
00472 gauss_points (gp.a,w.a,intordvlm[ri][ci]);
00473
00474 if (Cp->savemode==0)
00475 ipp=Ct->elements[eid].ippl[ri][ci];
00476 if (Cp->savemode==1)
00477 ipp=Ct->elements[eid].ippl[0][0];
00478
00479 for (i=0;i<intordvlm[ri][ci];i++){
00480 xi=gp[i];
00481 for (j=0;j<intordvlm[ri][ci];j++){
00482 eta=gp[j];
00483
00484 if (Cp->lbb==lin_lin){
00485
00486 Asymlq->geom_matrix (gm,x,y,xi,eta,jac);
00487
00488 Lqat->bf_matrix (n,xi,eta);
00489
00490 r = Asymlq->approx (xi,eta,x);
00491 }
00492 if (Cp->lbb==quad_lin){
00493
00494 Asymqq->geom_matrix (gm,x,y,xi,eta,jac);
00495
00496 Lqat->bf_matrix (n,xi,eta);
00497
00498 r = Asymqq->approx (xi,eta,x);
00499 }
00500 if (Cp->lbb==quad_quad){
00501
00502 Asymqq->geom_matrix (gm,x,y,xi,eta,jac);
00503
00504 Qqat->bf_matrix (n,xi,eta);
00505
00506 r = Asymqq->approx (xi,eta,x);
00507 }
00508
00509
00510
00511 Cml->matcond (d,ipp,ri,ci);
00512
00513 jac*=w[i]*w[j]*r;
00514
00515 bdbjac (vm,n,d,gm,jac);
00516
00517 ipp++;
00518 }
00519 }
00520 }
00521
00522
00523
00524
00525
00526
00527
00528
00529
00530
00531
00532 void axiquadc::upper_cap_coup_matrix (long eid,long ri,long ci,matrix &vm)
00533 {
00534 long i,j,ipp;
00535 double xi,eta,jac,r;
00536 vector x(nnemp),y(nnemp),w(intordvum[ri][ci]),gp(intordvum[ri][ci]);
00537 matrix gm(tnmcomp,mndofe),d(tnmcomp,1),n(1,nnetp);
00538
00539
00540 Mt->give_node_coord2d (x,y,eid);
00541
00542
00543 gauss_points (gp.a,w.a,intordvum[ri][ci]);
00544
00545 if (Cp->savemode==0)
00546 ipp=Ct->elements[eid].ippu[ri][ci];
00547 if (Cp->savemode==1)
00548 ipp=Ct->elements[eid].ippu[0][0];
00549
00550 for (i=0;i<intordvum[ri][ci];i++){
00551 xi=gp[i];
00552 for (j=0;j<intordvum[ri][ci];j++){
00553 eta=gp[j];
00554
00555 if (Cp->lbb==lin_lin){
00556
00557 Asymlq->geom_matrix (gm,x,y,xi,eta,jac);
00558
00559 Lqat->bf_matrix (n,xi,eta);
00560
00561 r = Asymlq->approx (xi,eta,x);
00562 }
00563 if (Cp->lbb==quad_lin){
00564
00565 Asymqq->geom_matrix (gm,x,y,xi,eta,jac);
00566
00567 Lqat->bf_matrix (n,xi,eta);
00568
00569 r = Asymqq->approx (xi,eta,x);
00570 }
00571 if (Cp->lbb==quad_quad){
00572
00573 Asymqq->geom_matrix (gm,x,y,xi,eta,jac);
00574
00575 Qqat->bf_matrix (n,xi,eta);
00576
00577 r = Asymqq->approx (xi,eta,x);
00578 }
00579
00580
00581
00582 Cmu->matcap (d,ipp,ri,ci);
00583
00584 jac*=w[i]*w[j]*r;
00585
00586 bdbjac (vm,gm,d,n,jac);
00587
00588 ipp++;
00589 }
00590 }
00591 }
00592
00593
00594
00595
00596
00597
00598
00599
00600
00601
00602 void axiquadc::lower_cap_coup_matrix (long eid,long ri,long ci,matrix &vm)
00603 {
00604 long i,j,ipp;
00605 double xi,eta,jac,r;
00606 vector x(nnemp),y(nnemp),w(intordvum[ri][ci]),gp(intordvum[ri][ci]);
00607 matrix gm(tnmcomp,mndofe),d(1,tnmcomp),n(1,nnetp);
00608
00609
00610 Mt->give_node_coord2d (x,y,eid);
00611
00612
00613 gauss_points (gp.a,w.a,intordvlm[ri][ci]);
00614
00615 if (Cp->savemode==0)
00616 ipp=Ct->elements[eid].ippl[ri][ci];
00617 if (Cp->savemode==1)
00618 ipp=Ct->elements[eid].ippl[0][0];
00619
00620 for (i=0;i<intordvlm[ri][ci];i++){
00621 xi=gp[i];
00622 for (j=0;j<intordvlm[ri][ci];j++){
00623 eta=gp[j];
00624
00625 if (Cp->lbb==lin_lin){
00626
00627 Asymlq->geom_matrix (gm,x,y,xi,eta,jac);
00628
00629 Lqat->bf_matrix (n,xi,eta);
00630
00631 r = Asymlq->approx (xi,eta,x);
00632 }
00633 if (Cp->lbb==quad_lin){
00634
00635 Asymqq->geom_matrix (gm,x,y,xi,eta,jac);
00636
00637 Lqat->bf_matrix (n,xi,eta);
00638
00639 r = Asymqq->approx (xi,eta,x);
00640 }
00641 if (Cp->lbb==quad_quad){
00642
00643 Asymqq->geom_matrix (gm,x,y,xi,eta,jac);
00644
00645 Qqat->bf_matrix (n,xi,eta);
00646
00647 r = Asymqq->approx (xi,eta,x);
00648 }
00649
00650
00651
00652 Cml->matcap (d,ipp,ri,ci);
00653
00654 jac*=w[i]*w[j]*r;
00655
00656 bdbjac (vm,n,d,gm,jac);
00657
00658 ipp++;
00659 }
00660 }
00661 }
00662
00663
00664
00665
00666
00667
00668
00669
00670
00671
00672
00673
00674
00675
00676
00677
00678
00679
00680
00681 void axiquadc::res_upper_cond_coup_matrix (long eid,matrix &vm)
00682 {
00683 long i,j,*ccn;
00684 matrix lvm(mndofe,nnetp);
00685
00686 ccn = new long [nnetp];
00687
00688 fillm(0.0,vm);
00689
00690 for (i=0;i<mnb;i++){
00691 for (j=0;j<ntm;j++){
00692
00693
00694 upper_cond_coup_matrix (eid,i,j,lvm);
00695
00696 if (Cp->lbb==lin_lin || Cp->lbb==quad_lin){
00697 Lqat->codnum (ccn,j);
00698 }
00699 if (Cp->lbb==quad_quad){
00700 Qqat->codnum (ccn,j);
00701 }
00702
00703 mat_localize (vm,lvm,mordering,ccn);
00704 }
00705 }
00706 delete [] ccn;
00707 }
00708
00709
00710
00711
00712
00713
00714
00715
00716
00717 void axiquadc::res_lower_cond_coup_matrix (long eid,matrix &vm)
00718 {
00719 long i,j,*rcn;
00720 matrix lvm(nnetp,mndofe);
00721
00722 rcn = new long [nnetp];
00723
00724 fillm(0.0,vm);
00725
00726 for (i=0;i<ntm;i++){
00727 for (j=0;j<mnb;j++){
00728
00729
00730 lower_cond_coup_matrix (eid,i,j,lvm);
00731
00732 if (Cp->lbb==lin_lin || Cp->lbb==quad_lin){
00733 Lqat->codnum (rcn,j);
00734 }
00735 if (Cp->lbb==quad_quad){
00736 Qqat->codnum (rcn,j);
00737 }
00738
00739 mat_localize (vm,lvm,rcn,mordering);
00740 }
00741 }
00742 delete [] rcn;
00743 }
00744
00745
00746
00747
00748
00749
00750
00751
00752
00753 void axiquadc::res_upper_cap_coup_matrix (long eid,matrix &vm)
00754 {
00755 long i,j,*ccn;
00756 matrix lvm(mndofe,nnetp);
00757
00758 ccn = new long [nnetp];
00759
00760 fillm(0.0,vm);
00761
00762 for (i=0;i<mnb;i++){
00763 for (j=0;j<ntm;j++){
00764
00765
00766 upper_cap_coup_matrix (eid,i,j,lvm);
00767
00768 if (Cp->lbb==lin_lin || Cp->lbb==quad_lin){
00769 Lqat->codnum (ccn,j);
00770 }
00771 if (Cp->lbb==quad_quad){
00772 Qqat->codnum (ccn,j);
00773 }
00774
00775 mat_localize (vm,lvm,mordering,ccn);
00776 }
00777 }
00778 delete [] ccn;
00779 }
00780
00781
00782
00783
00784
00785
00786
00787
00788
00789 void axiquadc::res_lower_cap_coup_matrix (long eid,matrix &vm)
00790 {
00791 long i,j,*rcn;
00792 matrix lvm(nnetp,mndofe);
00793
00794 rcn = new long [nnetp];
00795
00796 fillm(0.0,vm);
00797
00798 for (i=0;i<ntm;i++){
00799 for (j=0;j<mnb;j++){
00800
00801
00802 lower_cap_coup_matrix (eid,i,j,lvm);
00803
00804 if (Cp->lbb==lin_lin || Cp->lbb==quad_lin){
00805 Lqat->codnum (rcn,j);
00806 }
00807 if (Cp->lbb==quad_quad){
00808 Qqat->codnum (rcn,j);
00809 }
00810
00811 mat_localize (vm,lvm,rcn,mordering);
00812 }
00813 }
00814 delete [] rcn;
00815 }
00816
00817
00818
00819
00820
00821
00822
00823
00824
00825
00826
00827
00828
00829
00830
00831 void axiquadc::upper_cond_coup_vector (vector &tvm,vector &nodval,long eid,long ri,long ci)
00832 {
00833 long i,j,ipp;
00834 double xi,eta,jac,r;
00835 ivector nodes(nnetp);
00836 vector x(nnemp),y(nnemp),w(intordvum[ri][ci]),gp(intordvum[ri][ci]),av(mndofe),t(nnetp);
00837 matrix gm(mncomp[ri],mndofe),d(tnmcomp,1),dd(mncomp[ri],1),n(1,dofe[ci]),ucm(mndofe,dofe[ci]);
00838
00839
00840 Tt->give_elemnodes (eid,nodes);
00841
00842 Tc->give_thickness (eid,nodes,t);
00843
00844 Mt->give_node_coord2d (x,y,eid);
00845
00846 gauss_points (gp.a,w.a,intordvum[ri][ci]);
00847
00848 fillm (0.0,ucm);
00849
00850 if (Cp->savemode==0){
00851 ipp=Ct->elements[eid].ippu[ri][ci];
00852 }
00853 if (Cp->savemode==1){
00854 ipp=Ct->elements[eid].ippu[0][0];
00855 }
00856
00857 for (i=0;i<intordvum[ri][ci];i++){
00858 xi=gp[i];
00859 for (j=0;j<intordvum[ri][ci];j++){
00860 eta=gp[j];
00861
00862 if (Cp->lbb==lin_lin){
00863 Asymlq->geom_matrix_block (gm,ri,x,y,xi,eta,jac);
00864
00865 r = Asymlq->approx (xi,eta,x);
00866 Lqat->bf_matrix (n,xi,eta);
00867 }
00868 if (Cp->lbb==quad_lin){
00869 Asymqq->geom_matrix_block (gm,ri,x,y,xi,eta,jac);
00870
00871 r = Asymqq->approx (xi,eta,x);
00872 Lqat->bf_matrix (n,xi,eta);
00873 }
00874 if (Cp->lbb==quad_quad){
00875 Asymqq->geom_matrix_block (gm,ri,x,y,xi,eta,jac);
00876
00877 r = Asymqq->approx (xi,eta,x);
00878 Qqat->bf_matrix (n,xi,eta);
00879 }
00880
00881
00882
00883 Cmu->volume_rhs1 (d,ipp,ri,ci);
00884
00885
00886 dmatblockcol (ri,ci,d,dd);
00887
00888 jac*=r*w[i]*w[j];
00889
00890
00891 bdbjac (ucm,gm,dd,n,jac);
00892
00893 ipp++;
00894 }
00895 }
00896
00897 mxv (ucm,nodval,av);
00898 addv (av,tvm,tvm);
00899
00900 }
00901
00902
00903
00904
00905
00906
00907
00908
00909
00910
00911
00912
00913 void axiquadc::res_upper_cond_coup_vector (vector &f,long eid)
00914 {
00915 long i,j;
00916 long *cn,*lcn;
00917 double *r;
00918 vector lr;
00919
00920 cn = new long [tndofe];
00921 r = new double [tndofe];
00922
00923
00924 Tt->give_code_numbers (eid,cn);
00925
00926
00927 initialvalues (r,cn,tndofe);
00928
00929 fillv (0.0,f);
00930
00931 for (i=0;i<ntm;i++){
00932
00933 lcn = new long [dofe[i]];
00934
00935 allocv (dofe[i],lr);
00936
00937
00938 if (Cp->lbb==lin_lin || Cp->lbb==quad_lin)
00939 Lqat->codnum (lcn,i);
00940 else
00941 Qqat->codnum (lcn,i);
00942
00943
00944 globloc (r,lr.a,lcn,dofe[i]);
00945
00946 for (j=0;j<mnb;j++){
00947
00948 upper_cond_coup_vector (f,lr,eid,j,i);
00949 }
00950 delete [] lcn;
00951 destrv (lr);
00952 }
00953
00954 delete [] cn;
00955 delete [] r;
00956 }
00957
00958
00959
00960
00961
00962
00963
00964
00965
00966
00967
00968
00969
00970
00971
00972
00973
00974
00975
00976
00977
00978
00979
00980
00981
00982
00983
00984
00985
00986
00987
00988
00989
00990
00991
00992
00993
00994
00995
00996
00997
00998
00999
01000
01001
01002
01003
01004
01005
01006
01007
01008
01009
01010
01011
01012
01013
01014
01015
01016
01017
01018
01019
01020
01021
01022
01023
01024
01025
01026
01027
01028
01029
01030
01031
01032
01033
01034
01035
01036
01037
01038
01039
01040
01041
01042
01043
01044
01045
01046
01047
01048
01049
01050
01051
01052
01053
01054
01055
01056
01057
01058
01059
01060
01061
01062
01063
01064
01065
01066
01067
01068
01069
01070
01071
01072
01073
01074
01075
01076
01077
01078
01079
01080
01081
01082
01083
01084
01085
01086
01087
01088
01089
01090
01091
01092
01093
01094
01095
01096
01097
01098
01099
01100
01101
01102
01103
01104
01105
01106
01107
01108
01109
01110
01111
01112
01113
01114
01115
01116
01117
01118
01119
01120
01121
01122
01123
01124
01125
01126
01127
01128
01129
01130
01131
01132
01133
01134
01135
01136
01137
01138
01139
01140
01141
01142
01143
01144
01145
01146
01147
01148
01149
01150
01151
01152
01153
01154
01155
01156
01157
01158
01159
01160
01161
01162
01163
01164
01165
01166
01167
01168
01169
01170
01171
01172