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 )Cinclude 'vaba_param.inc'Cdimension curCoords(nBlock,ndim), velocity(nBlock,ndim),1 dirCos(nBlock,ndim,ndim), value(nBlock)snameCdo 100 km = 1, nBlockuser coding to define value100 continuereturnend
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 = 6L2 = 3a = 0.5!运动参数speed = 10000t = steptimedistance = speed*steptime!初始时刻车身距路面的距离xm = 3ym = 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) thenvalue(k) = pressure!第二个轮子else if (x >= xm+distance+a+L1 .and. x <= xm+distance+a*2+L1) thenvalue(k) = pressureelsevalue(k) = 0.0end if!下面一排车轮elif (y >= ym+a+L2 .and. y <= ym+a*2+L2) then!第一个轮子if (x >= xm+distance .and. x <= xm+distance+a) thenvalue(k) = pressure!第二个轮子else if (x >= xm+distance+a+L1 .and. x <= xm+distance+a*2+L1) thenvalue(k) = pressureelsevalue(k) = 0.0end ifelsevalue(k) = 0.0end ifend 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 )Cinclude 'vaba_param.inc'Cdimension curCoords(nBlock,ndim), velocity(nBlock,ndim),1 dirCos(nBlock,ndim,ndim), value(nBlock)snameC!车身参数L1 = 6L2 = 3a = 0.5!运动参数speed = 10000t = steptimedistance = speed*steptime!初始时刻车身距路面的距离xm = 3ym = 3pressure = 3000do 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= pressure!第二个轮子else if (x >= xm+distance+a+L1 .and. x <= xm+distance+a*2+L1) then= pressureelse= 0.0end if!下面一排车轮else if (y >= ym+a+L2 .and. y <= ym+a*2+L2) then!第一个轮子if (x >= xm+distance .and. x <= xm+distance+a) then= pressure!第二个轮子else if (x >= xm+distance+a+L1 .and. x <= xm+distance+a*2+L1) then= pressureelse= 0.0end ifelse= 0.0end ifend doreturnend
为了能够观察到塑性应变,对载荷进行更改为3000MPa。
下面是运行一段时间的结果。
