-
Notifications
You must be signed in to change notification settings - Fork 6
/
Taylor.f
133 lines (133 loc) · 5.55 KB
/
Taylor.f
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
!
! Subroutine Taylor
!
!-----------------------------------------------------------------------
!-----------------------------------------------------------------------
! Subroutines should be inlined by the compiler
!-----------------------------------------------------------------------
!DIR$ ATTRIBUTES FORCEINLINE :: Taylor
!-----------------------------------------------------------------------
! Preprocessor definitions
!-----------------------------------------------------------------------
#ifndef SCMM_HYPO_TAYLOR
#define SCMM_HYPO_TAYLOR
!-----------------------------------------------------------------------
! Subroutine Taylor
!-----------------------------------------------------------------------
subroutine Taylor(stressNew,stateNew,defgradNew,
+ stressOld,stateOld,defgradOld,dt,props,
+ nblock,nstatev,nprops,Dissipation)
!-----------------------------------------------------------------------
implicit none
!-----------------------------------------------------------------------
integer nblock, nstatev, nprops
real*8 dt
real*8 props(nprops),
+ stressOld(nblock,6), stateOld(nblock,nstatev),
+ defgradNew(nblock,9), defgradOld(nblock,9),
+ stressNew(nblock,6), stateNew(nblock,nstatev)
real*8 Dissipation(nblock)! The change in dissipated inelastic
!-----------------------------------------------------------------------
integer i, j, k
integer, parameter :: Nsdv = SCMM_HYPO_NSTATEV, Ngrain = 8
real*8, parameter :: fraction = 1.d0/Ngrain
real*8 tempDissipation(nblock,Ngrain)
#if SCMM_HYPO_DFLAG != 0
real*8 isActive ! Is the integration point active (0=deleted,
! 1=active)
#endif
!-----------------------------------------------------------------------
! First step
!-----------------------------------------------------------------------
if(stateold(1,13).lt.1.d-6)then ! First step
if(nstatev.lt.((Nsdv+6)*Ngrain))then
#if defined SCMM_HYPO_STANDARD
call STDB_ABQERR(-3,
. 'The number of SDVs must be equal to %I',(Nsdv+6)*Ngrain,,)
#elif defined SCMM_HYPO_EXPLICIT
call XPLB_ABQERR(-3,
. 'The number of SDVs must be equal to %I',(Nsdv+6)*Ngrain,,)
#else
write(*,*) 'The number of SDVs must be equal to ',
. (Nsdv+6)*Ngrain
error stop 'ERROR: wrong number of SDVs'
#endif
endif
do i=1,Ngrain
do j=1,6
do k=1,nblock
stateOld(k,j+Nsdv+(i-1)*(Nsdv+6)) = stressOld(k,j)
enddo
enddo
enddo
endif
!-----------------------------------------------------------------------
! FC-Taylor homogenization for an Ngrain polycrystal
!-----------------------------------------------------------------------
do i=1,Ngrain
!-----------------------------------------------------------------------
! Call the subroutine Hypo
!-----------------------------------------------------------------------
#if SCMM_HYPO_MODEL == 3
call Hypo(stateNew(:,1+Nsdv+(i-1)*(Nsdv+6):6+Nsdv+(i-1)*(Nsdv+6)),
+ stateNew(:,1+(i-1)*(Nsdv+6):Nsdv+(i-1)*(Nsdv+6)),
+ defgradNew,
+ stateOld(:,1+Nsdv+(i-1)*(Nsdv+6):6+Nsdv+(i-1)*(Nsdv+6)),
+ stateOld(:,1+(i-1)*(Nsdv+6):Nsdv+(i-1)*(Nsdv+6)),
+ defgradOld,
+ dt,props,nblock,Nsdv,nprops,tempDissipation(:,i))
#elif SCMM_HYPO_MODEL == 4
call CCCP(stateNew(:,1+Nsdv+(i-1)*(Nsdv+6):6+Nsdv+(i-1)*(Nsdv+6)),
+ stateNew(:,1+(i-1)*(Nsdv+6):Nsdv+(i-1)*(Nsdv+6)),
+ defgradNew,
+ stateOld(:,1+Nsdv+(i-1)*(Nsdv+6):6+Nsdv+(i-1)*(Nsdv+6)),
+ stateOld(:,1+(i-1)*(Nsdv+6):Nsdv+(i-1)*(Nsdv+6)),
+ defgradOld,
+ dt,props,nblock,Nsdv,nprops,tempDissipation(:,i))
#endif
enddo
!-----------------------------------------------------------------------
! FC-Taylor homogenization
!-----------------------------------------------------------------------
stressNew = 0.d0
do i=1,Ngrain
do j=1,6
do k=1,nblock
stressNew(k,j) = stressNew(k,j) +
+ fraction*stateNew(k,j+Nsdv+(i-1)*(Nsdv+6))
enddo
enddo
enddo
Dissipation = 0.d0
do i=1,Ngrain
do k=1,nblock
Dissipation(k) = Dissipation(k) +
+ fraction*tempDissipation(k,i)
enddo
enddo
!-----------------------------------------------------------------------
! An element is deleted once the first FC-Taylor grain fails
!-----------------------------------------------------------------------
#if SCMM_HYPO_DFLAG != 0
do k=1,nblock
isActive = 1
do i=1,Ngrain
isActive = min(stateNew(k,30+(i-1)*(Nsdv+6)), isActive)
enddo
do i=1,Ngrain
stateNew(k,30+(i-1)*(Nsdv+6)) = isActive
enddo
enddo
#endif
!-----------------------------------------------------------------------
! End Subroutine
!-----------------------------------------------------------------------
return
end subroutine Taylor
!-----------------------------------------------------------------------
! End preprocessor definitions
!-----------------------------------------------------------------------
#endif
!-----------------------------------------------------------------------