基于ABAQUS的VDLOAD子程序实现四轮小车移动

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。

下面是运行一段时间的结果。