FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
mechgauss.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 implicit none
10
11 ! ----------------------------------------------------------------------------
14 type(tmaterial), pointer :: pmaterial => null()
15 real(kind=kreal) :: strain(6)
16 real(kind=kreal) :: stress(6)
17 integer, pointer :: istatus(:) =>null()
18 real(kind=kreal), pointer :: fstatus(:) => null()
19 real(kind=kreal) :: plstrain
20 real(kind=kreal) :: strain_bak(6)
21 real(kind=kreal) :: stress_bak(6)
22 real(kind=kreal) :: nqm(12)
23 real(kind=kreal) :: strain_out(6)
24 real(kind=kreal) :: stress_out(6)
25 end type
26
27 ! ----------------------------------------------------------------------------
30 integer :: etype
31 integer :: iset
32 real(kind=kreal), pointer :: equiforces(:) => null()
33 type(tgaussstatus), pointer :: gausses(:) => null()
34 real(kind=kreal), pointer :: aux(:,:) => null()
35 end type
36
37contains
38
40 subroutine fstr_init_gauss( gauss )
41 type( tgaussstatus ), intent(inout) :: gauss
42 integer :: n
43 gauss%strain=0.d0; gauss%stress=0.d0
44 gauss%strain_bak=0.d0; gauss%stress_bak=0.d0
45 gauss%strain_out=0.d0; gauss%stress_out=0.d0
46 gauss%plstrain =0.d0
47 gauss%nqm =0.d0
48 if( gauss%pMaterial%mtype==usermaterial ) then
49 if( gauss%pMaterial%nfstatus> 0 ) then
50 allocate( gauss%fstatus(gauss%pMaterial%nfstatus) )
51 gauss%fstatus(:) = 0.d0
52 endif
53 else if( iselastoplastic(gauss%pMaterial%mtype) ) then
54 allocate( gauss%istatus(1) ) ! 0:elastic 1:plastic
55 if( iskinematicharden( gauss%pMaterial%mtype ) ) then
56 allocate( gauss%fstatus(7+6) ) ! plastic strain, back stress
57 else
58 allocate( gauss%fstatus(2) ) ! plastic strain
59 endif
60 gauss%istatus = 0
61 gauss%fstatus = 0.d0
62 else if( isviscoelastic(gauss%pMaterial%mtype) ) then
63 n = fetch_tablerow( mc_viscoelastic, gauss%pMaterial%dict )
64 if( n>0 ) then
65 allocate( gauss%fstatus(12*n+6) ) ! visco stress components
66 gauss%fstatus = 0.d0
67 else
68 stop "Viscoelastic properties not defined"
69 endif
70 else if( gauss%pMaterial%mtype==norton ) then
71 allocate( gauss%fstatus(2) ) ! effective stress, effective viscoplastic strain
72 gauss%fstatus = 0.d0
73 gauss%plstrain = 0.d0
74 endif
75 end subroutine fstr_init_gauss
76
78 subroutine fstr_finalize_gauss( gauss )
79 type( tgaussstatus ), intent(inout) :: gauss
80 if( associated( gauss%istatus ) ) deallocate( gauss%istatus )
81 if( associated( gauss%fstatus ) ) deallocate( gauss%fstatus )
82 end subroutine
83
85 subroutine fstr_copy_gauss( gauss1, gauss2 )
86 type( tgaussstatus ), intent(in) :: gauss1
87 type( tgaussstatus ), intent(inout) :: gauss2
88
89 gauss2%strain = gauss1%strain
90 gauss2%stress = gauss1%stress
91 gauss2%strain_bak = gauss1%strain_bak
92 gauss2%stress_bak = gauss1%stress_bak
93 gauss2%plstrain = gauss1%plstrain
94
95 if( associated(gauss1%istatus) .and. associated(gauss2%istatus) ) then
96 gauss2%istatus = gauss1%istatus
97 end if
98 if( associated(gauss1%fstatus) .and. associated(gauss2%fstatus) ) then
99 gauss2%fstatus = gauss1%fstatus
100 end if
101 end subroutine fstr_copy_gauss
102
103
104end module
105
106
107
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
This module summarizes all infomation of material properties.
Definition: material.f90:6
character(len=dict_key_length) mc_viscoelastic
Definition: material.f90:124
integer(kind=kint), parameter norton
Definition: material.f90:71
logical function iskinematicharden(mtype)
If it is a kinematic hardening material?
Definition: material.f90:318
integer(kind=kint), parameter usermaterial
Definition: material.f90:56
logical function isviscoelastic(mtype)
If it is an viscoelastic material?
Definition: material.f90:354
logical function iselastoplastic(mtype)
If it is an elastoplastic material?
Definition: material.f90:336
This modules defines a structure to record history dependent parameter in static analysis.
Definition: mechgauss.f90:6
subroutine fstr_init_gauss(gauss)
Initializer.
Definition: mechgauss.f90:41
subroutine fstr_finalize_gauss(gauss)
Finializer.
Definition: mechgauss.f90:79
subroutine fstr_copy_gauss(gauss1, gauss2)
Copy.
Definition: mechgauss.f90:86
Stucture to management all material relates data.
Definition: material.f90:144
All data should be recorded in every elements.
Definition: mechgauss.f90:29
All data should be recorded in every quadrature points.
Definition: mechgauss.f90:13