FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_ctrl_material.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
8 use mmaterial
9 use m_table
10 implicit none
11
12 private :: read_user_matl
13
14 include 'fstr_ctrl_util_f.inc'
15
16contains
17
18
19 !----------------------------------------------------------------------
21 integer function fstr_ctrl_get_material( ctrl, matname )
22 integer(kind=kint), intent(in) :: ctrl
23 character(len=*), intent(out) :: matname
24
25 matname=""
26 fstr_ctrl_get_material = fstr_ctrl_get_param_ex( ctrl, 'NAME ', '# ', 0, 'S', matname )
27 end function fstr_ctrl_get_material
28
29 !----------------------------------------------------------------------
31 integer function fstr_ctrl_get_usermaterial( ctrl, mattype, nlgeom, nstatus, matval )
32 integer(kind=kint), intent(in) :: ctrl
33 integer(kind=kint), intent(inout) :: mattype
34 integer(kind=kint), intent(out) :: nlgeom
35 integer(kind=kint), intent(out) :: nstatus
36 real(kind=kreal),intent(out) :: matval(:)
37
38 integer(kind=kint) :: ipt
39 character(len=HECMW_NAME_LEN) :: data_fmt
40 character(len=256) :: s, fname
41
43 mattype = usermaterial
44 nlgeom = updatelag !default value
45 nstatus = 1
46 if( fstr_ctrl_get_param_ex( ctrl, 'NSTATUS ', '# ', 0, 'I', nstatus )/= 0) return
47 if( fstr_ctrl_get_param_ex( ctrl, 'KIRCHHOFF ', '# ', 0, 'E', ipt )/= 0) return
48 if( ipt/=0 ) nlgeom = totallag
49
50 fstr_ctrl_get_usermaterial = read_user_matl( ctrl, matval )
52
53 !----------------------------------------------------------------------
55 integer function fstr_ctrl_get_elasticity( ctrl, mattype, nlgeom, matval, dict )
56 integer(kind=kint), intent(in) :: ctrl
57 integer(kind=kint), intent(inout) :: mattype
58 integer(kind=kint), intent(out) :: nlgeom
59 real(kind=kreal),intent(out) :: matval(:)
60 type(dict_struct), pointer :: dict
61
62 integer(kind=kint) :: i,j, rcode, depends, ipt, n
63 real(kind=kreal),pointer :: fval(:,:)
64 character(len=HECMW_NAME_LEN) :: data_fmt
65 type( ttable ) :: mattable
66 logical :: isok
67 character(len=256) :: s
68
70 depends = 0
71 rcode = fstr_ctrl_get_param_ex( ctrl, 'DEPENDENCIES ', '# ', 0, 'I', depends )
72 if( depends>1 ) depends=1 ! temperature depends only currently
73 if( depends > 3 ) stop "We cannot read dependencies>3 right now"
74 nlgeom = totallag !default value
75
76 if( fstr_ctrl_get_param_ex( ctrl, 'CAUCHY ', '# ', 0, 'E', ipt )/= 0) return
77 if( ipt/=0 ) nlgeom = updatelag
78 if( fstr_ctrl_get_param_ex( ctrl, 'INFINITESIMAL ', '# ', 0, 'E', ipt )/= 0) return
79 if( ipt/=0 ) nlgeom = infinitesimal
80
81 ! for backward compatibility
82 if( fstr_ctrl_get_param_ex( ctrl, 'INFINITE ', '# ', 0, 'E', ipt )/= 0) return
83 if( ipt/=0 ) then
84 write(*,*) "Warning : !ELASTIC : parameter 'INFINITE' is deprecated." &
85 & // " Pleaase use the replacement parameter 'INFINITESIMAL'"
86 nlgeom = infinitesimal
87 endif
88
89 ipt=1
90 s = 'ISOTROPIC,ORTHOTROPIC,USER '
91 if( fstr_ctrl_get_param_ex( ctrl, 'TYPE ', s, 0, 'P', ipt ) /= 0 ) return
92
94 ! ISOTROPIC
95 if( ipt==1 ) then
96 allocate( fval(2+depends,n) )
97 if( depends==0 ) then
98 data_fmt = "RR "
100 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:) )
101 endif
102 if( depends==1 ) then
103 data_fmt = "RRR "
105 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:) )
106 endif
107 if( fstr_ctrl_get_elasticity ==0 ) then
108 matval(m_youngs) = fval(1,1)
109 matval(m_poisson) = fval(2,1)
110 call init_table( mattable, depends, 2+depends,n, fval )
111 call dict_add_key( dict, mc_isoelastic, mattable )
112 ! call print_table( mattable, 6 ); pause
113 endif
114 mattype = elastic
115
116 ! ORTHOTROPIC
117 else if( ipt==2 ) then
118 allocate( fval(9+depends,n) )
119 if( depends==0 ) then
120 data_fmt = "RRRRRRRRR "
122 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, &
123 fval(1,:), fval(2,:), fval(3,:), fval(4,:), fval(5,:), fval(6,:), &
124 fval(7,:), fval(8,:), fval(9,:) )
125 else if( depends==1 ) then
126 data_fmt = "RRRRRRRRRR "
128 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, &
129 fval(1,:), fval(2,:), fval(3,:), fval(4,:), fval(5,:), fval(6,:), &
130 fval(7,:), fval(8,:), fval(9,:), fval(10,:) )
131 endif
132 if( fstr_ctrl_get_elasticity ==0 ) then
133 isok = .true.
134 do i=1,n
135 if( fval(1,i)<=0.d0 .or. fval(2,i)<=0.d0 .or. fval(3,i)<=0.d0 .or. &
136 fval(7,i)<=0.d0 .or. fval(8,i)<=0.d0 .or. fval(9,i)<=0.d0 ) then
137 isok = .false.; fstr_ctrl_get_elasticity=1; exit
138 endif
139 enddo
140 if( isok ) then
141 call init_table( mattable, depends, 9+depends,n, fval )
142 call dict_add_key( dict, mc_orthoelastic, mattable )
143 mattype = mn_orthoelastic
144 endif
145 endif
146
147 else if( ipt==3 ) then
148 allocate( fval(10,10) )
149 fval =0.d0
150 fstr_ctrl_get_elasticity = fstr_ctrl_get_data_ex( ctrl, 1, 'rrrrrrrrrr ', &
151 fval(1,:), fval(2,:), fval(3,:), fval(4,:), fval(5,:), fval(6,:), &
152 fval(7,:), fval(8,:), fval(9,:), fval(10,:) )
153 if( fstr_ctrl_get_elasticity ==0 ) then
154 do i=1,10
155 do j=1,10
156 matval(100+(i-1)*10+j)=fval(i,j)
157 enddo
158 enddo
159 endif
160 mattype = userelastic
161 nlgeom = infinitesimal
162
163 else
164 stop "ERROR: Material type not supported"
165
166 endif
167
168 call finalize_table( mattable )
169 if( associated(fval) ) deallocate(fval)
170 end function fstr_ctrl_get_elasticity
171
172
173 !----------------------------------------------------------------------
175 integer function fstr_ctrl_get_hyperelastic( ctrl, mattype, nlgeom, matval )
176 integer(kind=kint), intent(in) :: ctrl
177 integer(kind=kint), intent(inout) :: mattype
178 integer(kind=kint), intent(out) :: nlgeom
179 real(kind=kreal),intent(out) :: matval(:)
180
181 integer(kind=kint) :: i,j, rcode, depends, ipt
182 real(kind=kreal),pointer :: fval(:,:)
183 character(len=HECMW_NAME_LEN) :: data_fmt
184 character(len=256) :: s
185
187 depends = 0
188 rcode = fstr_ctrl_get_param_ex( ctrl, 'DEPENDENCIES ', '# ', 0, 'I', depends )
189 if( depends > 3 ) stop "We cannot read dependencies>3 right now"
190 nlgeom = totallag !default value
191 if( fstr_ctrl_get_param_ex( ctrl, 'CAUCHY ', '# ', 0, 'E', ipt )/= 0) return
192 if( ipt/=0 ) nlgeom = updatelag
193
194 ipt=1
195 s = 'NEOHOOKE,MOONEY-RIVLIN,ARRUDA-BOYCE,USER,MOONEY-RIVLIN-ANISO '
196 if( fstr_ctrl_get_param_ex( ctrl, 'TYPE ', s, 0, 'P', ipt ) /= 0 ) return
197
198 ! NEOHOOKE
199 if( ipt==1 ) then
200 allocate( fval(2,depends+1) )
201 fval =0.0d0
202 if( depends==0 ) then
203 data_fmt = "RR "
205 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:) )
206 endif
207 if( fstr_ctrl_get_hyperelastic ==0 ) then
208 if( fval(2,1)==0.d0 ) stop "We cannot deal with incompressible material currently"
209 matval(m_plconst1) = fval(1,1)
210 matval(m_plconst2) = 0.d0
211 matval(m_plconst3) = fval(2,1)
212 ! matval(M_YOUNGS) = fval(1,1)
213 ! matval(M_POISSON) = fval(2,1)
214 endif
215 mattype = neohooke
216
217 ! MOONEY
218 else if( ipt==2 ) then
219 allocate( fval(3,depends+1) )
220 fval =0.0d0
221 if( depends==0 ) then
222 data_fmt = "RRR "
224 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:) )
225 endif
226 if( fstr_ctrl_get_hyperelastic ==0 ) then
227 matval(m_plconst1) = fval(1,1)
228 matval(m_plconst2) = fval(2,1)
229 matval(m_plconst3) = fval(3,1)
230 endif
231 mattype = mooneyrivlin
232
233 ! ARRUDA
234 else if( ipt==3 ) then
235 allocate( fval(3,depends+1) )
236 fval =0.0d0
237 if( depends==0 ) then
238 data_fmt = "RRR "
240 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:) )
241 endif
242 if( fstr_ctrl_get_hyperelastic ==0 ) then
243 matval(m_plconst1) = fval(1,1)
244 matval(m_plconst2) = fval(2,1)
245 matval(m_plconst3) = fval(3,1)
246 endif
247 mattype = arrudaboyce
248
249 else if( ipt==4 ) then !User
250 allocate( fval(10,10) )
251 fval =0.0d0
252 fstr_ctrl_get_hyperelastic = fstr_ctrl_get_data_ex( ctrl, 1, 'rrrrrrrrrr ', &
253 fval(1,:), fval(2,:), fval(3,:), fval(4,:), fval(5,:), fval(6,:), &
254 fval(7,:), fval(8,:), fval(9,:), fval(10,:) )
255 if( fstr_ctrl_get_hyperelastic ==0 ) then
256 do i=1,10
257 do j=1,10
258 matval(100+(i-1)*10+j)=fval(i,j)
259 enddo
260 enddo
261 endif
262 mattype = userhyperelastic
263
264 ! MOONEY-ORTHO
265 else if( ipt==5 ) then
266 allocate( fval(10,depends+1) )
267 fval =0.0d0
268 if( depends==0 ) then
269 data_fmt = "RRRRRrrrrr "
271 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, &
272 fval(1,:), fval(2,:), fval(3,:), fval(4,:), fval(5,:), &
273 fval(6,:), fval(7,:), fval(8,:), fval(9,:), fval(10,:) )
274 endif
275 if( fstr_ctrl_get_hyperelastic ==0 ) then
276 matval(m_plconst1) = fval(1,1)
277 matval(m_plconst2) = fval(2,1)
278 matval(m_plconst3) = fval(3,1)
279 matval(m_plconst4) = fval(4,1)
280 matval(m_plconst5) = fval(5,1)
281 matval(m_plconst6) = fval(6,1)
282 matval(m_plconst7) = fval(7,1)
283 matval(m_plconst8) = fval(8,1)
284 matval(m_plconst9) = fval(9,1)
285 matval(m_plconst10) = fval(10,1)
286 endif
287 mattype = mooneyrivlin_aniso
288
289 endif
290
291 if( associated(fval) ) deallocate(fval)
292 end function fstr_ctrl_get_hyperelastic
293
294
295 !----------------------------------------------------------------------
297 integer function fstr_ctrl_get_viscoelasticity( ctrl, mattype, nlgeom, dict )
298 integer(kind=kint), intent(in) :: ctrl
299 integer(kind=kint), intent(inout) :: mattype
300 integer(kind=kint), intent(out) :: nlgeom
301 type(dict_struct), pointer :: dict
302
303 integer(kind=kint) :: i,j, rcode, depends, ipt, n
304 real(kind=kreal),pointer :: fval(:,:)
305 character(len=HECMW_NAME_LEN) :: data_fmt
306 type( ttable ) :: mattable
307 character(len=256) :: s
308
310 depends = 0
311 rcode = fstr_ctrl_get_param_ex( ctrl, 'DEPENDENCIES ', '# ', 0, 'I', depends )
312 if( depends>1 ) depends=1 ! temperature depends only currently
313 !depends = 0
314 nlgeom = totallag !default value
315 if( fstr_ctrl_get_param_ex( ctrl, 'INFINITESIMAL ', '# ', 0, 'E', ipt )/= 0) return
316 if( ipt/=0 ) nlgeom = infinitesimal
317
318 ! for backward compatibility
319 if( fstr_ctrl_get_param_ex( ctrl, 'INFINITE ', '# ', 0, 'E', ipt )/= 0) return
320 if( ipt/=0 ) then
321 write(*,*) "Warning : !VISCOELASTIC : parameter 'INFINITE' is deprecated." &
322 & // " Pleaase use the replacement parameter 'INFINITESIMAL'"
323 nlgeom = infinitesimal
324 endif
325
326 ipt=1
327 s = 'ISOTROPIC,USER '
328 if( fstr_ctrl_get_param_ex( ctrl, 'TYPE ', s, 0, 'P', ipt ) /= 0 ) return
329 ipt = 1
330
331 ! ISOTROPIC
332 if( ipt==1 ) then
333 n = fstr_ctrl_get_data_line_n( ctrl )
334 allocate( fval(2+depends,n) )
335 if( depends==0 ) then
336 data_fmt = "RR "
338 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:) )
339 if( fval(2,1)==0.d0 ) stop "Error in defining viscoelasticity: Relaxation time cannot be zero!"
340 endif
341 if( depends==1 ) then
342 data_fmt = "RRR "
344 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:) )
345 endif
346 if( fstr_ctrl_get_viscoelasticity ==0 ) then
347 call init_table( mattable, 1, 2+depends,n, fval )
348 call dict_add_key( dict, mc_viscoelastic, mattable )
349 ! call print_table( mattable, 6 ); pause
350 endif
351 mattype = viscoelastic
352
353 else
354 stop "ERROR: Material type not supported"
355
356 endif
357
358 call finalize_table( mattable )
359 if( associated(fval) ) deallocate(fval)
361
362
364 integer function fstr_ctrl_get_trs( ctrl, mattype, matval )
365 integer(kind=kint), intent(in) :: ctrl
366 integer(kind=kint), intent(inout) :: mattype
367 real(kind=kreal),intent(out) :: matval(:)
368
369 integer :: ipt
370 character(len=256) :: s
371
372 ipt=1
373 s = 'WLF,ARRHENIUS '
374 if( fstr_ctrl_get_param_ex( ctrl, 'DEFINITION ', s, 0, 'P', ipt ) /= 0 ) return
375
377 = fstr_ctrl_get_data_ex( ctrl, 1, "RRR ", matval(1), matval(2), matval(3) )
378 if( fstr_ctrl_get_trs/=0 ) return
379 mattype = mattype+ipt
380
381 end function fstr_ctrl_get_trs
382
383
384 !----------------------------------------------------------------------
386 integer function fstr_ctrl_get_plasticity( ctrl, mattype, nlgeom, matval, mattable, dict )
387 integer(kind=kint), intent(in) :: ctrl
388 integer(kind=kint), intent(inout) :: mattype
389 integer(kind=kint), intent(out) :: nlgeom
390 real(kind=kreal),intent(out) :: matval(:)
391 real(kind=kreal), pointer :: mattable(:)
392 type(dict_struct), pointer :: dict
393
394 integer(kind=kint) :: i, n, rcode, depends, ipt, hipt
395 real(kind=kreal),pointer :: fval(:,:)
396 real(kind=kreal) :: dum, fdum
397 character(len=HECMW_NAME_LEN) :: data_fmt
398 character(len=256) :: s
399 type( ttable ) :: mttable
400 real(kind=kreal), parameter :: pi=3.14159265358979d0
401
403 ipt = 0; hipt = 0
404
405 depends = 0
406 rcode = fstr_ctrl_get_param_ex( ctrl, 'DEPENDENCIES ', '# ', 0, 'I', depends )
407 if( depends>1 ) depends = 1 ! we consider temprature dependence only currently
408 if( depends > 3 ) stop "We cannot read dependencies>3 right now"
409 nlgeom = updatelag !default value
410 if( fstr_ctrl_get_param_ex( ctrl, 'KIRCHHOFF ', '# ', 0, 'E', ipt )/= 0) return
411 ! rcode = fstr_ctrl_get_param_ex( ctrl, 'FILE ', '# ', 0, 'S', fname )
412 if( ipt/=0 ) nlgeom = totallag
413 if( fstr_ctrl_get_param_ex( ctrl, 'INFINITESIMAL ', '# ', 0, 'E', ipt )/= 0) return
414 if( ipt/=0 ) nlgeom = infinitesimal
415
416 ! for backward compatibility
417 if( fstr_ctrl_get_param_ex( ctrl, 'INFINITE ', '# ', 0, 'E', ipt )/= 0) return
418 if( ipt/=0 ) then
419 write(*,*) "Warning : !PLASTIC : parameter 'INFINITE' is deprecated." &
420 & // " Pleaase use the replacement parameter 'INFINITESIMAL'"
421 nlgeom = infinitesimal
422 endif
423
424 call setdigit( 1, 1, mattype )
425 call setdigit( 2, 2, mattype )
426
427 ! hardening
428 s = 'BILINEAR,MULTILINEAR,SWIFT,RAMBERG-OSGOOD,KINEMATIC,COMBINED '
429 if( fstr_ctrl_get_param_ex( ctrl, 'HARDEN ', s , 0, 'P', hipt ) /= 0 ) return
430 if( hipt==0 ) hipt=1 ! default: linear hardening
431 call setdigit( 5, hipt-1, mattype )
432
433 ! yield function
434 s = 'MISES,MOHR-COULOMB,DRUCKER-PRAGER,USER '
435 call setdigit( 2, 2, mattype )
436 if( fstr_ctrl_get_param_ex( ctrl, 'YIELD ', s , 0, 'P', ipt ) /= 0 ) return
437 if( ipt==0 ) ipt=1 ! default: mises yield function
438 call setdigit( 4, ipt-1, mattype )
439
440 n = fstr_ctrl_get_data_line_n( ctrl )
441 if( n == 0 ) return ! fail in reading plastic
442 if( hipt==2 .and. n<2 ) return ! not enought data
443 if( ( ipt==3 .or. ipt==4 ) .and. hipt>2 ) hipt = 1
444
445 select case (ipt)
446 case (1) !Mises
447 select case (hipt)
448 case (1,5) ! linear hardening, kinematic hardening
449 allocate( fval(2,n) )
450 data_fmt = "RR "
452 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:) )
453 if( fstr_ctrl_get_plasticity ==0 ) then
454 matval(m_plconst1) = fval(1,1)
455 if(hipt==1) then
456 matval(m_plconst2) = fval(2,1)
457 else
458 matval(m_plconst2) = 0.d0
459 matval(m_plconst3) = fval(2,1)
460 endif
461 endif
462 case (2) ! multilinear approximation
463 allocate( fval(depends+2,n) )
464 if( depends==0 ) then
465 data_fmt = "RR "
467 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:) )
468 if( fstr_ctrl_get_plasticity ==0 ) then
469 if( fval(2,1)/=0.d0 ) then
470 print *, "Multilinear hardening: First plastic strain must be zero"
471 stop
472 endif
473 do i=1,n
474 if( fval(2,i)<0.0 ) &
475 stop "Multilinear hardening: Error in plastic strain definition"
476 enddo
477 call init_table( mttable,1, 2+depends, n, fval )
478 call dict_add_key( dict, mc_yield, mttable )
479
480 endif
481 else ! depends==1
482 data_fmt = "RRR "
484 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:) )
485 if( fstr_ctrl_get_plasticity ==0 ) then
486 call init_table( mttable,2, 2+depends,n, fval )
487 call dict_add_key( dict, mc_yield, mttable )
488 endif
489 endif
490 case (3, 4, 6) ! swift, Ramberg-Osgood, Combined
491 allocate( fval(3,1) )
492 data_fmt = "RRR "
494 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:) )
495 if( fstr_ctrl_get_plasticity ==0 ) then
496 matval(m_plconst1) = fval(1,1)
497 matval(m_plconst2) = fval(2,1)
498 matval(m_plconst3) = fval(3,1)
499 endif
500 case default
501 print *, "Error in hardening definition!"
502 stop
503 end select
504 case (2, 3) ! Mohr-Coulomb, Drucker-Prager
505 call setdigit( 5, 0, mattype )
506 allocate( fval(3,depends+1) )
507 data_fmt = "RRr "
509 = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:) )
510 if( fstr_ctrl_get_plasticity ==0 ) then
511 matval(m_plconst1) = fval(1,1) ! c
512 matval(m_plconst2) = fval(3,1) ! H
513 if( ipt==3 ) then ! Drucker-Prager
514 dum = fval(2,1)*pi/180.d0
515 fdum = 2.d0*sin(dum)/ ( sqrt(3.d0)*(3.d0+sin(dum)) )
516 matval(m_plconst3) = fdum
517 fdum = 6.d0* cos(dum)/ ( sqrt(3.d0)*(3.d0+sin(dum)) )
518 matval(m_plconst4) = fdum
519 else ! Mohr-Coulomb
520 matval(m_plconst3) = fval(2,1)*pi/180.d0
521 endif
522 endif
523
524 case(4)
525 fstr_ctrl_get_plasticity = read_user_matl( ctrl, matval )
526
527 case default
528 stop "Yield function not supported"
529 end select
530
531 if( associated(fval) ) deallocate(fval)
532 call finalize_table( mttable )
533 end function fstr_ctrl_get_plasticity
534
535
536 !----------------------------------------------------------------------
538 integer function fstr_ctrl_get_viscoplasticity( ctrl, mattype, nlgeom, dict )
539 integer(kind=kint), intent(in) :: ctrl
540 integer(kind=kint), intent(inout) :: mattype
541 integer(kind=kint), intent(out) :: nlgeom
542 type(dict_struct), pointer :: dict
543
544 integer(kind=kint) :: i,j, rcode, depends, ipt, n
545 real(kind=kreal),pointer :: fval(:,:)
546 character(len=HECMW_NAME_LEN) :: data_fmt
547 type( ttable ) :: mattable
548 character(len=256) :: s
549
551 depends = 0
552 rcode = fstr_ctrl_get_param_ex( ctrl, 'DEPENDENCIES ', '# ', 0, 'I', depends )
553 if( depends>1 ) depends=1 ! temperature depends only currently
554 depends = 0
555 nlgeom = updatelag !default value
556 if( fstr_ctrl_get_param_ex( ctrl, 'KIRCHHOFF ', '# ', 0, 'E', ipt )/= 0) return
557 if( ipt/=0 ) nlgeom = totallag
558
559 ipt=1
560 s = 'NORTON,USER '
561 if( fstr_ctrl_get_param_ex( ctrl, 'TYPE ', s, 0, 'P', ipt ) /= 0 ) return
562 ipt = 1
563
564 ! NORTON
565 if( ipt==1 ) then
566 n = fstr_ctrl_get_data_line_n( ctrl )
567 allocate( fval(3+depends,n) )
568 if( depends==0 ) then
569 data_fmt = "RRR "
571 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:) )
572 endif
573 if( depends==1 ) then
574 data_fmt = "RRRR "
576 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:), fval(4,:) )
577 endif
578 if( fstr_ctrl_get_viscoplasticity ==0 ) then
579 call init_table( mattable, depends, 3+depends,n, fval )
580 call dict_add_key( dict, mc_norton, mattable )
581 ! call print_table( mattable, 6 ); pause
582 endif
583 mattype = norton
584
585 else
586 stop "ERROR: Material type not supported"
587
588 endif
589
590 call finalize_table( mattable )
591 if( associated(fval) ) deallocate(fval)
593
594 !----------------------------------------------------------------------
596 integer function fstr_ctrl_get_density( ctrl, matval )
597 integer(kind=kint), intent(in) :: ctrl
598 real(kind=kreal),intent(out) :: matval(:)
599
600 integer(kind=kint) :: i, rcode, depends
601 real(kind=kreal),pointer :: fval(:,:)
602 character(len=HECMW_NAME_LEN) :: data_fmt
603
604 data_fmt = "R "
605
607
608 depends = 0
609 rcode = fstr_ctrl_get_param_ex( ctrl, 'DEPENDENCIES ', '# ', 0, 'I', depends )
610 if( depends>1 ) depends = 1 ! we consider temprature dependence only currently
611
612 allocate( fval(1,depends+1) )
613 do i=2,1+depends
614 data_fmt = data_fmt //"R "
615 enddo
617 = fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:) )
618 if( fstr_ctrl_get_density==0 ) matval(m_density) = fval(1,1)
619
620 if( associated(fval) ) deallocate(fval)
621
622 end function fstr_ctrl_get_density
623
624
625 !----------------------------------------------------------------------
627 integer function fstr_ctrl_get_expansion_coeff( ctrl, matval, dict )
628 integer(kind=kint), intent(in) :: ctrl
629 real(kind=kreal),intent(out) :: matval(:)
630 type(dict_struct), pointer :: dict
631
632 integer(kind=kint) :: i, n, rcode, depends, ipt
633 real(kind=kreal),pointer :: fval(:,:)
634 type( ttable ) :: mttable
635 character(len=HECMW_NAME_LEN) :: data_fmt, ss
636
637 data_fmt = "R "
638
640 n = fstr_ctrl_get_data_line_n( ctrl )
641 if( n == 0 ) return ! fail in reading plastic
642
643 ss = 'ISOTROPIC,ORTHOTROPIC '
644 ipt = 1 !default
645 if( fstr_ctrl_get_param_ex( ctrl, 'TYPE ', ss, 0, 'P', ipt ) /= 0 ) return
646
647 depends = 0
648 rcode = fstr_ctrl_get_param_ex( ctrl, 'DEPENDENCIES ', '# ', 0, 'I', depends )
649 if( depends>1 ) depends = 1 ! we consider temprature dependence only currently
650
651 if( ipt==1 ) then
652 allocate( fval(depends+1, n) )
653 do i=2,1+depends
654 data_fmt = data_fmt //"R "
655 enddo
656 if( depends==0 ) then
658 fstr_ctrl_get_data_array_ex( ctrl, "R ", fval(1,:) )
659 else
661 fstr_ctrl_get_data_array_ex( ctrl, "RR ", fval(1,:), fval(2,:) )
662 endif
664 matval(m_exapnsion) = fval(1,1)
665 call init_table( mttable,depends, 1+depends, n, fval )
666 call dict_add_key( dict, mc_themoexp, mttable )
667 endif
668 else
669 allocate( fval(3+depends,n) )
670 do i=2,3+depends
671 data_fmt = trim(data_fmt) //"R "
672 enddo
673 if( depends==0 ) then
675 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:) )
676 elseif( depends==1 ) then
678 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:), fval(3,:), fval(4,:) )
679 endif
681 call init_table( mttable, depends, 3+depends,n, fval )
682 if( fstr_ctrl_get_expansion_coeff==0 ) call dict_add_key( dict, mc_orthoexp, mttable )
683 endif
684 endif
685
686 call finalize_table( mttable )
687 if( associated(fval) ) deallocate(fval)
689
690
691 integer function read_user_matl( ctrl, matval )
692 integer(kind=kint), intent(in) :: ctrl
693 real(kind=kreal),intent(out) :: matval(:)
694
695 integer(kind=kint) :: i, j
696 real(kind=kreal) :: fval(10,10)
697
698 read_user_matl = -1
699
700 fval =0.d0
701 if( fstr_ctrl_get_data_array_ex( ctrl, 'rrrrrrrrrr ', fval(1,:), fval(2,:), fval(3,:), &
702 fval(4,:), fval(5,:), fval(6,:), fval(7,:), fval(8,:), fval(9,:), fval(10,:) ) /= 0 ) return
703 do i=1,10
704 do j=1,10
705 matval((i-1)*10+j)=fval(i,j)
706 enddo
707 enddo
708
709 read_user_matl = 0
710 end function read_user_matl
711
712 !----------------------------------------------------------------------
714 integer function fstr_ctrl_get_fluid( ctrl, mattype, nlgeom, matval, dict )
715 integer(kind=kint), intent(in) :: ctrl
716 integer(kind=kint), intent(inout) :: mattype
717 integer(kind=kint), intent(out) :: nlgeom
718 real(kind=kreal),intent(out) :: matval(:)
719 type(dict_struct), pointer :: dict
720
721 integer(kind=kint) :: i,j, rcode, depends, ipt, n
722 real(kind=kreal),pointer :: fval(:,:)
723 character(len=HECMW_NAME_LEN) :: data_fmt
724 type( ttable ) :: mattable
725 logical :: isok
726 character(len=256) :: s
727
729 depends = 0
730 rcode = fstr_ctrl_get_param_ex( ctrl, 'DEPENDENCIES ', '# ', 0, 'I', depends )
731 if( depends>1 ) depends=1 ! temperature depends only currently
732 if( depends > 3 ) stop "We cannot read dependencies>3 right now"
733 nlgeom = totallag !default value
734
735 ipt=1
736 s = 'INCOMP_NEWTONIAN '
737 if( fstr_ctrl_get_param_ex( ctrl, 'TYPE ', s, 0, 'P', ipt ) /= 0 ) return
738
739 n = fstr_ctrl_get_data_line_n( ctrl )
740 ! ISOTROPIC
741 if( ipt==1 ) then
742 allocate( fval(1+depends,n) )
743 if( depends==0 ) then
744 data_fmt = "R "
746 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:) )
747 endif
748 if( depends==1 ) then
749 data_fmt = "RR "
751 fstr_ctrl_get_data_array_ex( ctrl, data_fmt, fval(1,:), fval(2,:) )
752 endif
753 if( fstr_ctrl_get_fluid ==0 ) then
754 matval(m_viscocity) = fval(1,1)
755 call init_table( mattable, depends, 1+depends,n, fval )
756 call dict_add_key( dict, mc_incomp_newtonian, mattable )
757 ! call print_table( mattable, 6 ); pause
758 endif
759 mattype = incomp_newtonian
760
761 else
762 stop "ERROR: Material type not supported"
763
764 endif
765
766 call finalize_table( mattable )
767 if( associated(fval) ) deallocate(fval)
768 end function fstr_ctrl_get_fluid
769
770end module fstr_ctrl_material
int fstr_ctrl_get_param_ex(int *ctrl, const char *param_name, const char *value_list, int *necessity, char *type, void *val)
int fstr_ctrl_get_data_line_n(int *ctrl)
int fstr_ctrl_get_data_array_ex(int *ctrl, const char *format,...)
int fstr_ctrl_get_data_ex(int *ctrl, int *line_no, const char *format,...)
This module manages read in of various material properties.
integer function fstr_ctrl_get_hyperelastic(ctrl, mattype, nlgeom, matval)
Read in !HYPERELASTIC.
integer function fstr_ctrl_get_viscoelasticity(ctrl, mattype, nlgeom, dict)
Read in !VISCOELASTIC.
integer function fstr_ctrl_get_viscoplasticity(ctrl, mattype, nlgeom, dict)
Read in !VISCOELASTIC.
integer function fstr_ctrl_get_usermaterial(ctrl, mattype, nlgeom, nstatus, matval)
Read in !USER_MATERIAL.
integer function fstr_ctrl_get_expansion_coeff(ctrl, matval, dict)
Read in !EXPANSION_COEFF.
integer function fstr_ctrl_get_trs(ctrl, mattype, matval)
Read in !TRS.
integer function fstr_ctrl_get_elasticity(ctrl, mattype, nlgeom, matval, dict)
Read in !ELASTIC.
integer function fstr_ctrl_get_plasticity(ctrl, mattype, nlgeom, matval, mattable, dict)
Read in !PLASTIC.
integer function fstr_ctrl_get_material(ctrl, matname)
Read in !MATERIAL.
integer function fstr_ctrl_get_density(ctrl, matval)
Read in !DENSITY.
integer function fstr_ctrl_get_fluid(ctrl, mattype, nlgeom, matval, dict)
Read in !FLUID.
Definition: hecmw.f90:6
This module provides data structure table which would be dictionaried afterwards.
Definition: ttable.f90:7
subroutine finalize_table(table)
Definition: ttable.f90:81
subroutine init_table(table, ndp, col, row, tbval)
Definition: ttable.f90:43
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 m_plconst5
Definition: material.f90:94
character(len=dict_key_length) mc_incomp_newtonian
Definition: material.f90:126
integer(kind=kint), parameter m_plconst6
Definition: material.f90:112
integer(kind=kint), parameter mooneyrivlin
Definition: material.f90:65
integer(kind=kint), parameter mooneyrivlin_aniso
Definition: material.f90:68
character(len=dict_key_length) mc_viscoelastic
Definition: material.f90:124
integer(kind=kint), parameter m_plconst4
Definition: material.f90:93
integer(kind=kint), parameter viscoelastic
Definition: material.f90:70
integer(kind=kint), parameter m_exapnsion
Definition: material.f90:97
integer(kind=kint), parameter m_plconst1
Definition: material.f90:90
character(len=dict_key_length) mc_themoexp
Definition: material.f90:122
character(len=dict_key_length) mc_norton
Definition: material.f90:125
integer(kind=kint), parameter m_plconst10
Definition: material.f90:116
integer(kind=kint), parameter elastic
Definition: material.f90:58
integer(kind=kint), parameter m_plconst2
Definition: material.f90:91
integer(kind=kint), parameter arrudaboyce
Definition: material.f90:66
integer(kind=kint), parameter mn_orthoelastic
Definition: material.f90:59
integer(kind=kint), parameter incomp_newtonian
Definition: material.f90:73
integer(kind=kint), parameter totallag
Definition: material.f90:14
integer(kind=kint), parameter m_density
Definition: material.f90:86
subroutine setdigit(npos, ival, mtype)
Modify material type.
Definition: material.f90:252
integer(kind=kint), parameter norton
Definition: material.f90:71
integer(kind=kint), parameter m_plconst9
Definition: material.f90:115
integer(kind=kint), parameter m_poisson
Definition: material.f90:85
integer(kind=kint), parameter m_plconst8
Definition: material.f90:114
character(len=dict_key_length) mc_orthoexp
Definition: material.f90:123
integer(kind=kint), parameter infinitesimal
Definition: material.f90:13
character(len=dict_key_length) mc_yield
Definition: material.f90:121
integer(kind=kint), parameter userelastic
Definition: material.f90:60
integer(kind=kint), parameter m_viscocity
Definition: material.f90:109
integer(kind=kint), parameter neohooke
Definition: material.f90:64
integer(kind=kint), parameter m_plconst7
Definition: material.f90:113
integer(kind=kint), parameter usermaterial
Definition: material.f90:56
character(len=dict_key_length) mc_isoelastic
Definition: material.f90:119
integer(kind=kint), parameter userhyperelastic
Definition: material.f90:67
integer(kind=kint), parameter m_plconst3
Definition: material.f90:92
character(len=dict_key_length) mc_orthoelastic
Definition: material.f90:120
integer(kind=kint), parameter updatelag
Definition: material.f90:15