00001 #include <stdlib.h>
00002 #include <math.h>
00003 #include "linwedge.h"
00004 #include "gadaptivity.h"
00005 #include "global.h"
00006 #include "globmat.h"
00007 #include "genfile.h"
00008 #include "intpoints.h"
00009 #include "node.h"
00010 #include "element.h"
00011 #include "loadcase.h"
00012
00013
00014
00015
00016 linwedge::linwedge (void)
00017 {
00018 long i;
00019
00020
00021 nne=6;
00022
00023 ndofe=18;
00024
00025 tncomp=6;
00026
00027 napfun=3;
00028
00029 intordmm=1;
00030
00031 ned=9;
00032
00033 nned=2;
00034
00035 nsurf=5;
00036
00037 nnsurf=4;
00038
00039 ssst=spacestress;
00040
00041
00042 nb=1;
00043
00044
00045 ncomp = new long [nb];
00046 ncomp[0]=6;
00047
00048
00049 cncomp = new long [nb];
00050 cncomp[0]=0;
00051
00052
00053
00054 nip = new long* [nb];
00055 intordsmt = new long* [nb];
00056 intordsmz = new long* [nb];
00057 for (i=0;i<nb;i++){
00058 nip[i] = new long [nb];
00059 intordsmt[i] = new long [nb];
00060 intordsmz[i] = new long [nb];
00061 }
00062
00063 nip[0][0]=6;
00064
00065
00066 tnip=6;
00067
00068
00069 intordsmt[0][0]=3;
00070 intordsmz[0][0]=2;
00071
00072 }
00073
00074 linwedge::~linwedge (void)
00075 {
00076 long i;
00077
00078 for (i=0;i<nb;i++){
00079 delete [] nip[i];
00080 delete [] intordsmt[i];
00081 delete [] intordsmz[i];
00082 }
00083 delete [] nip;
00084 delete [] intordsmt;
00085 delete [] intordsmz;
00086
00087 delete [] cncomp;
00088 delete [] ncomp;
00089 }
00090
00091
00092
00093
00094
00095
00096
00097
00098
00099
00100
00101 double linwedge::approx (double xi,double eta,double zeta,vector &nodval)
00102 {
00103 double f;
00104 vector bf(nne);
00105
00106 bf_lin_wed_3d (bf.a,xi,eta,zeta);
00107
00108 scprd (bf,nodval,f);
00109
00110 return f;
00111 }
00112
00113
00114
00115
00116
00117
00118
00119
00120
00121 void linwedge::bf_matrix (matrix &n,double xi,double eta,double zeta)
00122 {
00123 long i,j,k,l;
00124 vector bf(nne);
00125
00126 fillm (0.0,n);
00127
00128 bf_lin_wed_3d (bf.a,xi,eta,zeta);
00129
00130 j=0; k=1; l=2;
00131 for (i=0;i<nne;i++){
00132 n[0][j]=bf[i]; j+=3;
00133 n[1][k]=bf[i]; k+=3;
00134 n[2][l]=bf[i]; l+=3;
00135 }
00136
00137 }
00138
00139
00140
00141
00142
00143
00144
00145
00146
00147
00148
00149
00150
00151 void linwedge::geom_matrix (matrix &gm,vector &x,vector &y,vector &z,
00152 double xi,double eta,double zeta,double &jac)
00153 {
00154 long i,j,k,l;
00155 vector dx(nne),dy(nne),dz(nne);
00156
00157 dx_bf_lin_wed_3d (dx.a,zeta);
00158 dy_bf_lin_wed_3d (dy.a,zeta);
00159 dz_bf_lin_wed_3d (dz.a,xi,eta);
00160
00161 derivatives_3d (dx,dy,dz,jac,x,y,z,xi,eta,zeta);
00162
00163 fillm (0.0,gm);
00164
00165 j=0; k=1; l=2;
00166 for (i=0;i<nne;i++){
00167 gm[0][j]=dx[i];
00168 gm[1][k]=dy[i];
00169 gm[2][l]=dz[i];
00170
00171 gm[3][k]=dz[i];
00172 gm[3][l]=dy[i];
00173
00174 gm[4][j]=dz[i];
00175 gm[4][l]=dx[i];
00176
00177 gm[5][j]=dy[i];
00178 gm[5][k]=dx[i];
00179
00180 j+=3; k+=3; l+=3;
00181 }
00182
00183 }
00184
00185
00186
00187
00188
00189
00190
00191
00192
00193
00194 void linwedge::transf_matrix (ivector &nodes,matrix &tmat)
00195 {
00196 long i,n,m;
00197
00198 fillm (0.0,tmat);
00199
00200 n=nodes.n;
00201 m=tmat.m;
00202 for (i=0;i<m;i++){
00203 tmat[i][i]=1.0;
00204 }
00205
00206 for (i=0;i<n;i++){
00207 if (Mt->nodes[nodes[i]].transf>0){
00208 tmat[i*3+0][i*3]=Mt->nodes[nodes[i]].e1[0];
00209 tmat[i*3+1][i*3]=Mt->nodes[nodes[i]].e1[1];
00210 tmat[i*3+2][i*3]=Mt->nodes[nodes[i]].e1[2];
00211
00212 tmat[i*3+0][i*3+1]=Mt->nodes[nodes[i]].e2[0];
00213 tmat[i*3+1][i*3+1]=Mt->nodes[nodes[i]].e2[1];
00214 tmat[i*3+2][i*3+1]=Mt->nodes[nodes[i]].e2[2];
00215
00216 tmat[i*3+0][i*3+2]=Mt->nodes[nodes[i]].e3[0];
00217 tmat[i*3+1][i*3+2]=Mt->nodes[nodes[i]].e3[1];
00218 tmat[i*3+2][i*3+2]=Mt->nodes[nodes[i]].e3[2];
00219 }
00220 }
00221 }
00222
00223
00224
00225
00226
00227
00228
00229
00230
00231
00232
00233 void linwedge::stiffness_matrix (long eid,long ri,long ci,matrix &sm)
00234 {
00235 long i,k,ii,jj,ipp,transf;
00236 double xi,eta,zeta,jac;
00237 vector x(nne),y(nne),z(nne),w,wt,gp,gp1,gp2;
00238 matrix gm,d(tncomp,tncomp);
00239
00240 Mt->give_node_coord3d (x,y,z,eid);
00241
00242 fillm (0.0,sm);
00243
00244 for (ii=0;ii<nb;ii++){
00245 allocm (ncomp[ii],ndofe,gm);
00246 for (jj=0;jj<nb;jj++){
00247 if (intordsmt[ii][jj]==0) continue;
00248
00249 allocv (intordsmz[ii][jj],w);
00250 allocv (intordsmz[ii][jj],gp);
00251 allocv (intordsmt[ii][jj],wt);
00252 allocv (intordsmt[ii][jj],gp1);
00253 allocv (intordsmt[ii][jj],gp2);
00254
00255 gauss_points (gp.a,w.a,intordsmz[ii][jj]);
00256 gauss_points_tr (gp1.a,gp2.a,wt.a,intordsmt[ii][jj]);
00257
00258 ipp=Mt->elements[eid].ipp[ri+ii][ci+jj];
00259
00260 for (i=0;i<intordsmt[ii][jj];i++){
00261 xi=gp1[i]; eta=gp2[i];
00262 for (k=0;k<intordsmz[ii][jj];k++){
00263 zeta=gp[k];
00264
00265
00266 geom_matrix (gm,x,y,z,xi,eta,zeta,jac);
00267
00268 Mm->matstiff (d,ipp); ipp++;
00269
00270 jac*=wt[i]*w[k];
00271
00272
00273 bdbjac (sm,gm,d,gm,jac);
00274
00275 }
00276 }
00277 destrv (gp); destrv (w);
00278 }
00279 destrm (gm);
00280 }
00281
00282
00283 ivector nodes (nne);
00284 Mt->give_elemnodes (eid,nodes);
00285 transf = Mt->locsystems (nodes);
00286 if (transf>0){
00287 matrix tmat (ndofe,ndofe);
00288 transf_matrix (nodes,tmat);
00289 glmatrixtransf (sm,tmat);
00290 }
00291 }
00292
00293
00294
00295
00296
00297
00298
00299
00300
00301 void linwedge::res_stiffness_matrix (long eid,matrix &sm)
00302 {
00303 stiffness_matrix (eid,0,0,sm);
00304 }
00305