黄永刚单晶塑性有限元umat子程序

C IWKDIR -- all possible slip directions (integer) C IWKNOR -- all possible slip planes (integer)

C----- Use single precision on cray C

IMPLICIT REAL*8 (A-H,O-Z)

DIMENSION ISPDIR(3), ISPNOR(3), SLPDIR(3,50), SLPNOR(3,50), * ROTATE(3,3), IWKDIR(3,24), IWKNOR(3,24), TERM(3)

NSLIP=0 NSPDIR=0 NSPNOR=0

C----- Generating all possible slip directions in this set C

C Denote the slip direction by [lmn]. I1 is the minimum of the C absolute value of l, m and n, I3 is the maximum and I2 is the C mode, e.g. (1 -3 2), I1=1, I2=2 and I3=3. I1<=I2<=I3.

I1=MIN(IABS(ISPDIR(1)),IABS(ISPDIR(2)),IABS(ISPDIR(3))) I3=MAX(IABS(ISPDIR(1)),IABS(ISPDIR(2)),IABS(ISPDIR(3))) I2=IABS(ISPDIR(1))+IABS(ISPDIR(2))+IABS(ISPDIR(3))-I1-I3

RMODIR=SQRT(FLOAT(I1*I1+I2*I2+I3*I3))

C I1=I2=I3=0

IF (I3.EQ.0) THEN

WRITE (6,*) '***ERROR - slip direction is [000]' STOP

C I1=I2=0, I3>0 --- [001] type ELSE IF (I2.EQ.0) THEN NSPDIR=3 DO J=1,3 DO I=1,3

IWKDIR(I,J)=0

IF (I.EQ.J) IWKDIR(I,J)=I3 END DO END DO

C I1=0, I3>=I2>0

ELSE IF (I1.EQ.0) THEN

C I1=0, I3=I2>0 --- [011] type IF (I2.EQ.I3) THEN NSPDIR=6 DO J=1,6 DO I=1,3

IWKDIR(I,J)=I2

IF (I.EQ.J.OR.J-I.EQ.3) IWKDIR(I,J)=0 IWKDIR(1,6)=-I2 IWKDIR(2,4)=-I2 IWKDIR(3,5)=-I2 END DO END DO

C I1=0, I3>I2>0 --- [012] type ELSE

NSPDIR=12

CALL LINE1 (I2, I3, IWKDIR(1,1), 1) CALL LINE1 (I3, I2, IWKDIR(1,3), 1) CALL LINE1 (I2, I3, IWKDIR(1,5), 2) CALL LINE1 (I3, I2, IWKDIR(1,7), 2) CALL LINE1 (I2, I3, IWKDIR(1,9), 3) CALL LINE1 (I3, I2, IWKDIR(1,11), 3)

END IF

C I1=I2=I3>0 --- [111] type ELSE IF (I1.EQ.I3) THEN NSPDIR=4

CALL LINE (I1, I1, I1, IWKDIR)

C I3>I2=I1>0 --- [112] type ELSE IF (I1.EQ.I2) THEN NSPDIR=12

CALL LINE (I1, I1, I3, IWKDIR(1,1)) CALL LINE (I1, I3, I1, IWKDIR(1,5)) CALL LINE (I3, I1, I1, IWKDIR(1,9))

C I3=I2>I1>0 --- [122] type ELSE IF (I2.EQ.I3) THEN NSPDIR=12

CALL LINE (I1, I2, I2, IWKDIR(1,1)) CALL LINE (I2, I1, I2, IWKDIR(1,5)) CALL LINE (I2, I2, I1, IWKDIR(1,9))

C I3>I2>I1>0 --- [123] type ELSE

NSPDIR=24

CALL LINE (I1, I2, I3, IWKDIR(1,1)) CALL LINE (I3, I1, I2, IWKDIR(1,5)) CALL LINE (I2, I3, I1, IWKDIR(1,9)) CALL LINE (I1, I3, I2, IWKDIR(1,13)) CALL LINE (I2, I1, I3, IWKDIR(1,17)) CALL LINE (I3, I2, I1, IWKDIR(1,21))

END IF

C----- Generating all possible slip planes in this set C

C Denote the normal to slip plane by (pqr). J1 is the minimum of C the absolute value of p, q and r, J3 is the maximum and J2 is the C mode, e.g. (1 -2 1), J1=1, J2=1 and J3=2. J1<=J2<=J3.

J1=MIN(IABS(ISPNOR(1)),IABS(ISPNOR(2)),IABS(ISPNOR(3))) J3=MAX(IABS(ISPNOR(1)),IABS(ISPNOR(2)),IABS(ISPNOR(3))) J2=IABS(ISPNOR(1))+IABS(ISPNOR(2))+IABS(ISPNOR(3))-J1-J3

RMONOR=SQRT(FLOAT(J1*J1+J2*J2+J3*J3))

IF (J3.EQ.0) THEN

WRITE (6,*) '***ERROR - slip plane is [000]' STOP

C (001) type

ELSE IF (J2.EQ.0) THEN NSPNOR=3 DO J=1,3 DO I=1,3

IWKNOR(I,J)=0

IF (I.EQ.J) IWKNOR(I,J)=J3 END DO END DO

ELSE IF (J1.EQ.0) THEN

C (011) type

IF (J2.EQ.J3) THEN NSPNOR=6 DO J=1,6

DO I=1,3

IWKNOR(I,J)=J2

IF (I.EQ.J.OR.J-I.EQ.3) IWKNOR(I,J)=0 IWKNOR(1,6)=-J2 IWKNOR(2,4)=-J2 IWKNOR(3,5)=-J2 END DO END DO

C (012) type ELSE

NSPNOR=12

CALL LINE1 (J2, J3, IWKNOR(1,1), 1) CALL LINE1 (J3, J2, IWKNOR(1,3), 1) CALL LINE1 (J2, J3, IWKNOR(1,5), 2) CALL LINE1 (J3, J2, IWKNOR(1,7), 2) CALL LINE1 (J2, J3, IWKNOR(1,9), 3) CALL LINE1 (J3, J2, IWKNOR(1,11), 3)

END IF

C (111) type

ELSE IF (J1.EQ.J3) THEN NSPNOR=4

CALL LINE (J1, J1, J1, IWKNOR)

C (112) type

ELSE IF (J1.EQ.J2) THEN NSPNOR=12

CALL LINE (J1, J1, J3, IWKNOR(1,1)) CALL LINE (J1, J3, J1, IWKNOR(1,5)) CALL LINE (J3, J1, J1, IWKNOR(1,9))

C (122) type

ELSE IF (J2.EQ.J3) THEN NSPNOR=12

CALL LINE (J1, J2, J2, IWKNOR(1,1)) CALL LINE (J2, J1, J2, IWKNOR(1,5)) CALL LINE (J2, J2, J1, IWKNOR(1,9))

C (123) type ELSE

NSPNOR=24

CALL LINE (J1, J2, J3, IWKNOR(1,1))