shp3dv.cppGo to the documentation of this file.00001 // Calculate Shape functions and parametric derivatives for a 00002 00003 // 3-D finite element with varied nen between 8 and 27. 00004 00005 00006 00007 // By Jinchi Lu, 04/15/2004 00008 00009 // (Based on a Fortran code from R. Brockman 30 May 1986) 00010 00011 00012 00013 #include <stdio.h> 00014 00015 #include <math.h> 00016 00017 00018 00019 void shap3dv(double *R, int *NP, double Q[27][4]){ 00020 00021 00022 00023 // Formal parameters: 00024 00025 // 00026 00027 // R (Input) = Parametric coordinates, in the interval [-1,1], 00028 00029 // at which shape functions are to be evaluated 00030 00031 // NP (Input) = List of flags for each local node of the finite 00032 00033 // element. (Usually this is just the connectivity list) 00034 00035 // = 0 if node is absent 00036 00037 // > 0 if node is present 00038 00039 // * Q (Output) = Shape functions and derivatives, as follows: 00040 00041 // * Q[i][3] is the shape function for local node 'i' 00042 00043 // * Q[i][0] is the r-derivative 00044 00045 // * Q[i][1] is the s-derivative 00046 00047 // * Q[i][2] is the t-derivative 00048 00049 // 00050 00051 // Local Node Pattern for Three-Dimensional Elements: 00052 00053 // 00054 00055 // * Nodes 1 - 4 Lower surface, counterclockwise 00056 00057 // * Nodes 5 - 8 Upper surface, counterclockwise 00058 00059 // * Nodes 9 - 12 Midsides of edges 1-2, 2-3, 3-4, 4-1 00060 00061 // * Nodes 13 - 16 Midsides of edges 5-6, 6-7, 7-8, 8-5 00062 00063 // * Nodes 17 - 20 Midsides of edges 1-5, 2-6, 3-7, 4-8 00064 00065 // * Nodes 21 - 26 Mid-face nodes on +r, +s, +t, -r, -s, -t 00066 00067 // * Node 27 Centroid node 00068 00069 // 00070 00071 00072 00073 double G[3][3], D[3][3], C; 00074 00075 int L[27] = {3,1,1,3,3,1,1,3,2,1,2,3,2,1,2,3,3,1,1,3,1,2,2,3,2,2,2}, 00076 00077 M[27] = {3,3,1,1,3,3,1,1,3,2,1,2,3,2,1,2,3,3,1,1,2,1,2,2,3,2,2}, 00078 00079 N[27] = {3,3,3,3,1,1,1,1,3,3,3,3,1,1,1,1,2,2,2,2,2,2,1,2,2,3,2}; 00080 00081 00082 00083 int i, j, LR, LS, LT; 00084 00085 for( i = 0; i < 3; i++ ) { 00086 00087 G[0][i] = 0.5 + 0.5*R[i]; 00088 00089 G[1][i] = 1.0 - R[i]*R[i]; 00090 00091 G[2][i] = 0.5 - 0.5*R[i]; 00092 00093 D[0][i] = 0.5 ; 00094 00095 D[1][i] = -2.0*R[i]; 00096 00097 D[2][i] = -0.5; 00098 00099 } 00100 00101 00102 00103 // Construct basic three-dimensional quadratic shape functions 00104 00105 00106 00107 for( i = 0; i < 27; i++ ) { 00108 00109 LR = L[i]-1; 00110 00111 LS = M[i]-1; 00112 00113 LT = N[i]-1; 00114 00115 Q[i][0] = D[LR][0] * G[LS][1] * G[LT][2]; 00116 00117 Q[i][1] = G[LR][0] * D[LS][1] * G[LT][2]; 00118 00119 Q[i][2] = G[LR][0] * G[LS][1] * D[LT][2]; 00120 00121 Q[i][3] = G[LR][0] * G[LS][1] * G[LT][2]; 00122 00123 // printf("%d %15.6e %15.6e %15.6e %15.6e\n", i, Q[i][0], Q[i][1], Q[i][2], Q[i][3]); 00124 00125 } 00126 00127 00128 00129 // Modify basic shape functions to account for omitted nodes 00130 00131 00132 00133 for(j = 0; j < 4; j++ ) { 00134 00135 if( NP[26] == 0 ) Q[26][j] = 0.; 00136 00137 C = -0.5*Q[26][j]; 00138 00139 for( i = 20; i < 26; i++ ) { 00140 00141 Q[i][j] += C ; 00142 00143 if( NP[i] == 0 ) Q[i][j] = 0.; 00144 00145 } 00146 00147 00148 00149 C = 0.5*C; 00150 00151 Q[ 8][j] += - 0.5* (Q[25][j] + Q[24][j]) + C; 00152 00153 Q[ 9][j] += - 0.5* (Q[25][j] + Q[20][j]) + C; 00154 00155 Q[10][j] += - 0.5* (Q[25][j] + Q[21][j]) + C; 00156 00157 Q[11][j] += - 0.5* (Q[25][j] + Q[23][j]) + C; 00158 00159 Q[12][j] += - 0.5* (Q[24][j] + Q[22][j]) + C; 00160 00161 Q[13][j] += - 0.5* (Q[20][j] + Q[22][j]) + C; 00162 00163 Q[14][j] += - 0.5* (Q[21][j] + Q[22][j]) + C; 00164 00165 Q[15][j] += - 0.5* (Q[23][j] + Q[22][j]) + C; 00166 00167 Q[16][j] += - 0.5* (Q[23][j] + Q[24][j]) + C; 00168 00169 Q[17][j] += - 0.5* (Q[24][j] + Q[20][j]) + C; 00170 00171 Q[18][j] += - 0.5* (Q[20][j] + Q[21][j]) + C; 00172 00173 Q[19][j] += - 0.5* (Q[21][j] + Q[23][j]) + C; 00174 00175 00176 00177 for( i = 8; i < 20; i++ ) { 00178 00179 if( NP[i] == 0 ) Q[i][j] = 0.; 00180 00181 } 00182 00183 00184 00185 C = 0.5*C; 00186 00187 for( i = 0; i < 8; i++ ) { 00188 00189 Q[i][j] += C; 00190 00191 } 00192 00193 00194 00195 Q[0][j] += - 0.50* (Q[16][j] + Q[11][j] + Q[ 8][j]) - 0.25* (Q[23][j] + Q[24][j] + Q[25][j]); 00196 00197 Q[1][j] += - 0.50* (Q[17][j] + Q[ 9][j] + Q[ 8][j]) - 0.25* (Q[20][j] + Q[24][j] + Q[25][j]); 00198 00199 Q[2][j] += - 0.50* (Q[10][j] + Q[ 9][j] + Q[18][j]) - 0.25* (Q[25][j] + Q[20][j] + Q[21][j]); 00200 00201 Q[3][j] += - 0.50* (Q[19][j] + Q[11][j] + Q[10][j]) - 0.25* (Q[25][j] + Q[21][j] + Q[23][j]); 00202 00203 Q[4][j] += - 0.50* (Q[16][j] + Q[15][j] + Q[12][j]) - 0.25* (Q[22][j] + Q[23][j] + Q[24][j]); 00204 00205 Q[5][j] += - 0.50* (Q[17][j] + Q[12][j] + Q[13][j]) - 0.25* (Q[24][j] + Q[20][j] + Q[22][j]); 00206 00207 Q[6][j] += - 0.50* (Q[18][j] + Q[13][j] + Q[14][j]) - 0.25* (Q[20][j] + Q[21][j] + Q[22][j]); 00208 00209 Q[7][j] += - 0.50* (Q[19][j] + Q[15][j] + Q[14][j]) - 0.25* (Q[21][j] + Q[22][j] + Q[23][j]); 00210 00211 } 00212 00213 00214 00215 } 00216 00217 00218 00219 int brcshl(double shl[4][20][27], double w[27], int nint, int nen) { 00220 00221 /* 00222 00223 PROGRAM TO CALCULATE INTEGRATION-RULE WEIGHTS, SHAPE FUNCTIONS 00224 00225 AND LOCAL DERIVATIVES FOR A EIGHT-NODE BRICK ELEMENT 00226 00227 00228 00229 R,S,T = LOCAL ELEMENT COORD ("XI", "ETA", "ZETA" RESP.) 00230 00231 SHL(1,I,L) = LOCAL ("XI") DERIVATIVE OF SHAPE FUNCTION 00232 00233 SHL(2,I,L) = LOCAL ("ETA") DERIVATIVE OF SHAPE FUNCTION 00234 00235 SHL(3,I,L) = LOCAL ("ZETA") DERIVATIVE OF SHAPE FUNCTION 00236 00237 SHL(4,I,L) = LOCAL SHAPE FUNCTION 00238 00239 W(L) = INTEGRATION-RULE WEIGHT 00240 00241 I = LOCAL NODE NUMBER 00242 00243 L = INTEGRATION POINT NUMBER 00244 00245 NINT = NUMBER OF INTEGRATION POINTS, EQ. 1, 8 OR 27 00246 00247 */ 00248 00249 00250 00251 double RA[27] = {-0.50, 0.50, 0.50,-0.50,-0.50, 0.50, 0.50,-0.50, 00252 00253 0.00, 0.50, 0.00,-0.50, 0.00, 0.50, 0.00,-0.50, 00254 00255 -0.50, 0.50, 0.50,-0.50, 0.50, 0.00, 0.00,-0.50, 00256 00257 0.00, 0.00, 0.00}; 00258 00259 00260 00261 double SA[27] = {-0.50,-0.50, 0.50, 0.50,-0.50,-0.50, 0.50, 0.50, 00262 00263 -0.50, 0.00, 0.50, 0.00,-0.50, 0.00, 0.50, 0.00, 00264 00265 -0.50,-0.50, 0.50, 0.50, 0.00, 0.50, 0.00, 0.00, 00266 00267 -0.50, 0.00, 0.00}; 00268 00269 00270 00271 double TA[27] = {-0.50,-0.50,-0.50,-0.50, 0.50, 0.50, 0.50, 0.50, 00272 00273 -0.50,-0.50,-0.50,-0.50, 0.50, 0.50, 0.50, 0.50, 00274 00275 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.50, 0.00, 00276 00277 0.00,-0.50, 0.00}; 00278 00279 double five9 = 0.5555555555555556, eight9 = 0.8888888888888889; 00280 00281 double G = 0., R[3], Q[27][4]; 00282 00283 int i, j, L, NP[27]; 00284 00285 w[0] = 8; 00286 00287 if( nint == 8) { 00288 00289 G = 2/sqrt(3.0); 00290 00291 for(i = 0; i < nint; i++ ) 00292 00293 w[i] = 1.; 00294 00295 } 00296 00297 else if ( nint == 27 ) { 00298 00299 G = 2.0*sqrt(3./5.); 00300 00301 w[0] = five9 * five9 * five9; 00302 00303 for( i = 1; i < 8; i++ ) 00304 00305 w[i] = w[0]; 00306 00307 w[8] = five9 * five9 * eight9; 00308 00309 for( i = 9; i < 20; i++ ) 00310 00311 w[i] = w[8]; 00312 00313 w[20] = five9 * eight9 * eight9; 00314 00315 for( i = 21; i < 26; i++ ) 00316 00317 w[i] = w[20]; 00318 00319 w[26] = eight9 * eight9 * eight9; 00320 00321 } 00322 00323 else { 00324 00325 //printf("invalid nint %d\n", nint); 00326 00327 //exit(-1); 00328 00329 return -1; 00330 00331 } 00332 00333 00334 00335 for( i = 0; i < 27; i ++ ) 00336 00337 NP[i] = 1; 00338 00339 if( nen < 27) { 00340 00341 for( i = nen; i < 27; i++ ) 00342 00343 NP[i] = 0; 00344 00345 } 00346 00347 else if( nen < 8 ) { 00348 00349 //printf("invalid nen %d\n", nen); 00350 00351 //exit(-1); 00352 00353 return -1; 00354 00355 } 00356 00357 00358 00359 for( L = 0; L < nint; L++ ) { 00360 00361 R[0] = G * RA[L]; 00362 00363 R[1] = G * SA[L]; 00364 00365 R[2] = G * TA[L]; 00366 00367 00368 00369 shap3dv(R, NP, Q); 00370 00371 00372 00373 for( j = 0; j < nen; j++) { 00374 00375 for( i = 0; i < 4; i++ ) { 00376 00377 shl[i][j][L] = Q[j][i]; 00378 00379 } 00380 00381 //printf("%5d %5d %15.6e %15.6e %15.6e %15.6e\n", L+1, j+1, 00382 00383 // shl[0][j][L],shl[1][j][L],shl[2][j][L],shl[3][j][L]); 00384 00385 } 00386 00387 00388 00389 } 00390 00391 return 0; 00392 00393 } 00394 00395 00396 00397 00398 00399 00400 00401 00402 00403 00404 00405 00406 00407 00408 00409 00410 00411 00412 00413 00414 00415 00416 00417 00418 00419 00420 00421 00422 00423 00424 00425 00426 00427 00428 00429 00430 00431 00432 00433 00434 00435 00436
Generated on Mon Oct 23 15:05:11 2006 for OpenSees by |
opensees-support @ berkeley.edu
Supported by the National Science Foundation
©2006, UC Regents