ABAQUS因其具有強大的子程序二次開發功能而倍受研究人員歡迎。VDLOAD子程序可以定義隨時間和空間變化的載荷分布，可用于模擬汽車移動和激光沖擊過程。

其子程序模板如下：

      subroutine vdload (
C Read only (unmodifiable)variables -
     1 nBlock, ndim, stepTime, totalTime,
     2 amplitude, curCoords, velocity, dirCos, jltyp, sname,
C Write only (modifiable) variable -
     1 value )
C
      include 'vaba_param.inc'
C
      dimension curCoords(nBlock,ndim), velocity(nBlock,ndim),
     1  dirCos(nBlock,ndim,ndim), value(nBlock)
      character*80 sname
C


      do 100 km = 1, nBlock
        user coding to define value


  100 continue


      return
      end

value(nblock)是用戶自定義的載荷大小，是程序中最重要的參數。

ndim:坐標的方向數，1、2、和3分別表示三維空間中、全局坐標下的X、Y和Z三個坐標分量。

steptime:表示當前分析步時間值，若有多個分析步存在，則其變化為0—>step1time—>0—>step2time—>…….—>0—>stepntime

totaltime:表示當前時間值，從0開始逐漸增大。

curCoords(nblock, ndim):當前節點坐標。

以上是用戶常用的參數及其含義。

本貼就以VDLOAD子程序實現生活中常見的四輪小車移動，程序相對簡單，主要難點在于車輪的坐標控制，代碼開源。

對于四輪小車而言，假設車輪與地面的接觸區域為正方形，示意圖如下。

其中前后車輪間距為L1，并排車輪距離為L2。藍色區域為加載區，加載區長度為a*a。在ABAQUS中建模時，路面為XY平面，坐標原點位于路面模型的頂點。本模型默認整車位于路面的正中央，路面的長和寬分別為30mm和10mm，路面長的方向為X方向，寬為Y方向。初始時刻車輪距路面邊緣的距離分別為xm和ym。

先定義車身參數和運動參數。

!車身參數
L1 = 6
L2 = 3
a = 0.5
!運動參數
speed = 10000
t = steptime
distance = speed*steptime
!初始時刻車身距路面的距離
xm = 3
ym = 3

每個車輪對路面的壓強為：

pressure = 2000

主程序：

do k=1, nblock
  !節點坐標
  x = curCoords(k,1)
  y = curCoords(k,2)
  !上面一排車輪施加載荷
  if (y >= ym .and. y <= ym+a)then
    !第一個輪子
    if (x >= xm+distance .and. x <= xm+distance+a) then
      value(k) = pressure
    !第二個輪子
    else if (x >= xm+distance+a+L1 .and. x <= xm+distance+a*2+L1) then
      value(k) = pressure
    else
      value(k) = 0.0
    end if
  !下面一排車輪
  elif (y >= ym+a+L2 .and. y <= ym+a*2+L2) then
    !第一個輪子
    if (x >= xm+distance .and. x <= xm+distance+a) then
      value(k) = pressure
    !第二個輪子
    else if (x >= xm+distance+a+L1 .and. x <= xm+distance+a*2+L1) then
      value(k) = pressure
    else
      value(k) = 0.0
    end if
  else
    value(k) = 0.0
  end if
end do

程序編寫完畢，接下來就是在ABAQUS軟件建模，模型相對簡單，就是一個板材。構建的模型如下。

完整的代碼如下

      subroutine vdload (
C Read only (unmodifiable)variables -
     1 nBlock, ndim, stepTime, totalTime,
     2 amplitude, curCoords, velocity, dirCos, jltyp, sname,
C Write only (modifiable) variable -
     1 value )
C
      include 'vaba_param.inc'
C
      dimension curCoords(nBlock,ndim), velocity(nBlock,ndim),
     1  dirCos(nBlock,ndim,ndim), value(nBlock)
      character*80 sname
C


      !車身參數
      L1 = 6
      L2 = 3
      a = 0.5
      !運動參數
      speed = 10000
      t = steptime
      distance = speed*steptime
      !初始時刻車身距路面的距離
      xm = 3
      ym = 3
      pressure = 3000
      
      do k=1, nblock
      !節點坐標
      x = curCoords(k,1)
      y = curCoords(k,2)
      !上面一排車輪施加載荷
      if (y >= ym .and. y <= ym+a)then
        !第一個輪子
        if (x >= xm+distance .and. x <= xm+distance+a) then
          value(k) = pressure
        !第二個輪子
        else if (x >= xm+distance+a+L1 .and. x <= xm+distance+a*2+L1) then
          value(k) = pressure
        else
          value(k) = 0.0
        end if
      !下面一排車輪
      else if (y >= ym+a+L2 .and. y <= ym+a*2+L2) then
        !第一個輪子
        if (x >= xm+distance .and. x <= xm+distance+a) then
          value(k) = pressure
        !第二個輪子
        else if (x >= xm+distance+a+L1 .and. x <= xm+distance+a*2+L1) then
          value(k) = pressure
        else
          value(k) = 0.0
        end if
      else
        value(k) = 0.0
      end if
      
      end do


      return
      end

為了能夠觀察到塑性應變，對載荷進行更改為3000MPa。

下面是運行一段時間的結果。