FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
ElasticLinear.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
7 use hecmw_util
8 use mmaterial
9
10 implicit none
11
12contains
13
15 subroutine calelasticmatrix( matl, sectType, D, temp )
16 type( tmaterial ), intent(in) :: matl
17 integer, intent(in) :: sectType
18 real(kind=kreal), intent(out) :: d(:,:)
19 real(kind=kreal), optional :: temp
20 real(kind=kreal) :: ee, pp, coef1, coef2, ina(1), outa(2)
21 logical :: ierr
22
23 d(:,:)=0.d0
24 if( present(temp) ) then
25 ina(1) = temp
26 call fetch_tabledata( mc_isoelastic, matl%dict, outa, ierr, ina )
27 if( ierr ) then
28 ee = matl%variables(m_youngs)
29 pp = matl%variables(m_poisson)
30 else
31 ee = outa(1)
32 pp = outa(2)
33 endif
34 else
35 call fetch_tabledata( mc_isoelastic, matl%dict, outa, ierr )
36 if( ierr ) then
37 ee = matl%variables(m_youngs)
38 pp = matl%variables(m_poisson)
39 else
40 ee = outa(1)
41 pp = outa(2)
42 endif
43 endif
44
45
46 select case (secttype)
47 case (d3)
48 d(1,1)=ee*(1.d0-pp)/(1.d0-2.d0*pp)/(1.d0+pp)
49 d(1,2)=ee*pp/(1.d0-2.d0*pp)/(1.d0+pp)
50 d(1,3)=d(1,2)
51 d(2,1)=d(1,2)
52 d(2,2)=d(1,1)
53 d(2,3)=d(1,2)
54 d(3,1)=d(1,3)
55 d(3,2)=d(2,3)
56 d(3,3)=d(1,1)
57 d(4,4)=ee/(1.d0+pp)*0.5d0
58 d(5,5)=ee/(1.d0+pp)*0.5d0
59 d(6,6)=ee/(1.d0+pp)*0.5d0
60 case (planestress)
61 coef1=ee/(1.d0-pp*pp)
62 coef2=0.5d0*(1.d0-pp)
63 d(1,1)=coef1
64 d(1,2)=coef1*pp
65 d(1,3)=0.d0
66 d(2,1)=d(1,2)
67 d(2,2)=d(1,1)
68 d(2,3)=0.d0
69 d(3,1)=0.d0
70 d(3,2)=0.d0
71 d(3,3)=coef1*coef2
72 case (planestrain)
73 coef1=ee/((1.d0+pp)*(1.d0-2.d0*pp))
74 coef2=ee/(2.d0*(1.d0+pp))
75 d(1,1)=coef1*(1.d0-pp)
76 d(1,2)=coef1*pp
77 d(1,3)=0.d0
78 d(2,1)=d(1,2)
79 d(2,2)=d(1,1)
80 d(2,3)=0.d0
81 d(3,1)=0.d0
82 d(3,2)=0.d0
83 d(3,3)=coef2
84 case (axissymetric)
85 coef1=ee*(1.d0-pp)/((1.d0+pp)*(1.d0-2.d0*pp))
86 coef2=(1.d0-2.d0*pp)/(2.d0*(1.d0-pp))
87 d(1,1)=coef1
88 d(1,2)=coef1*pp/(1.d0-pp)
89 d(1,3)=0.d0
90 d(1,4)=d(1,2)
91 d(2,1)=d(1,2)
92 d(2,2)=d(1,1)
93 d(2,3)=0.d0
94 d(2,4)=d(1,2)
95 d(3,1)=0.d0
96 d(3,2)=0.d0
97 d(3,3)=coef1*coef2
98 d(3,4)=0.d0
99 d(4,1)=d(1,4)
100 d(4,2)=d(2,4)
101 d(4,3)=0.d0
102 d(4,4)=d(1,1)
103 case default
104 stop "Section type not defined"
105 end select
106
107 end subroutine
108
109
111 subroutine calelasticmatrix_ortho( matl, sectType, bij, DMAT, temp )
112 use m_utilities
113 type( tmaterial ), intent(in) :: matl
114 integer, intent(in) :: sectType
115 real(kind=kreal), intent(in) :: bij(3,3)
116 real(kind=kreal), intent(out) :: dmat(:,:)
117 real(kind=kreal), optional :: temp
118 real(kind=kreal) :: e1, e2, e3, g12, g23, g13, nyu12, nyu23,nyu13
119 real(kind=kreal) :: nyu21,nyu32,nyu31, delta1, ina(1), outa(9)
120 real(kind=kreal) :: tm(6,6)
121 logical :: ierr
122
123 if( present(temp) ) then
124 ina(1) = temp
125 call fetch_tabledata( mc_orthoelastic, matl%dict, outa, ierr, ina )
126 if( ierr ) then
127 stop "Fails in fetching orthotropic elastic constants!"
128 endif
129 else
130 call fetch_tabledata( mc_orthoelastic, matl%dict, outa, ierr )
131 if( ierr ) then
132 stop "Fails in fetching orthotropic elastic constants!"
133 endif
134 endif
135
136 e1 = outa(1)
137 e2 = outa(2)
138 e3 = outa(3)
139 nyu12 = outa(4)
140 nyu13 = outa(5)
141 nyu23 = outa(6)
142 g12 = outa(7)
143 g13 = outa(8)
144 g23 = outa(9)
145 nyu21 = e2/e1*nyu12
146 nyu32 = e3/e2*nyu23
147 nyu31 = e3/e1*nyu13
148 delta1 = 1.d0/(1.d0 -nyu12*nyu21 -nyu23*nyu32 -nyu31*nyu13 -2.d0*nyu21*nyu32*nyu13)
149
150 dmat(:,:)=0.d0
151 dmat(1,1) = e1*(1.d0-nyu23*nyu32)*delta1
152 dmat(2,2) = e2*(1.d0-nyu13*nyu31)*delta1
153 dmat(3,3) = e3*(1.d0-nyu12*nyu21)*delta1
154 dmat(1,2) = e1*(nyu21+nyu31*nyu23)*delta1
155 dmat(1,3) = e1*(nyu31+nyu21*nyu32)*delta1
156 dmat(2,3) = e2*(nyu32+nyu12*nyu31)*delta1
157 dmat(4,4) = g12
158 dmat(5,5) = g23
159 dmat(6,6) = g13
160
161 dmat(2,1) = dmat(1,2)
162 dmat(3,2) = dmat(2,3)
163 dmat(3,1) = dmat(1,3)
164
165 call transformation(bij, tm)
166
167 dmat = matmul( transpose(tm), dmat)
168 dmat = matmul( dmat, (tm) )
169
170 end subroutine
171
172
173 !####################################################################
175 (matl, secttype, c, &
176 e1_hat, e2_hat, e3_hat, cg1, cg2, cg3, &
177 alpha, n_layer)
178 !####################################################################
179
180 type( tmaterial ), intent(in) :: matl
181 integer, intent(in) :: sectType, n_layer
182 real(kind = kreal), intent(out) :: c(:, :, :, :)
183 real(kind = kreal), intent(in) :: e1_hat(3), e2_hat(3), e3_hat(3)
184 real(kind = kreal), intent(in) :: cg1(3), cg2(3), cg3(3)
185 real(kind = kreal), intent(out) :: alpha
186
187 !--------------------------------------------------------------------
188
189 real(kind = kreal) :: ee, pp, ee2, g12, g23, g31, theta, pp2
190 real(kind = kreal) :: outa(2)
191 real(kind = kreal) :: lambda1, lambda2, mu, k_correction
192 real(kind = kreal) :: c_hat(3, 3, 3, 3), d(5,5), d_hat(5,5), d_temp(5,5), t(5,5)
193 real(kind = kreal) :: e_hat_dot_cg(3, 3)
194 real(kind = kreal) :: alpha_over_mu
195
196 integer :: index_i(5), index_j(5), index_k(5), index_l(5)
197 integer :: is, js, ii, ij, ik, il
198 integer :: i, j, k, l, total_layer, matl_type
199
200 logical :: ierr
201
202 !--------------------------------------------------------------------
203
204 call fetch_tabledata(mc_isoelastic, matl%dict, outa, ierr)
205
206 !--------------------------------------------------------------------
207
208 matl_type = matl%shell_var(n_layer)%ortho
209 if(matl_type == 0)then
210 if( ierr ) then
211
212 ee = matl%shell_var(n_layer)%ee
213 pp = matl%shell_var(n_layer)%pp
214 alpha_over_mu = matl%variables(m_alpha_over_mu)
215
216 else
217
218 ee = outa(1)
219 pp = outa(2)
220 alpha_over_mu = matl%variables(m_alpha_over_mu)
221
222 end if
223
224 !--------------------------------------------------------------------
225
226 ! Elastic constant
227 lambda1 = ee/( 1.0d0-pp*pp )
228 lambda2 = pp*lambda1
229 mu = 0.5d0*ee/( 1.0d0+pp )
230
231 !--------------------------------------------------------------------
232
233 ! Shear correction factor
234 k_correction = 5.0d0/6.0d0
235 !k_correction = 1.0D0
236
237 !--------------------------------------------------------------------
238
239 alpha = alpha_over_mu*mu
240
241 !--------------------------------------------------------------------
242
243 ! Constitutive tensor
244 c_hat(:, :, :, :) = 0.0d0
245
246 !--------------------------------------------------------
247
248 c_hat(1, 1, 1, 1) = lambda1
249 c_hat(1, 1, 2, 2) = lambda2
250 c_hat(2, 2, 1, 1) = lambda2
251 c_hat(2, 2, 2, 2) = lambda1
252 c_hat(1, 2, 1, 2) = mu
253 c_hat(1, 2, 2, 1) = mu
254 c_hat(2, 1, 1, 2) = mu
255 c_hat(2, 1, 2, 1) = mu
256 c_hat(1, 3, 1, 3) = k_correction*mu
257 c_hat(1, 3, 3, 1) = k_correction*mu
258 c_hat(2, 3, 2, 3) = k_correction*mu
259 c_hat(2, 3, 3, 2) = k_correction*mu
260 c_hat(3, 1, 3, 1) = k_correction*mu
261 c_hat(3, 1, 1, 3) = k_correction*mu
262 c_hat(3, 2, 3, 2) = k_correction*mu
263 c_hat(3, 2, 2, 3) = k_correction*mu
264
265 !--------------------------------------------------------
266
267 elseif(matl_type == 1)then
268 total_layer = matl%totallyr
269
270 ee = matl%shell_var(n_layer)%ee
271 pp = matl%shell_var(n_layer)%pp
272 ee2 = matl%shell_var(n_layer)%ee2
273 g12 = matl%shell_var(n_layer)%g12
274 g23 = matl%shell_var(n_layer)%g23
275 g31 = matl%shell_var(n_layer)%g31
276 theta = matl%shell_var(n_layer)%angle
277
278 alpha_over_mu = matl%variables(m_alpha_over_mu)
279
280 !--------------------------------------------------------------------
281
282 ! Shear correction factor
283 k_correction = 5.0d0/6.0d0
284
285 !--------------------------------------------------------------------
286
287 alpha = alpha_over_mu * 0.5d0 * ee / ( 1.0d0+pp )
288
289 !--------------------------------------------------------------------
290
291 ! write(*,*) ee,pp,ee2,g12,g23,g31,theta,alpha_over_mu,n_layer
292
293 d(:,:) = 0.0d0
294 d_hat(:,:) = 0.0d0
295 d_temp(:,:) = 0.0d0
296 t(:,:) = 0.0d0
297
298 pp2 = pp * ee2 / ee
299
300 d(1,1) = ee / (1.0d0 - pp * pp2)
301 d(1,2) = pp2 * ee / (1.0d0- pp * pp2)
302 d(2,1) = pp2 * ee / (1.0d0- pp * pp2)
303 d(2,2) = ee2 / (1.0d0 - pp * pp2)
304 d(3,3) = g12
305 d(4,4) = g23
306 d(5,5) = g31
307
308 t(1,1) = cos(theta) * cos(theta)
309 t(1,2) = sin(theta) * sin(theta)
310 t(2,1) = sin(theta) * sin(theta)
311 t(2,2) = cos(theta) * cos(theta)
312 t(3,3) = cos(theta) * cos(theta) - sin(theta) * sin(theta)
313 t(1,3) = sin(theta) * cos(theta)
314 t(2,3) = -sin(theta) * cos(theta)
315 t(3,1) = -2.0d0 * sin(theta) * cos(theta)
316 t(3,2) = 2.0d0 * sin(theta) * cos(theta)
317 t(4,4) = cos(theta)
318 t(4,5) = sin(theta)
319 t(5,4) = -sin(theta)
320 t(5,5) = cos(theta)
321
322 !--------------------- D_temp = [D]*[T]
323
324 d_temp(1,1) = d(1,1)*t(1,1)+d(1,2)*t(2,1)
325 d_temp(1,2) = d(1,1)*t(1,2)+d(1,2)*t(2,2)
326 d_temp(2,1) = d(2,1)*t(1,1)+d(2,2)*t(2,1)
327 d_temp(2,2) = d(2,1)*t(1,2)+d(2,2)*t(2,2)
328 d_temp(3,1) = d(3,3)*t(3,1)
329 d_temp(3,2) = d(3,3)*t(3,2)
330 d_temp(1,3) = d(1,1)*t(1,3)+d(1,2)*t(2,3)
331 d_temp(2,3) = d(2,1)*t(1,3)+d(2,2)*t(2,3)
332 d_temp(3,3) = d(3,3)*t(3,3)
333 d_temp(4,4) = d(4,4)*t(4,4)
334 d_temp(4,5) = d(4,4)*t(4,5)
335 d_temp(5,4) = d(5,5)*t(5,4)
336 d_temp(5,5) = d(5,5)*t(5,5)
337
338 !--------------------- D_hat = [trans_T]*[D_temp]
339
340 d_hat(1,1) = t(1,1)*d_temp(1,1)+t(1,2)*d_temp(2,1)+t(3,1)*d_temp(3,1)
341 d_hat(1,2) = t(1,1)*d_temp(1,2)+t(1,2)*d_temp(2,2)+t(3,1)*d_temp(3,2)
342 d_hat(2,1) = t(2,1)*d_temp(1,1)+t(2,2)*d_temp(2,1)+t(3,2)*d_temp(3,1)
343 d_hat(2,2) = t(2,1)*d_temp(1,2)+t(2,2)*d_temp(2,2)+t(3,2)*d_temp(3,2)
344 d_hat(3,1) = t(1,3)*d_temp(1,1)+t(2,3)*d_temp(2,1)+t(3,3)*d_temp(3,1)
345 d_hat(3,2) = t(1,3)*d_temp(1,2)+t(2,3)*d_temp(2,2)+t(3,3)*d_temp(3,2)
346 d_hat(1,3) = t(1,1)*d_temp(1,3)+t(1,2)*d_temp(2,3)+t(3,1)*d_temp(3,3)
347 d_hat(2,3) = t(2,1)*d_temp(1,3)+t(2,2)*d_temp(2,3)+t(3,2)*d_temp(3,3)
348 d_hat(3,3) = t(1,3)*d_temp(1,3)+t(2,3)*d_temp(2,3)+t(3,3)*d_temp(3,3)
349 d_hat(4,4) = t(4,4)*d_temp(4,4)+t(5,4)*d_temp(5,4)
350 d_hat(4,5) = t(4,4)*d_temp(4,5)+t(5,4)*d_temp(5,5)
351 d_hat(5,4) = t(4,5)*d_temp(4,4)+t(5,5)*d_temp(5,4)
352 d_hat(5,5) = t(4,5)*d_temp(4,5)+t(5,5)*d_temp(5,5)
353
354 !--------------------------------------------------------------------
355
356 ! Constitutive tensor
357 c_hat(:, :, :, :) = 0.0d0
358 !write(*,*) 'Elastic linear.f90 make c_hat'
359
360 !-------D_hat to c_hat
361
362 index_i(1) = 1
363 index_i(2) = 2
364 index_i(3) = 1
365 index_i(4) = 2
366 index_i(5) = 3
367
368 index_j(1) = 1
369 index_j(2) = 2
370 index_j(3) = 2
371 index_j(4) = 3
372 index_j(5) = 1
373
374 index_k(1) = 1
375 index_k(2) = 2
376 index_k(3) = 1
377 index_k(4) = 2
378 index_k(5) = 3
379
380 index_l(1) = 1
381 index_l(2) = 2
382 index_l(3) = 2
383 index_l(4) = 3
384 index_l(5) = 1
385
386 !--------------------------------------------------------------------
387
388 do js = 1, 5
389 do is = 1, 5
390
391 ii = index_i(is)
392 ij = index_j(is)
393 ik = index_k(js)
394 il = index_l(js)
395
396 c_hat(ii, ij, ik, il) = d_hat(is, js)
397
398 end do
399 end do
400
401 !--------------------------------------------------------
402
403 else
404 write(*,*) 'shell matl type isnot collect'
405 stop
406 endif
407
408 select case( secttype )
409 case( shell )
410
411 e_hat_dot_cg(1, 1) &
412 = e1_hat(1)*cg1(1)+e1_hat(2)*cg1(2) &
413 +e1_hat(3)*cg1(3)
414 e_hat_dot_cg(2, 1) &
415 = e2_hat(1)*cg1(1)+e2_hat(2)*cg1(2) &
416 +e2_hat(3)*cg1(3)
417 e_hat_dot_cg(3, 1) = 0.0d0
418 e_hat_dot_cg(1, 2) &
419 = e1_hat(1)*cg2(1)+e1_hat(2)*cg2(2) &
420 +e1_hat(3)*cg2(3)
421 e_hat_dot_cg(2, 2) &
422 = e2_hat(1)*cg2(1)+e2_hat(2)*cg2(2) &
423 +e2_hat(3)*cg2(3)
424 e_hat_dot_cg(3, 2) = 0.0d0
425 e_hat_dot_cg(1, 3) &
426 = e1_hat(1)*cg3(1)+e1_hat(2)*cg3(2) &
427 +e1_hat(3)*cg3(3)
428 e_hat_dot_cg(2, 3) &
429 = e2_hat(1)*cg3(1)+e2_hat(2)*cg3(2) &
430 +e2_hat(3)*cg3(3)
431 e_hat_dot_cg(3, 3) &
432 = e3_hat(1)*cg3(1)+e3_hat(2)*cg3(2) &
433 +e3_hat(3)*cg3(3)
434
435 !--------------------------------------------------------
436
437 ! Constitutive tensor
438
439 c(1, 1, 1, 1) = 0.0d0
440 c(2, 2, 1, 1) = 0.0d0
441 c(1, 2, 1, 1) = 0.0d0
442 c(2, 2, 2, 2) = 0.0d0
443 c(1, 2, 2, 2) = 0.0d0
444 c(1, 2, 1, 2) = 0.0d0
445 c(3, 1, 1, 1) = 0.0d0
446 c(3, 1, 2, 2) = 0.0d0
447 c(3, 1, 1, 2) = 0.0d0
448 c(2, 3, 1, 1) = 0.0d0
449 c(2, 3, 2, 2) = 0.0d0
450 c(2, 3, 1, 2) = 0.0d0
451 c(3, 1, 3, 1) = 0.0d0
452 c(3, 1, 2, 3) = 0.0d0
453 c(2, 3, 2, 3) = 0.0d0
454
455 do l = 1, 2
456
457 do k = 1, 2
458
459 do j = 1, 2
460
461 do i = 1, 2
462
463 c(1, 1, 1, 1) &
464 = c(1, 1, 1, 1) &
465 +c_hat(i, j, k, l) &
466 *e_hat_dot_cg(i, 1)*e_hat_dot_cg(j ,1) &
467 *e_hat_dot_cg(k, 1)*e_hat_dot_cg(l, 1)
468 c(2, 2, 1, 1) &
469 = c(2, 2, 1, 1) &
470 +c_hat(i, j, k, l) &
471 *e_hat_dot_cg(i, 2)*e_hat_dot_cg(j, 2) &
472 *e_hat_dot_cg(k, 1)*e_hat_dot_cg(l, 1)
473 c(1, 2, 1, 1) &
474 = c(1, 2, 1, 1) &
475 +c_hat(i, j, k, l) &
476 *e_hat_dot_cg(i, 1)*e_hat_dot_cg(j, 2) &
477 *e_hat_dot_cg(k, 1)*e_hat_dot_cg(l, 1)
478 c(2, 2, 2, 2) &
479 = c(2, 2, 2, 2) &
480 +c_hat(i, j, k, l) &
481 *e_hat_dot_cg(i, 2)*e_hat_dot_cg(j, 2) &
482 *e_hat_dot_cg(k, 2)*e_hat_dot_cg(l, 2)
483
484 c(1, 2, 2, 2) &
485 = c(1, 2, 2, 2) &
486 +c_hat(i, j, k, l) &
487 *e_hat_dot_cg(i, 1)*e_hat_dot_cg(j, 2) &
488 *e_hat_dot_cg(k, 2)*e_hat_dot_cg(l, 2)
489 c(1, 2, 1, 2) &
490 = c(1, 2, 1, 2) &
491 +c_hat(i, j, k, l) &
492 *e_hat_dot_cg(i, 1)*e_hat_dot_cg(j, 2) &
493 *e_hat_dot_cg(k, 1)*e_hat_dot_cg(l, 2)
494
495 end do
496
497 do i = 1, 3
498
499 c(3, 1, 1, 1) &
500 = c(3, 1, 1, 1) &
501 +c_hat(i, j, k, l) &
502 *e_hat_dot_cg(i, 3)*e_hat_dot_cg(j, 1) &
503 *e_hat_dot_cg(k, 1)*e_hat_dot_cg(l, 1)
504 c(3, 1, 2, 2) &
505 = c(3, 1, 2, 2) &
506 +c_hat(i, j, k, l) &
507 *e_hat_dot_cg(i, 3)*e_hat_dot_cg(j, 1) &
508 *e_hat_dot_cg(k, 2)*e_hat_dot_cg(l, 2)
509 c(3, 1, 1, 2) &
510 = c(3, 1, 1, 2) &
511 +c_hat(i, j, k, l) &
512 *e_hat_dot_cg(i, 3)*e_hat_dot_cg(j, 1) &
513 *e_hat_dot_cg(k, 1)*e_hat_dot_cg(l, 2)
514
515 end do
516
517 end do
518
519 do j = 1, 3
520
521 do i = 1, 2
522
523 c(2, 3, 1, 1) &
524 = c(2, 3, 1, 1) &
525 +c_hat(i, j, k, l) &
526 *e_hat_dot_cg(i, 2)*e_hat_dot_cg(j, 3) &
527 *e_hat_dot_cg(k, 1)*e_hat_dot_cg(l, 1)
528 c(2, 3, 2, 2) &
529 = c(2, 3, 2, 2) &
530 +c_hat(i, j, k, l) &
531 *e_hat_dot_cg(i, 2)*e_hat_dot_cg(j, 3) &
532 *e_hat_dot_cg(k, 2)*e_hat_dot_cg(l, 2)
533 c(2, 3, 1, 2) &
534 = c(2, 3, 1, 2) &
535 +c_hat(i, j, k, l) &
536 *e_hat_dot_cg(i, 2)*e_hat_dot_cg(j, 3) &
537 *e_hat_dot_cg(k, 1)*e_hat_dot_cg(l, 2)
538
539 end do
540
541 end do
542
543 end do
544
545 do k = 1, 3
546
547 do j = 1, 2
548
549 do i = 1, 3
550
551 c(3, 1, 3, 1) &
552 = c(3, 1, 3, 1) &
553 +c_hat(i, j, k, l) &
554 *e_hat_dot_cg(i, 3)*e_hat_dot_cg(j, 1) &
555 *e_hat_dot_cg(k, 3)*e_hat_dot_cg(l, 1)
556
557 end do
558
559 end do
560
561 end do
562
563 end do
564
565 do l = 1, 3
566
567 do k = 1, 2
568
569 do j = 1, 2
570
571 do i = 1, 3
572
573 c(3, 1, 2, 3) &
574 = c(3, 1, 2, 3) &
575 +c_hat(i, j, k, l) &
576 *e_hat_dot_cg(i, 3)*e_hat_dot_cg(j, 1) &
577 *e_hat_dot_cg(k, 2)*e_hat_dot_cg(l, 3)
578
579 end do
580
581 end do
582
583 do j = 1, 3
584
585 do i = 1, 2
586
587 c(2, 3, 2, 3) &
588 = c(2, 3, 2, 3) &
589 +c_hat(i, j, k, l) &
590 *e_hat_dot_cg(i, 2)*e_hat_dot_cg(j, 3) &
591 *e_hat_dot_cg(k, 2)*e_hat_dot_cg(l, 3)
592
593 end do
594
595 end do
596
597 end do
598
599 end do
600
601 c(1, 1, 2, 2) = c(2, 2, 1, 1)
602 c(1, 1, 1, 2) = c(1, 2, 1, 1)
603 c(1, 1, 2, 3) = c(2, 3, 1, 1)
604 c(1, 1, 3, 1) = c(3, 1, 1, 1)
605 c(2, 2, 1, 2) = c(1, 2, 2, 2)
606 c(2, 2, 2, 3) = c(2, 3, 2, 2)
607 c(2, 2, 3, 1) = c(3, 1, 2, 2)
608 c(1, 2, 2, 3) = c(2, 3, 1, 2)
609 c(1, 2, 3, 1) = c(3, 1, 1, 2)
610 c(2, 3, 3, 1) = c(3, 1, 2, 3)
611
612 !--------------------------------------------------------
613
614 ! DO l = 1, 3
615 !
616 ! DO k = 1, 3
617 !
618 ! DO j = 1, 3
619 !
620 ! DO i = 1, 3
621 !
622 ! c(i, j, k, l) = 0.0D0
623 !
624 ! DO ll = 1, 3
625 !
626 ! DO kk = 1, 3
627 !
628 ! DO jj = 1, 3
629 !
630 ! DO ii = 1, 3
631 !
632 ! c(i, j, k, l) &
633 ! = c(i, j, k, l) &
634 ! +c_hat(ii, jj, kk, ll) &
635 ! *e_hat_dot_cg(ii, i)*e_hat_dot_cg(jj, j) &
636 ! *e_hat_dot_cg(kk, k)*e_hat_dot_cg(ll, l)
637 !
638 ! END DO
639 !
640 ! END DO
641 !
642 ! END DO
643 !
644 ! END DO
645 !
646 ! END DO
647 !
648 ! END DO
649 !
650 ! END DO
651 !
652 ! END DO
653
654 !--------------------------------------------------------
655
656 end select
657
658 !--------------------------------------------------------------------
659
660 return
661
662 !####################################################################
663 end subroutine linearelastic_shell
664 !####################################################################
665 ! > (Gaku Hashimoto, The University of Tokyo, 2012/11/15)
666
667
668
669end module m_elasticlinear
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
This module provides functions for elastic material.
subroutine calelasticmatrix_ortho(matl, secttype, bij, dmat, temp)
Calculate orthotropic elastic matrix.
subroutine linearelastic_shell(matl, secttype, c, e1_hat, e2_hat, e3_hat, cg1, cg2, cg3, alpha, n_layer)
subroutine calelasticmatrix(matl, secttype, d, temp)
Calculate isotropic elastic matrix.
This module provides aux functions.
Definition: utilities.f90:6
subroutine transformation(jacob, tm)
Definition: utilities.f90:339
This module summarizes all infomation of material properties.
Definition: material.f90:6
integer(kind=kint), parameter m_youngs
Definition: material.f90:84
integer(kind=kint), parameter planestress
Definition: material.f90:77
integer(kind=kint), parameter d3
Definition: material.f90:76
integer(kind=kint), parameter planestrain
Definition: material.f90:78
integer(kind=kint), parameter shell
Definition: material.f90:80
integer(kind=kint), parameter m_poisson
Definition: material.f90:85
integer(kind=kint), parameter axissymetric
Definition: material.f90:79
character(len=dict_key_length) mc_isoelastic
Definition: material.f90:119
character(len=dict_key_length) mc_orthoelastic
Definition: material.f90:120
integer(kind=kint), parameter m_alpha_over_mu
Definition: material.f90:99
Stucture to management all material relates data.
Definition: material.f90:144