40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62#include "implicit_f.inc"
63#include "comlock.inc"
64
65
66
67#include "param_c.inc"
68#include "com04_c.inc"
69#include "com06_c.inc"
70#include "com08_c.inc"
71#include "vect01_c.inc"
72#include "scr06_c.inc"
73
74
75
76 INTEGER MAT(NEL), IFLAG, NEL
78 . off(nel) ,eint(nel) ,mu(nel) ,
79 . mu2(nel) ,espe(nel) ,dvol(nel) ,df(nel) ,
80 . vnew(nel) ,pnew(nel) ,dpdm(nel),
81 . dpde(nel)
82 my_real,
INTENT(INOUT) :: psh(nel)
83
84
85
86 INTEGER I,MX,IBFRAC
88 . rho0 , aa , bb , r1,
89 . r2, w1, bhe, p0, bulk,
90 . r1df,r2df,er1df,er2df,dpdmu
91
92
93
94
95 IF(iflag == 0) THEN
96
97 ELSEIF(iflag == 1) THEN
98
99 ELSEIF(iflag == 2) THEN
100 mx = mat(1)
101 rho0 = pm( 1,mx)
102 aa = pm(33,mx)
103 bb = pm(34,mx)
104 r1 = pm(35,mx)
105 r2 = pm(36,mx)
106 w1 = pm(45,mx)
107 vdet = pm(38,mx)
108 bhe = pm(40,mx)
109 p0 = pm(31,mx)
110 bulk = pm(44,mx)
111 ibfrac = nint(pm(41,mx))
112 psh(1:nel) = pm(88,mx)
113 DO i=1,nel
114 bfrac(i)=one
115 ENDDO
116 DO i=1, nel
117 IF (vnew(i) > zero) THEN
118 r1df = r1*df(i)
119 r2df = r2*df(i)
120 er1df = exp(-r1df)
121 er2df = exp(-r2df)
122 pnew(i) = - psh(i) + aa*(one-w1/r1df)*er1df + bb*(one-w1/r2df)*er2df + w1*espe(i)/df(i)
123 pnew(i) =
max(zero - psh(i), pnew(i))
124
125
126
127 dpde(i) = w1/df(i)
128 dpdmu =
129 . -aa*w1*er1df/r1 + aa*(one-w1/(r1df))*r1df*r1df*er1df
130 . -bb*w1*er2df/r2 + bb*(one-w1/(r2df))*r2df*r2df*er2df
131 . +w1*espe(i)
132 dpdm(i) = dpdmu + (pnew(i)+psh(i))*df(i)*df(i)*dpde(i)
133 ENDIF
134 ENDDO
135 ENDIF
136
137
138 RETURN