UMAT中是否一定要给出雅克比矩阵？

漂泊星6

UMAT中是否一定要给出雅克比矩阵？谢谢大侠们指导！

changbing

有助于收敛，你可以自己试试，如果不写雅克比矩阵也能收敛，那就可以不用加。雅克比矩阵只影响计算的收敛性，不影响最终的计算结果。只要能收敛，加不加雅克比矩阵算出来的结果都一样。但是如果稍微复杂的一点的本构方程，不加入雅克比矩阵，基本都计算不收敛。
不过需要提醒的是，有时候本构关系相对较为复杂，雅克比矩阵可能很难推导，给你一个小技巧，那就是直接把刚度矩阵作为雅克比矩阵，很多时候也能够达到收敛。因为雅克比矩阵只是一个为了加快收敛的作用，所以即使你的雅克比矩阵不完全是按照它定义的方法得到的，但是差别只要不是很大，也同样能够达到收敛的目的，只是收敛速度会慢一些。

漂泊星6

changbing 发表于 2016-7-16 12:58
有助于收敛，你可以自己试试，如果不写雅克比矩阵也能收敛，那就可以不用加。雅克比矩阵只影响计算的收敛性 ...

非常感谢，我想写一个拉压不同模量的弹性本钩子程序，拉压条件下雅克比矩阵是不同。所以就不知道怎么给了

漂泊星6

changbing 发表于 2016-7-16 12:58
有助于收敛，你可以自己试试，如果不写雅克比矩阵也能收敛，那就可以不用加。雅克比矩阵只影响计算的收敛性 ...

您好：下面是我自己写的双模量材料本构，计算过程中没有用到雅克比矩阵，在UMAT程序结束后给了雅克比矩阵，提交后出现如下错误且计算并没真正开始：The executable d:\SIMULIA\Abaqus\6.10-1\exec\standard.exe aborted with system error code 5. Please check the .dat, .msg, and .sta files for error messages if the files exist.  If there are no error messages and you cannot resolve the problem, please run the command "abaqus job=support information=support" to report and save your system information.  Use the same command to run Abaqus that you used when the problem occurred.  Please contact your local Abaqus support office and send them the input file, the file support.log which you just created, the executable name, and the error code.将雅克比矩阵删除后仍然出现这个错误，在相关文件中也没找到错误。
UMAT程序：
      SUBROUTINE UMAT(STRESS,STATEV,DDSDDE,SSE,SPD,SCD,
     1     RPL,DDSDDT,DRPLDE,DRPLDT,
     2     STRAN,DSTRAN,TIME,DTIME,TEMP,DTEMP,PREDEF,DPRED,CMNAME,
     3     NDI,NSHR,NTENS,NSTATV,PROPS,NPROPS,COORDS,DROT,PNEWDT,
     4     CELENT,DFGRD0,DFGRD1,NOEL,NPT,LAYER,KSPT,KSTEP,KINC)
C     
      INCLUDE 'ABA_PARAM.INC'
C     
      CHARACTER*80 CMNAME
      DIMENSION STRESS(NTENS),STATEV(NSTATV),
     1     DDSDDE(NTENS,NTENS),
     2     DDSDDT(NTENS),DRPLDE(NTENS),
     3     STRAN(NTENS),DSTRAN(NTENS),TIME(2),PREDEF(1),DPRED(1),
     4     PROPS(NPROPS)
C --------------------------------------------------------------------
C     UMAT FOR ISOTROPIC WITH DIFFERENT MODULUS
C     CAN ONLY USED FOR PLANE STRESS
C--------------------------------------------------------------------
C     PROPS(1) - E1 TENSION
C     PROPS(2) - NU1 TENSION POISSON   
C     PROPS(3) - E2 COMPRESS
C     PROPS(4) - NU2 COMPRESS POISSON
C---------------------------------------------------------------------
C     LOCAL ARRAYS
C---------------------------------------------------------------------
C     EELASP - PRINCIPAL ELASTIC STRAINS
C      
      DIMENSION STRANT(NTENS),STRESSP(NDI),QTEN(NDI,NDI),
     1 QCOM(NDI,NDI),QTENC(NDI,NDI),DSTFF(NDI,NDI)
      DIMENSION STRANP(NDI),STRAND(NDI,NDI)
      DIMENSION L(2,3),LT(3,2),DD(3,2)
      PARAMETER (ZERO = 0.0,ONE = 1.0,TWO = 2.0, HALF = 0.50)

      IF (NDI.NE.2) THEN
         WRITE (7,*) 'THIS UMAT MAY ONLY USED FOR PLANE STRESS'
         CALL XIT
      ENDIF
