工作环境是ubuntu20.04，python3.9.12，taichi1.1.3，cuda11.3，llvm 10.0.0，commit 1262a70a
我是一个刚接触计算机编程半年左右的初学者，最近在使用taichi构建一个布料对折的物理模拟，然后通过点的位置误差，向后传播来更新其他参数（比如弹簧质点模型里弹簧的K值）。
运行的时候会出现如下错误提示
发生异常: RuntimeError
[jit_cuda.h:lookup_function@58] Cannot look up function simulation_c79_0_reverse_grad_reverse_grad_kernel_0_range_for
File “/home/ais/difftaichi/examples/test14.py”, line 166, in main
forward()
File “/home/ais/difftaichi/examples/test14.py”, line 183, in
main()
完整源代码:test14文件
代码大致思路
定义需要grad的变量:
初始化布料
布料物理模拟
计算误差、清除grade数据
主程序
大致代码内容
import matplotlib.pyplot as plt
ti.init(arch=ti.cuda,device_memory_fraction=0.3)
scalar = lambda: ti.field(dtype=ti.f32)
vec = lambda: ti.Vector.field(3, dtype=ti.f32)
pos=vec()
vel=vec()
F=vec()
acc=vec()
KStruct=scalar()
KShear=scalar()
KBend=scalar()
loss_n=scalar()
mass=scalar()
damping=scalar()
ti.root.dense(ti.ij,(ClothResX+1,ClothResX+1)).place(pos,vel,F,acc)
ti.root.dense(ti.i,1).place(KStruct,KShear,KBend,mass,damping)
ti.root.place(loss_n)
ti.root.lazy_grad()
@ti.func
def Compute_Force(coord):
p1=pos[coord]
F[coord]=gravity*mass[0]-vel[coord]*damping[0]
for n in range(0,12):
#计算弹簧力。
Sping_Type=ti.Vector([Get_X(n),Get_Y(n)])
#struct，shear，bend，各四个弹簧，总共十二个弹簧。
Coord_Neigh=coord+Sping_Type
#弹簧的另一端点坐标
if (Coord_Neigh.x>=0) and (Coord_Neigh.x<=ClothResX) and (Coord_Neigh.y>=0) and (Coord_Neigh.y<=ClothResX):
#保证弹簧在布料范围内
Sping_Vector=Sping_Type*ti.Vector([ClothWid/ClothResX,ClothHgt/ClothResX])
Rest_Length=ti.sqrt(Sping_Vector.x**2+Sping_Vector.y**2)
#计算弹簧原始长度
p2=pos[Coord_Neigh]
deltaP=p1-p2
dist=ti.sqrt(deltaP.x**2+deltaP.y**2+deltaP.z**2)
#计算弹簧现在的长度
Sping_Force=ti.Vector([0.0,0.0,0.0])
if (n<4):
Sping_Force=(KStruct[0]*(dist-Rest_Length)*((deltaP)/dist))
elif (n>=4 and n<8):
Sping_Force=(KShear[0]*(dist-Rest_Length)*((deltaP)/dist))
else:
Sping_Force=(KBend[0]*(dist-Rest_Length)*((deltaP)/dist))
#胡克定律计算弹簧力
F[coord]+=Sping_Force
@ti.func
def collision(coord):
#将y=0的平面定为桌面，检查是否碰撞
if(pos[coord].y<0):
pos[coord].y=0
@ti.kernel
def Reset_Cloth():
#初始化布料模型，布料正中心点为[0.0,0.0,0.0]
for i,j in pos:
pos[i,j] = ti.Vector([ClothWid*(i/ClothResX)-ClothWid/2.,0.0,ClothHgt*(j/ClothResX)-ClothHgt/2.0])
vel[i,j]=ti.Vector([0.0,0.0,0.0])
F[i,j]=ti.Vector([0.0,0.0,0.0])
acc[i,j]=ti.Vector([0.0,0.0,0.0])
@ti.kernel
def simulation(t:ti.i32):
for i,j in pos:
coord=ti.Vector([i,j])
Compute_Force(coord)
if (i == 0 and j == 0 ) or (i == 35 and j == 0 ):
#抓取布料两个顶点，改变它们的坐标来折叠布料，坐标点为（-2.0，0.0,-2.0)和(2.0,0.0,-2.0)
z=(ClothWid/step)*t
y=a*z**2+b
x=0.0
if pos[coord].x <0.0:
x=-ClothHgt/2.0
else:
x=ClothHgt/2.0
pos[coord]=[x,y,z]
#顶点移动路径为f(z)=-0.3z**2+1.2
else:
#除布料顶点外其他点更新加速度、速度和位置。
acc[coord]=F[coord]/mass[0]
vel[coord]=acc[coord]*dt
pos[coord]+=vel[coord]*dt
collision(coord)
#将y=0的平面定为桌面，检查是否碰撞
@ti.kernel
def Compute_Loss():
current=pos[18,0]
target=ti.Vector([0.0,0.0,2.0])
Loss_Vector=current-target
loss_n[None]=0.5*ti.sqrt(Loss_Vector.x**2+Loss_Vector.y**2+Loss_Vector.z**2)
#将两个抓取点所在的边的正中间的点定为误差计算对象，布料初始化时，位置为[0.0,0.0,-2.0]对折后目标位置为[0.0,0.0,2.0]
@ti.kernel
def Vec_Clear():
#初始化所有向量的grad
for i,j in pos:
pos.grad[i,j]=ti.Vector([0.0,0.0,0.0])
vel.grad[i,j]=ti.Vector([0.0,0.0,0.0])
F.grad[i,j]=ti.Vector([0.0,0.0,0.0])
acc.grad[i,j]=ti.Vector([0.0,0.0,0.0])
@ti.kernel
def Scalar_Clear():
#初始化所有标量的grad
KStruct.grad[0]=0.0
KShear.grad[0]=0.0
KBend.grad[0]=0.0
loss_n.grad[None]=0.0
mass.grad[0]=0.0
damping.grad[0]=0.0
def Grad_Clear():
#初始化所有可微变量
Vec_Clear()
Scalar_Clear()
def forward():
for t in range(step):
simulation(t)
Compute_Loss()
def main():
Reset_Cloth()
forward()
Grad_Clear()
for i in range(step):
with ti.ad.Tape(loss=loss_n,clear_gradients=True):
forward(）
KStruct[0]-=learning_rate * KStruct.grad[0]
KShear[0]-=learning_rate * KShear.grad[0]
KBend[0]-=learning_rate * KBend.grad[0]
Reset_Cloth()