OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
no_asso_qplas76c.F File Reference
#include "implicit_f.inc"
#include "tabsiz_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine no_asso_qplas76c (nel, nuparam, nuvar, nfunc, ifunc, npf, tf, numtabl, table, time, timestep, uparam, uvar, rho, depsxx, depsyy, depsxy, depsyz, depszx, sigoxx, sigoyy, sigoxy, sigoyz, sigozx, signxx, signyy, signxy, signyz, signzx, pla, dpla, epsd, off, gs, yld, soundsp, dezz, inloc, dplanl, nvartmp, vartmp, loff)

Function/Subroutine Documentation

◆ no_asso_qplas76c()

subroutine no_asso_qplas76c ( integer, intent(in) nel,
integer, intent(in) nuparam,
integer, intent(in) nuvar,
integer, intent(in) nfunc,
integer, dimension(nfunc), intent(in) ifunc,
integer, dimension(snpc) npf,
tf,
integer, intent(in) numtabl,
type(table_4d_), dimension(numtabl), intent(in) table,
intent(in) time,
intent(in) timestep,
intent(in) uparam,
intent(inout) uvar,
intent(in) rho,
intent(in) depsxx,
intent(in) depsyy,
intent(in) depsxy,
intent(in) depsyz,
intent(in) depszx,
intent(in) sigoxx,
intent(in) sigoyy,
intent(in) sigoxy,
intent(in) sigoyz,
intent(in) sigozx,
intent(out) signxx,
intent(out) signyy,
intent(out) signxy,
intent(out) signyz,
intent(out) signzx,
intent(inout) pla,
intent(out) dpla,
intent(inout) epsd,
intent(inout) off,
intent(in) gs,
intent(inout) yld,
intent(inout) soundsp,
intent(out) dezz,
integer, intent(in) inloc,
intent(in) dplanl,
integer, intent(in) nvartmp,
integer, dimension(nel,nvartmp), intent(inout) vartmp,
intent(in) loff )

Definition at line 34 of file no_asso_qplas76c.F.