C
C     ELASTIC PROPERTIES
      EMODT=PROPS(1)
      ENUT=PROPS(2)
      EMODP=PROPS(3)
      ENUP=PROPS(4)
C
C     CALCULATE STIFFNESS COEFFICIENT
      DO K1=1, NDI
        DO K2=1, NDI
            QTEN(K2,K1)=EMODT*ENUT/(ONE-ENUT**2)  ! 拉
            QCOM(K2,K1)=EMODP*ENUP/(ONE-ENUP**2)  ! 压
          END DO
            QTEN(K1,K1)=EMODT/(ONE-ENUT**2)  ! 拉
            QCOM(K1,K1)=EMODP/(ONE-ENUP**2)  ! 压
      END DO
            QTENC(1,1)=EMODT/(ONE-ENUT*ENUP)  ! 拉-压
            QTENC(2,2)=EMODP/(ONE-ENUT*ENUP)  ! 拉-压
            QTENC(1,2)=EMODP*ENUT/(ONE-ENUT*ENUP)  ! 拉-压
            QTENC(2,1)=EMODT*ENUP/(ONE-ENUT*ENUP)  ! 拉-压               
C     CALCULATE PRINCIPLE STRAIN PRINCIPLE DIRECTION
      DO K1=1,NTENS   ! 计算总应变张量
         STRANT=STRAN(K1)+DSTRAN(K1)
      END DO
C      ！计算主应变、主方向  
      CALL SPRIND ( STRANT,STRANP,STRAND, 2, NDI, NSHR)
C
C     CALCULATE PRINCIPLE STRESS & 雅克比矩阵
C     DETERMINE STRESS STATUS AND GIVE STIFFNESS
        DO K1=1,NDI
          DO K2=1,NDI
             DSTFF(K1,K2)=ZERO
          END DO
        END DO
C
      IF (STRANP(1).LE.ZERO) THEN
        DO K1=1,NDI
          DO K2=1,NDI
             DSTFF(K1,K2)=QCOM(K1,K2)
          END DO
        END DO              
      else if (STRANP(2).GE.ZERO) THEN
         DO K1=1,NDI
           DO K2=1,NDI
             DSTFF(K1,K2)=QTEN(K1,K2)
           END DO
         END DO      
      ELSE
        DO K1=1,NDI
          DO K2=1,NDI
             DSTFF(K1,K2)=QTENC(K1,K2)
          END DO
        END DO         
      END IF
C     CALCULATE PRINCIPLE STRESS
C
      DO K1=1,NDI
         STRESSP(K1)=0
        DO K2=1,NDI
         STRESSP(K2)=STRESSP(K2)+DSTFF(K2,K1)*STRANP(K1)
        END DO
      END DO
C   
C     CALCULATE  STRESS           
      STRESS(1)=STRESSP(1)*STRAND(1,1)**2+STRESSP(2)*STRAND(2,1)**2
      STRESS(2)=STRESSP(1)*STRAND(2,1)**2+STRESSP(2)*STRAND(2,1)**2     
      STRESS(3)=(STRESSP(1)-STRESSP(2))*STRAND(1,1)*STRAND(2,1)      
C     CALCULATE  JACABIN MATRIX  
C
      L(1,1)=STRAND(1,1)**2
      L(1,2)=STRAND(2,1)**2
      L(1,3)=STRAND(1,1)*STRAND(2,1)     
      L(2,1)=STRAND(2,1)**2
      L(2,2)=STRAND(1,1)**2
      L(2,3)=-STRAND(1,1)*STRAND(2,1)
C
      LT(1,1)=STRAND(1,1)**2
      LT(2,1)=STRAND(2,1)**2
      LT(3,1)=STRAND(1,1)*STRAND(2,1)     
      LT(1,2)=STRAND(2,1)**2
      LT(2,2)=STRAND(1,1)**2
      LT(3,2)=-STRAND(1,1)*STRAND(2,1)   
C
      DD=MATMUL(LT,DSTFF)
      DDSDDE=MATMUL(DD,L)   
      RETURN
      END
请大侠指导，不胜感激！

qqcyhappy

这个是不是因为方向余弦的问题呀？

小仙女的尾巴

我也是编写的混凝土本构，可以收敛一点点，但是计算得到的应力全部都是零，挠度和支反力倒是有，我也不知道为啥

小仙女的尾巴

我建议你去查一下code 5到底是什么错误，我倒是没这个错误，就是计算一点就不收敛，且应力皆为零。

liufenjie_7

我最近在编写超弹性本构模型的UMAT程序遇到一个问题，总说我的单元没有赋予弹性特性，一直运行不了，有没有哪位前辈遇到过这种情况，帮忙指点一下

newmusic

如果能够直接写出应力与应变（随时间变化）的关系式，还要写雅克比矩阵吗？