umat子程序编写常用的fortran函数分享（四）
四阶张量双点乘：
subroutine aaaa_dot_dot_bbbb(n,a,b,c)
include 'ABA_PARAM.INC'
dimension a(n,n,n,n), b(n,n,n,n), c(n,n,n,n)
do i = 1,n
do j = 1,n
do k = 1,n
do l = 1,n
c(i,j,k,l) = 0
do m1 = 1,n
do m2 = 1,n
c(i,j,k,l) = c(i,j,k,l) + a(i,j,m1,m2) * b(m1,m2,k,l)
end do !m2
end do !m1
end do !l
end do !k
end do !j
end do !i
return
end
四阶张量与二阶张量双点乘：
subroutine aaaa_dot_dot_bb(n,a,b,c)
include 'ABA_PARAM.INC'
dimension a(n,n,n,n), b(n,n), c(n,n)
do i = 1,n
do j = 1,n
c(i,j) = 0
do k = 1,n
do l = 1,n
c(i,j) = c(i,j) + a(i,j,k,l) * b(k,l)
end do !l
end do !k
end do !j
end do !i
return
end
二阶张量与四阶张量双点乘：
subroutine aa_dot_dot_bbbb(n,a,b,c)
include 'ABA_PARAM.INC'
dimension a(n,n), b(n,n,n,n), c(n,n)
do i = 1,n
do j = 1,n
c(i,j) = 0
do k = 1,n
do l = 1,n
c(i,j) = c(i,j) + a(k,l) * b(k,l,i,j)
end do !l
end do !k
end do !j
end do !i
return
end
对称的四阶张量Voigt表示形式：
subroutine Voigt_to_forth(b,a)
include 'ABA_PARAM.INC'
dimension a(3,3,3,3), b(6,6)
do i = 1,3
do j = 1,3
ia = i
if (i.ne.j) ia=9-i-j
do k = 1,3
do l = 1,3
ib = k
if (k.ne.l) ib=9-k-l
a(i,j,k,l) = b(ia,ib)
if (ia.gt.3) a(i,j,k,l) = a(i,j,k,l) / 2
if (ib.gt.3) a(i,j,k,l) = a(i,j,k,l) / 2
end do
end do
end do
end do
return
end
sgn函数：
function sgn(a)
include 'ABA_PARAM.INC'
sgn = 1.0
if (a .lt. 0.0) sgn = -1.0
if (a .eq. 0.0) sgn = 0.0
return
end
LU分解：
SUBROUTINE LUDCMP (A, N, NP, INDX, D)
IMPLICIT REAL*8 (A-H,O-Z)
PARAMETER (NMAX=200, TINY=1.0E-20)
DIMENSION A(NP,NP), INDX(N), VV(NMAX)
D=1.
DO I=1,N
AAMAX=0.
DO J=1,N
IF (ABS(A(I,J)).GT.AAMAX) AAMAX=ABS(A(I,J))
END DO
IF (AAMAX.EQ.0.) PAUSE 'Singular matrix.'
VV(I)=1./AAMAX
END DO
C
DO J=1,N
DO I=1,J-1
SUM=A(I,J)
DO K=1,I-1
SUM=SUM-A(I,K)*A(K,J)
END DO
A(I,J)=SUM
END DO
AAMAX=0.
DO I=J,N
SUM=A(I,J)
DO K=1,J-1
SUM=SUM-A(I,K)*A(K,J)
END DO
A(I,J)=SUM
DUM=VV(I)*ABS(SUM)
IF (DUM.GE.AAMAX) THEN
IMAX=I
AAMAX=DUM
END IF
END DO
IF (J.NE.IMAX) THEN
DO K=1,N
DUM=A(IMAX,K)
A(IMAX,K)=A(J,K)
A(J,K)=DUM
END DO
D=-D
VV(IMAX)=VV(J)
END IF
INDX(J)=IMAX
IF (A(J,J).EQ.0.) A(J,J)=TINY
IF (J.NE.N) THEN
DUM=1./A(J,J)
DO I=J+1,N
A(I,J)=A(I,J)*DUM
END DO
END IF
END DO
RETURN
END
SUBROUTINE LUBKSB (A, N, NP, INDX, B)
IMPLICIT REAL*8 (A-H,O-Z)
DIMENSION A(NP,NP), INDX(N), B(N)
II=0
DO I=1,N
LL=INDX(I)
SUM=B(LL)
B(LL)=B(I)
IF (II.NE.0) THEN
DO J=II,I-1
SUM=SUM-A(I,J)*B(J)
END DO
ELSE IF (SUM.NE.0.) THEN
II=I
END IF
B(I)=SUM
END DO
DO I=N,1,-1
SUM=B(I)
IF (I.LT.N) THEN
DO J=I+1,N
SUM=SUM-A(I,J)*B(J)
END DO
END IF
B(I)=SUM/A(I,I)
END DO
RETURN
END