44C-----------------------------------------------
45C M o d u l e s
46C-----------------------------------------------
47 USE table4d_mod
49C-----------------------------------------------
50C I m p l i c i t T y p e s
51C-----------------------------------------------
52#include "implicit_f.inc"
53C-----------------------------------------------
54C G l o b a l P a r a m e t e r s
55C-----------------------------------------------
56#include "tabsiz_c.inc"
57C-----------------------------------------------
58C I N P U T A r g u m e n t s
59C-----------------------------------------------
60 INTEGER, INTENT(IN) :: NEL,NUPARAM,NUVAR,NFUNC,NUMTABL,INLOC
61 INTEGER, INTENT(IN) :: NVARTMP
62 INTEGER, DIMENSION(NFUNC), INTENT(IN) :: IFUNC
63 my_real, INTENT(IN) :: time,timestep
64 my_real, DIMENSION(NUPARAM), INTENT(IN) :: uparam
65 my_real, DIMENSION(NEL), INTENT(IN) :: gs,rho,dplanl,
66 . depsxx,depsyy,depsxy,depsyz,depszx,
67 . sigoxx,sigoyy,sigoxy,sigoyz,sigozx
68 TYPE(TABLE_4D_), DIMENSION(NUMTABL) ,INTENT(IN) :: TABLE
69 my_real, DIMENSION(NEL), INTENT(IN) :: loff
70C-----------------------------------------------
71C O U T P U T A r g u m e n t s
72C-----------------------------------------------
73 my_real, DIMENSION(NEL), INTENT(OUT) ::
74 . signxx,signyy,signxy,signyz,signzx,dezz,dpla
75C-----------------------------------------------
76C I N P U T O U T P U T A r g u m e n t s
77C-----------------------------------------------
78 INTEGER, DIMENSION(NEL,NVARTMP), INTENT(INOUT) :: VARTMP
79 my_real, DIMENSION(NEL,NUVAR) , INTENT(INOUT) :: uvar
80 my_real, DIMENSION(NEL), INTENT(INOUT) :: off,yld,pla,epsd,soundsp
81C-----------------------------------------------
82C VARIABLES FOR FUNCTION INTERPOLATION
83C-----------------------------------------------
84 INTEGER :: NPF(SNPC)
85 my_real :: tf(stf)
86C-----------------------------------------------
87C L o c a l V a r i a b l e s
88C-----------------------------------------------
89 INTEGER I,II,ITER,NITER,ICONV,ICAS,NINDX
90 INTEGER ,PARAMETER :: FUNC_TRAC = 1
91 INTEGER ,PARAMETER :: FUNC_COMP = 2
92 INTEGER ,PARAMETER :: FUNC_SHEAR = 3
93 INTEGER, DIMENSION(NEL) :: INDX,IAD,ILEN
94 my_real :: lam,dlam,epsdot,da0,da1,da2,
95 . ca,cb,aa,bb,cc,a1s,a1c,a1t,a2s,a2c,a2t,e,nu,nu1,nupc,xfac,
96 . yy,dx2,dfdsigdlam,yld_norm,epd_min,epd_max,
97 . normxx,normyy,normxy,normzz,normyz,normzx,alfa,alfi,dtinv,
98 . normgxx,normgyy,normgxy,normgzz,sig_dfdsig,
99 . epdt_min,epdt_max,epdc_min,epdc_max,epds_min,epds_max
100 my_real, DIMENSION(NEL) :: ff,p,svm,svm2,ylds,sxx,syy,sxy,szz,
101 . dpxx,dpyy,dpxy,dpzz,sigt,sigc,sigs,dft,dfc,dfs,a11_2d,a12_2d,g,
102 . a0,a1,a2,nup,epspt,epspc,epsps,epdt,epdc,epds,dydx,
103 . epdt_f,epdc_f,epds_f,dplat,dplac,dplas,gf,alpha,dlam_nl
104 my_real, DIMENSION(NEL,2) :: xvec
105 my_real, DIMENSION(NUMTABL) :: tfac
106 my_real, DIMENSION(NFUNC) :: yfac
107 LOGICAL :: CONV(NEL)
108 my_real, PARAMETER :: sfac = 1.05d0 ! Security factor of ICONV
109c-----------------------------------------------
110c associated plasticity with quadratic yield function
111c-----------------------------------------
112c icas ifunt | ifunc | ifuncs
113c -1 1 | 1 | 1
114c 0 1 | 0 | 0
115c 1 1 | 1 | 0
116c 2 1 | 0 | 1
117c-----------------------------------------
118c UVAR(1) : EPSPT
119c UVAR(2) : EPSPC
120c UVAR(3) : EPSPS
121c UVAR(4) : EPDT
122c UVAR(5) : EPDC
123c UVAR(6) : EPDS
124c UVAR(7) : NUP
125C=======================================================================
126c
127 !=======================================================================
128 ! - INITIALISATION OF COMPUTATION ON TIME STEP
129 !=======================================================================
130 ! Recovering model parameters
131 ! -> Elastic parameters
132 e = uparam(1)
133 nu = uparam(5)
134 ! -> Plastic parameters
135 nupc = uparam(9)
136 ! -> Flags
137 iconv = nint(uparam(15))
138 icas = uparam(17)
139 xfac = uparam(18)
140 alfa = min(one, uparam(16)*timestep)! filtering coefficient for plastic strain rate
141 epdt_min = uparam(19)
142 epdt_max = uparam(20)
143 epdc_min = uparam(21)
144 epdc_max = uparam(22)
145 epds_min = uparam(23)
146 epds_max = uparam(24)
147 tfac(1) = uparam(25)
148 tfac(2) = uparam(26)
149 tfac(3) = uparam(27)
150 yfac(1) = uparam(28)
151 yfac(2) = uparam(29)
152 a11_2d(1:nel) = uparam(2) ! E / (ONE - NU*NU)
153 a12_2d(1:nel) = uparam(3) ! AA2 * NU
154 g(1:nel) = uparam(4)
155 nu1 = nu/(one - nu) ! aa1/aa2
156 alfi = one-alfa
157 dtinv = one / max(em20, timestep)
158c
159 ! Initialize plastic Poisson ratio
160 IF (time == zero) THEN
161 nup(1:nel) = nupc
162 uvar(1:nel,7) = nupc
163 ELSE
164 nup(1:nel) = uvar(1:nel,7)
165 ENDIF
166c
167 ! Recovering internal variables
168 epspt(1:nel) = uvar(1:nel,1)
169 epspc(1:nel) = uvar(1:nel,2)
170 epsps(1:nel) = uvar(1:nel,3)
171 epdt_f(1:nel) = uvar(1:nel,4)
172 epdc_f(1:nel) = uvar(1:nel,5)
173 epds_f(1:nel) = uvar(1:nel,6)
174 dplat(1:nel) = zero
175 dplac(1:nel) = zero
176 dplas(1:nel) = zero
177 DO i=1,nel
178 epdt(i) = min(epdt_max, max(epdt_f(i),epdt_min)) * xfac
179 epdc(i) = min(epdc_max, max(epdc_f(i),epdc_min)) * xfac
180 epds(i) = min(epds_max, max(epds_f(i),epds_min)) * xfac
181 ENDDO
182c
183 ! Computation of yield stresses
184 xvec(1:nel,1) = epspt(1:nel)
185 xvec(1:nel,2) = epdt(1:nel)
186 CALL table_mat_vinterp(table(func_trac),nel,nel,vartmp(1,1),xvec,sigt,dft)
187 sigt(1:nel) = sigt(1:nel) * tfac(1)
188 dft(1:nel) = dft(1:nel) * tfac(1)
189 IF (table(func_comp)%NOTABLE > 0) THEN
190 xvec(1:nel,1) = epspc(1:nel)
191 xvec(1:nel,2) = epdc(1:nel)
192 CALL table_mat_vinterp(table(func_comp),nel,nel,vartmp(1,3),xvec,sigc,dfc)
193 sigc(1:nel) = sigc(1:nel) * tfac(2)
194 dfc(1:nel) = dfc(1:nel) * tfac(2)
195 END IF
196 IF (table(func_shear)%NOTABLE > 0) THEN
197 xvec(1:nel,1) = epsps(1:nel)
198 xvec(1:nel,2) = epds(1:nel)
199 CALL table_mat_vinterp(table(func_shear),nel,nel,vartmp(1,5),xvec,sigs,dfs)
200 sigs(1:nel) = sigs(1:nel) * tfac(3)
201 dfs(1:nel) = dfs(1:nel) * tfac(3)
202 END IF
203 IF (icas == 0) THEN
204 sigc(1:nel) = sigt(1:nel)
205 sigs(1:nel) = sigt(1:nel)/sqr3
206 ELSEIF (icas == 1) THEN
207 DO i=1,nel
208 sigs(i) = sqrt(sigt(i)*sigc(i)/three)
209 END DO
210 ENDIF
211 ! Ensured convexity
212 IF (iconv == 1) THEN ! Ensure convexity
213 DO i = 1,nel
214 conv(i) = .false.
215 IF (sigs(i) < sfac*sqrt(sigc(i)*sigt(i)/three)) THEN
216 sigs(i) = sfac*sqrt(sigc(i)*sigt(i)/three)
217 conv(i) = .true.
218 ENDIF
219 ENDDO
220 ENDIF
221 DO i=1,nel
222 aa = one/sigc(i)/sigt(i)
223 a0(i) = three*sigs(i)**2
224 a1(i) = nine*sigs(i)**2*(sigc(i) - sigt(i))*aa
225 a2(i) = nine*(sigc(i)*sigt(i) - three*sigs(i)**2)*aa
226 END DO
227c
228 !========================================================================
229 ! - COMPUTATION OF TRIAL VALUES
230 !========================================================================
231 DO i=1,nel
232 ! Computation of the trial stress tensor
233 signxx(i) = sigoxx(i) + a11_2d(i)*depsxx(i) + a12_2d(i)*depsyy(i)
234 signyy(i) = sigoyy(i) + a11_2d(i)*depsyy(i) + a12_2d(i)*depsxx(i)
235 signxy(i) = sigoxy(i) + depsxy(i)*g(i)
236 signyz(i) = sigoyz(i) + depsyz(i)*gs(i)
237 signzx(i) = sigozx(i) + depszx(i)*gs(i)
238 p(i) = -(signxx(i) + signyy(i)) * third
239 ! Computation of the deviatoric trial stress tensor
240 sxx(i) = signxx(i) + p(i)
241 syy(i) = signyy(i) + p(i)
242 szz(i) = p(i)
243 dezz(i)= -nu1 * (depsxx(i) + depsyy(i))
244 soundsp(i) = sqrt(a11_2d(i)/rho(i))
245 ENDDO
246c
247 !========================================================================
248 ! - COMPUTATION OF YIELD FONCTION
249 !========================================================================
250 ! Alpha coefficient for non associated function
251 DO i = 1,nel
252 alpha(i) = (nine/two)*((one-two*nup(i))/(one+nup(i)))
253 ENDDO
254 niter = 4
255 nindx = 0
256 DO i=1,nel
257 svm2(i) = three_half*(sxx(i)**2 + syy(i)**2 + szz(i)**2) + three*signxy(i)**2
258 svm(i) = sqrt(svm2(i))
259 ylds(i) = svm2(i) - a0(i) - a1(i)*p(i) - a2(i)*p(i)*p(i)
260 IF (ylds(i) > 0 .AND. off(i) == one) THEN
261 nindx = nindx + 1
262 indx(nindx) = i
263 ENDIF
264 ENDDO
265c
266 !====================================================================
267 ! - PLASTIC RETURN MAPPING WITH CUTTING PLANE METHOD
268 !====================================================================
269 IF (nindx > 0) THEN
270c
271 ! Loop over the iterations
272 DO iter = 1,niter
273c
274 ! Loop over yielding elements
275 DO ii = 1,nindx
276 i = indx(ii)
277c
278 ! Note : in this part, the purpose is to compute for each iteration
279 ! a plastic multiplier allowing to update internal variables to satisfy
280 ! the consistency condition using the cutting plane algorithm
281 ! Its expression at each iteration is : DLAMBDA = - PHI/DPHI_DLAMBDA
282 ! -> PHI : current value of yield function (known)
283 ! -> DPHI_DLAMBDA : derivative of PHI with respect to DLAMBDA by taking
284 ! into account of internal variables kinetic :
285 ! plasticity ...
286c
287 ! function g for non associated flow rule
288 gf(i) = sqrt(max(em20,(svm(i)**2)+alpha(i)*p(i)**2))
289c
290 ! dgf/dsig for non-associated plastic flow
291 normgxx = (three_half*sxx(i) - third*alpha(i)*p(i) ) /gf(i)
292 normgyy = (three_half*syy(i) - third*alpha(i)*p(i) ) /gf(i)
293 normgzz = (three_half*szz(i) - third*alpha(i)*p(i) ) /gf(i)
294 normgxy = three*signxy(i)/gf(i)
295c
296 ! df/dsig
297 cb = a1(i) + two*a2(i)*p(i)
298 normxx = three * sxx(i) + cb /three ! DF/DSIG
299 normyy = three * syy(i) + cb /three
300 normzz = three * szz(i) + cb /three
301 normxy = two *three * signxy(i)
302c
303 ! DF/DSIG * DSIG/DDLAM
304 dfdsigdlam = normxx * (a11_2d(i)*normgxx + a12_2d(i)*normgyy)
305 . + normyy * (a11_2d(i)*normgyy + a12_2d(i)*normgxx)
306 . + normxy * normgxy * g(i)
307c
308 yld_norm = svm(i)/gf(i)
309 bb = three_half/(one + nup(i))
310 dft(i) = dft(i) * yld_norm * bb
311 IF (table(func_comp)%NOTABLE > 0) dfc(i) = dfc(i) * yld_norm * bb
312 IF (table(func_shear)%NOTABLE > 0) dfs(i) = dfs(i) * yld_norm * sqr3/two
313 IF (icas == 0) THEN
314 dfc(i) = dft(i)
315 dfs(i) = (one/sqr3)*dft(i)
316 ELSEIF (icas == 1) THEN
317 dfs(i) = (one/sqr3)*(one/(two*sqrt(sigt(i)*sigc(i))))*
318 . (dfc(i)*sigt(i) + sigc(i)*dft(i))
319 ENDIF
320 IF (iconv == 1) THEN
321 IF (conv(i)) THEN
322 dfs(i) = sfac*(one/sqr3)*(one/(two*sqrt(sigt(i)*sigc(i))))*
323 . (dfc(i)*sigt(i) + sigc(i)*dft(i))
324 ENDIF
325 ENDIF
326c
327 ! derivatives dAi/dlam
328 cc = sigs(i)/sigc(i)/sigt(i)
329C
330 a1s = eighteen*(sigc(i) - sigt(i))*cc
331 a1c = nine*(sigs(i)/sigc(i))**2
332 a1t = -nine*(sigs(i)/sigt(i))**2
333C
334 a2s = -cinquante4*cc
335 a2c = twenty7*cc*sigs(i)/sigc(i)
336 a2t = twenty7*cc*sigs(i)/sigt(i)
337c
338 da0 = six*sigs(i)*dfs(i)
339 da1 = a1s*dfs(i) + a1t*dft(i) + a1c*dfc(i)
340 da2 = a2s*dfs(i) + a2t*dft(i) + a2c*dfc(i)
341C
342 ff(i) = dfdsigdlam + da0 + p(i)*da1 + p(i)**2 * da2
343 ff(i) = sign(max(abs(ff(i)),em20) ,ff(i))
344c
345 dlam = ylds(i)/ff(i)
346 ! Plastic strains tensor update
347 dpxx(i) = dlam * normgxx
348 dpyy(i) = dlam * normgyy
349 dpzz(i) = dlam * normgzz
350 dpxy(i) = dlam * normgxy
351c
352 ! Elasto-plastic stresses update
353 signxx(i) = signxx(i) - (a11_2d(i)*dpxx(i) + a12_2d(i)*dpyy(i))
354 signyy(i) = signyy(i) - (a11_2d(i)*dpyy(i) + a12_2d(i)*dpxx(i))
355 signxy(i) = signxy(i) - dpxy(i)*g(i)
356c
357 ! compute EPSPC(I), EPSPT(I), EPSPS(I)
358 epspt(i) = epspt(i) + dlam* yld_norm * bb
359 epspc(i) = epspt(i)
360 epsps(i) = epsps(i) + dlam* yld_norm*sqr3/two
361c
362 pla(i) = pla(i) + dlam*yld_norm*off(i)
363 dpla(i) = dpla(i) + dlam *yld_norm
364 dplat(i) = dplat(i) + dlam* yld_norm*bb
365 dplac(i) = dplat(i)
366 dplas(i) = dplas(i) + dlam* yld_norm*sqr3/two
367c
368 ENDDO
369c
370 ! Update Yld values and criterion with new plastic strain and strain rate
371 xvec(1:nel,1) = epspt(1:nel)
372 xvec(1:nel,2) = epdt(1:nel)
373 CALL table_mat_vinterp(table(func_trac),nel,nel,vartmp(1,1),xvec,sigt,dft)
374 sigt(1:nel) = sigt(1:nel) * tfac(1)
375 dft(1:nel) = dft(1:nel) * tfac(1)
376 IF (table(func_comp)%NOTABLE > 0) THEN
377 xvec(1:nel,1) = epspc(1:nel)
378 xvec(1:nel,2) = epdc(1:nel)
379 CALL table_mat_vinterp(table(func_comp),nel,nel,vartmp(1,3),xvec,sigc,dfc)
380 sigc(1:nel) = sigc(1:nel) * tfac(2)
381 dfc(1:nel) = dfc(1:nel) * tfac(2)
382 END IF
383 IF (table(func_shear)%NOTABLE > 0) THEN
384 xvec(1:nel,1) = epsps(1:nel)
385 xvec(1:nel,2) = epds(1:nel)
386 CALL table_mat_vinterp(table(func_shear),nel,nel,vartmp(1,5),xvec,sigs,dfs)
387 sigs(1:nel) = sigs(1:nel) * tfac(3)
388 dfs(1:nel) = dfs(1:nel) * tfac(3)
389 END IF
390 IF (icas == 0) THEN
391 DO ii = 1,nindx
392 i = indx(ii)
393 sigc(i) = sigt(i)
394 sigs(i) = sigt(i)/sqr3
395 ENDDO
396 ELSEIF (icas == 1) THEN
397 DO ii = 1,nindx
398 i = indx(ii)
399 sigs(i) = sqrt(sigt(i)*sigc(i)/three)
400 ENDDO
401 ENDIF
402 IF (iconv == 1) THEN
403 DO ii = 1,nindx
404 i = indx(ii)
405 conv(i) = .false.
406 IF (sigs(i) < sfac*sqrt(sigc(i)*sigt(i)/three)) THEN
407 sigs(i) = sfac*sqrt(sigc(i)*sigt(i)/three)
408 conv(i) = .true.
409 ENDIF
410 END DO
411 ENDIF
412 DO ii = 1,nindx
413 i = indx(ii)
414 aa = one/sigc(i)/sigt(i)
415 a0(i) = three*sigs(i)**2
416 a1(i) = nine*sigs(i)**2*(sigc(i) - sigt(i))*aa
417 a2(i) = nine*(sigc(i)*sigt(i) - three*sigs(i)**2)*aa
418 END DO
419c
420 ! update yield function value
421 DO ii = 1,nindx
422 i = indx(ii)
423 p(i) = -third*(signxx(i) + signyy(i) )
424 sxx(i) = signxx(i) + p(i)
425 syy(i) = signyy(i) + p(i)
426 szz(i) = p(i)
427 svm2(i)= three_half*(sxx(i)**2 + syy(i)**2 + szz(i)**2) + three*signxy(i)**2
428 svm(i) = sqrt(svm2(i))
429 ylds(i) = svm2(i) - a0(i) - a1(i)*p(i) - a2(i)*p(i)*p(i)
430 IF (inloc == 0) THEN
431 dezz(i) = dezz(i) + nu1*(dpxx(i) + dpyy(i)) + dpzz(i)
432 ENDIF
433 ENDDO
434 ENDDO ! End Newton iterations
435 ENDIF ! Plasticity
436c
437 !====================================================================
438 ! - UPDATE PLASTIC POISSON RATIO
439 !====================================================================
440 IF (ifunc(1) > 0) THEN
441 iad(1:nel) = npf(ifunc(1)) / 2 + 1
442 ilen(1:nel) = npf(ifunc(1)+1) / 2 - iad(1:nel) - vartmp(1:nel,8)
443!
444 CALL vinter(tf,iad,vartmp(1:nel,8),ilen,nel,pla,dydx,nup)
445!
446 uvar(1:nel,7) = yfac(1) * nup(1:nel)
447 uvar(1:nel,7) = max(zero, min(nup(1:nel), half))
448 END IF
449c
450 !====================================================================
451 ! - NON-LOCAL THICKNESS VARIATION
452 !====================================================================
453 IF (inloc > 0) THEN
454 DO i = 1,nel
455 IF (loff(i) == one) THEN
456 alpha(i)= (nine/two)*((one-two*nup(i))/(one+nup(i)))
457 gf(i) = sqrt(max(em20,(svm(i)**2)+alpha(i)*p(i)**2))
458 normgxx = (three_half*sxx(i) - third*alpha(i)*p(i))/gf(i)
459 normgyy = (three_half*syy(i) - third*alpha(i)*p(i))/gf(i)
460 normgzz = (three_half*szz(i) - third*alpha(i)*p(i))/gf(i)
461 yld_norm = svm(i)/gf(i)
462 IF (yld_norm /= zero) THEN
463 dlam_nl(i) = (one/yld_norm)*max(dplanl(i),zero)
464 dezz(i) = dezz(i) + nu1*(dlam_nl(i)*normgxx)
465 . + nu1*(dlam_nl(i)*normgyy)
466 . + dlam_nl(i)*normgzz
467 ENDIF
468 ENDIF
469 ENDDO
470 ENDIF
471c
472 !====================================================================
473 ! - STORING NEW VALUES
474 !====================================================================
475 uvar(1:nel,1) = epspt(1:nel)
476 uvar(1:nel,2) = epspc(1:nel)
477 uvar(1:nel,3) = epsps(1:nel)
478 uvar(1:nel,4) = alfa*dplat(1:nel)*dtinv + alfi*epdt_f(1:nel)
479 uvar(1:nel,5) = alfa*dplac(1:nel)*dtinv + alfi*epdc_f(1:nel)
480 uvar(1:nel,6) = alfa*dplas(1:nel)*dtinv + alfi*epds_f(1:nel)
481 epsd(1:nel) = alfa*dpla(1:nel)*dtinv + alfi*epsd(1:nel) ! save for strain rate output
482c
#define my_real
Definition cppsort.cpp:32
#define alpha
Definition eval.h:35
subroutine func_comp(table, ntable, nptmax, npt_trac, npt_shear, npt_comp, x_comp, y_comp, nup)
Definition law76_upd.F:344
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
subroutine table_mat_vinterp(table, dimx, nel, ipos, xx, yy, dydx)
subroutine vinter(tf, iad, ipos, ilen, nel, x, dydx, y)
Definition vinter.F:72