本发明属于空间操作机构运动学领域，尤其涉及一种空间可变形桁架逆运动学计算方法及系统。
背景技术：
可变形桁架(variablegeometrytruss,vgt)是指桁架的空间结构确定而几何参数可变，通过调节部分杆件的长度实现对末端负载位置和姿态的移动。它是空间操作机构的一种，能够实现在轨操作、在轨加注等任务。
在vgt的逆运动学问题上，miura和furuya提出了一组用三个角度值作为变量描述vgt非线性几何模型的方程，furuya提出了一组用三个角度值作为变量描述vgt非线性几何模型的方程，hughes给出了基于vgt-bay的逆运动学，实现了以vgt为执行机构的空间机器人操作。在miura和furuya工作的基础上，huang等人给出了一种空间八面体型vgt的定义和正逆运动学模型，并给出了自由漂浮系统的动力学模型和控制方法，通过物理仿真验证了有效性。上述方法能够解决特定的单级可变形桁架逆运动学问题。此外，也有通过蛮力搜索的方式，即对所有关节变量同时求解，寻找一组能够使末端到达期望位置的关节变量值，该方法思想简单，但运算效率低。现有算法在解决多级桁架逆运动学问题时，矩阵维数提高，运算量过大，难以得到逆运动学结果。
技术实现要素：
本发明解决的技术问题是：克服现有技术的不足，提供了一种空间可变形桁架逆运动学计算方法及系统，解决了多级桁架逆运动学运算量大、不易求解的问题。
本发明目的通过以下技术方案予以实现：根据本发明的一个方面，提供了一种空间可变形桁架逆运动学计算方法，所述方法包括如下步骤：
(1)将n个桁架单元级联构成空间可变形桁架，其中，n为偶数；
(2)根据第1级桁架单元建立基坐标系oxyz；
(3)建立各级桁架单元坐标系，将第n+1级桁架单元的固定杆面中心记为可变形桁架的工作点e，第n/2+1级桁架单元的固定杆面中心记为空间可变形桁架工作点e1；
(4)获取oxyz坐标系中的工作目标点d(xd,yd,zd)，计算第一阶段目标点d1(xd,0,zd/2)，其中，第一阶段为计算第1级到第n/2级桁架单元逆运动学；
(5)计算从第1级到第n/2级的桁架单元的每一级的工作点在其对应桁架单元坐标系的平行于主动杆面的坐标系平面的投影点，这些投影点构成坐标集合，根据坐标集合得到单级工作空间密度函数；
(6)在步骤(5)中的每一级桁架单元的坐标集合中使得e1在d1处概率密度函数最大的一个工作点作为该级工作点；
(7)由步骤(6)中的工作点得到第1至n/2级桁架单元各杆长度；
(8)计算从第n/2+1级到第n级的桁架单元的每一级的工作点在其对应桁架单元坐标系的平行于主动杆面的坐标系平面的投影点，这些投影点构成坐标集合，根据坐标集合得到单级工作空间密度函数；
(9)在步骤(8)中的每一级桁架单元的坐标集合中使得e在d处概率密度函数最大的一个工作点作为该级工作点；
(10)由步骤(9)中的工作点得到第n/2+1至n级桁架单元主动杆长度。
上述空间可变形桁架逆运动学计算方法中，在步骤(1)中，每级桁架单元为六面体，包括底面、四个侧面和顶面，底面和顶面为正方形，各边长度相同，且为常量w，自下而上，将各级桁架单元依次记为第1,…,n级桁架单元；对第1至n/2级桁架单元，前后两个面为主动杆面，左右两个面为被动杆面；对第n/2+1至n级桁架单元，左右两个面为主动杆面，前后两个面为被动杆面；主动杆面中，侧棱和面对角线均为主动调节长度的杆，且每组相对的杆的杆长同步运动；被动杆面中，侧棱与主动杆面的侧棱共用，面对角线为长度被动变化的杆；第k级桁架单元中，底面和顶面分别记为第k级和第k+1级固定杆面，主动杆面、被动杆面分别记为第k级主动杆面和第k级被动杆面，其中，k＝1,…,n。
上述空间可变形桁架逆运动学计算方法中，在步骤(2)中，原点o位于第1级桁架单元的固定杆面正方形中心，ox轴平行于第1级桁架单元的主动杆面与底面的交线，oy轴平行于第1级桁架单元的被动杆面与底面的交线，oz轴沿底面垂线向上，ox、oy、oz轴构成右手坐标系。
上述空间可变形桁架逆运动学计算方法中，在步骤(3)中，对第k级桁架单元，建立坐标系okxkykzk，原点ok位于第k级固定杆平面正方形中心，xkokyk平面与第k级固定杆面重合，okxk轴平行于第k级主动杆面与固定杆面的交线，okyk轴平行于第k级被动杆面与固定杆面的交线，okzk轴沿第k级平面的垂线向上，okxk、okyk、okzk轴构成右手坐标系；坐标系o1x1y1z1与基坐标系oxyz重合；在坐标系okxkykzk中，落在xkokyk平面第一至四象限内的固定杆面正方形顶点分别记为点pak、pbk、pck、pdk，记第k+1级固定杆面中心为第k级桁架单元的工作点ek，即ek与ok+1重合。
上述空间可变形桁架逆运动学计算方法中，在步骤(5)中，在第u级桁架单元中，保持y轴坐标不变，投影在xvovzv平面内，以长度步长rs遍历主动杆杆长，计算第u级桁架单元的工作点eu在平面xvovzv中投影的坐标集合，得到第u级桁架单元在平面xvovzv中的单级工作空间密度函数
上述空间可变形桁架逆运动学计算方法中，在步骤(8)中，在第k级桁架单元中，保持x轴坐标不变，投影在ykokzk平面内，以长度步长rs遍历主动杆杆长，计算第k级桁架单元的工作点ek在平面ykokzk中投影的坐标集合，得到第k级桁架单元在平面ykokzk中的单级工作空间密度函数
上述空间可变形桁架逆运动学计算方法中，在步骤(5)或步骤(8)中，单级桁架单元投影正运动学计算步骤为：
已知第k级桁架单元在平面xoz上的投影，已知第k级固定杆面投影产生的左、右顶点的坐标分别为[xkl,ykl]t和[xkr,ykr]t，主动杆的投影杆长分别为rk1，rk2和rk3，求取第k+1级固定杆面投影产生的左、右顶点坐标[xk+1,l,yk+1,l]t和[xk+1,r,yk+1,r]t：
当xkl＝xkr时，
当xkl≠ykl时，
式中，pa、pb、pa、pb、pc为中间变量，表达式为
取xkl≠ykl时存在两组解
η1、η2为用于判断增根的中间变量，选取z1、z2中的正值对应的解即为所求解点
类似地，当xkl＝xkr时，
当xkl≠ykl时，
式中，pa、pb、pa、pb、pc为中间变量，表达式为
取xkl≠ykl时存在两组解
η3、η4为用于判断增根的中间变量，选取z3、z4中的正值对应的解即为所求解点
即得到第k+1级固定杆面投影产生的左、右顶点坐标。
上述空间可变形桁架逆运动学计算方法中，在步骤(5)中，单级工作空间密度函数
(1)遍历主动杆杆长，可得到第u级桁架单元的工作点eu在平面xvovzv中的可达点集合，即eu的工作空间d，其中包含的可达点数量为σ，根据d中x的变化范围[xmin,xmax]和y的变化范围[ymin,ymax]，可以得到覆盖d的矩形v：
v＝[xmin,xmax]×[ymin,ymax]；
(2)沿x、y轴以步长s0将v划分为ωx×ωy个小矩形网格，ωx、ωy分别为沿x轴和y轴的网格数，并沿坐标轴的方向将网格小矩形编号为(ωi,ωj)(i＝1,…,ωx；j＝1,…,ωy)，对应的覆盖范围为：
式中，δx、δy表示网格小矩形的覆盖区域，δxe、δye为边界处小网格的覆盖区域，表达式为：
落在网格小矩形(ωi,ωj)内的可达点数量为σij，则网格(ωi,ωj)的工作空间密度值gij为：
(3)所有网格的工作空间密度构成了eu在坐标系xvovzv内的工作空间密度函数：
其中，下标u表示第u级桁架单元，上标v表示相对于平面xvovzv。
上述空间可变形桁架逆运动学计算方法中，在步骤(6)中，e1在d1处概率密度函数最大的一个工作点计算步骤为：
(1)令循环变量k＝1；
(2)在第k级桁架单元中，遍历ek在平面xkokzk中的工作点，对每个工作点，计算其后各级桁架单元投影工作空间概率密度函数的卷积为：
式中，*为二维卷积，计算方法如下：
设有两个二维离散函数：
式中，fa、fb为函数f的坐标范围，ga、gb为函数g的坐标范围，f、g的二维离散卷积为：
式中，fs与gs分别为f与g的周期化函数，表达式为
式中，
m＝fa+ga-1,n＝fb+gb-1；
m×n是二维离散卷积的一个周期；
(2)计算
(3)若k≤n/2，则转步骤(1)，否则转步骤(4)；
(4)选取其中使得e1在d1处概率密度函数最大的一个工作点作为第k级工作点。
上述空间可变形桁架逆运动学计算方法中，在步骤(9)中，e在d处概率密度函数最大的一个工作点计算步骤为：
(1)令循环变量k＝n/2+1；
(2)在第k级桁架单元中，遍历ek在平面xkokzk中的工作点，对每个工作点，计算其后各级桁架单元投影工作空间概率密度函数的卷积为：
式中，*为二维卷积，计算方法如下：
设有两个二维离散函数：
它们的二维离散卷积为：
式中，*与*分别为f与g的周期化函数；
式中，
m＝fa+ga-1,n＝fb+gb-1
m×n是二维离散卷积的一个周期；
(2)计算
(3)若k≤n，则转步骤(1)，否则转步骤(4)；
(4)选取其中使得e在d处概率密度函数最大的一个工作点作为第k级工作点。
上述空间可变形桁架逆运动学计算方法中，在步骤(7)或步骤(10)中，主动杆杆长计算步骤为：
已知投影中，第k级桁架单元的固定杆面四个顶点坐标分别为pak、pbk、pck和pdk，第k+1级桁架单元的固定杆面投影产生的四个顶点坐标pak+1、pbk+1、pck+1和pdk+1，则第k级桁架单元所有杆的长度分别为：
式中，||·||2为2-范数，即两点的距离。
根据本发明的另一方面，还提供了一种空间可变形桁架逆运动学计算系统，包括：第一模块，用于将n个桁架单元级联构成空间可变形桁架，其中，n为偶数；第二模块，用于根据第1级桁架单元建立基坐标系oxyz；第三模块，用于建立各级桁架单元坐标系，将第n+1级桁架单元的固定杆面中心记为可变形桁架的工作点e，第n/2+1级桁架单元的固定杆面中心记为空间可变形桁架工作点e1；第四模块，用于获取oxyz坐标系中的工作目标点d(xd,yd,zd)，计算第一阶段目标点d1(xd,0,zd/2)，其中，第一阶段为计算第1级到第n/2级桁架单元逆运动学；第五模块，用于计算从第1级到第n/2级的桁架单元的每一级的工作点在其对应桁架单元坐标系的平行于主动杆面的坐标系平面的投影点，这些投影点构成坐标集合，根据坐标集合得到单级工作空间密度函数；第六模块，用于在第五模块中的每一级桁架单元的坐标集合中使得e1在d1处概率密度函数最大的一个工作点作为该级工作点；第七模块，用于由第六模块中的工作点得到第1至n/2级桁架单元各杆长度；第八模块，用于计算从第n/2+1级到第n级的桁架单元的每一级的工作点在其对应桁架单元坐标系的平行于主动杆面的坐标系平面的投影点，这些投影点构成坐标集合，根据坐标集合得到单级工作空间密度函数；第九模块，用于在第八模块中的每一级桁架单元的坐标集合中使得e在d处概率密度函数最大的一个工作点作为该级工作点；第十模块，用于由第九模块中的工作点得到第n/2+1至n级桁架单元主动杆长度。
根据本发明的又一个方面，提供了一个或多个机器可读介质，其上存储有指令，当由一个或多个处理器执行时，使得设备执行本发明的一个方面的一个或多个的方法。
本发明与现有技术相比具有如下有益效果：
(1)相比现有的蛮力搜索方法在解决多级桁架逆运动学问题时运算量大，计算时间长，本发明利用工作空间概率密度函数的扩散过程，将多级桁架逆运动学的求解转化为递推形式，显著减小了运算量，缩短了计算时间；
(2)本发明通过选取使得末端工作点在期望位置工作空间密度函数值最大的一组工作点作为各级桁架单元的工作点，解决了现有的单级桁架逆运动学算法无法解决多级桁架的冗余性问题；
(3)本发明能够推广到同类型的多级重复结构逆运动学的求解中。
附图说明
通过阅读下文优选实施方式的详细描述，各种其他的优点和益处对于本领域普通技术人员将变得清楚明了。附图仅用于示出优选实施方式的目的，而并不认为是对本发明的限制。而且在整个附图中，用相同的参考符号表示相同的部件。在附图中：
图1为本发明的空间可变形桁架逆运动学计算方法的流程图；
图2为本发明的空间可变形桁架结构示意图；
图3为本发明的空间可变形桁架坐标系示意图；
图4为本发明的空间可变形桁架目标点示意图；
图5为本发明的单级投影正运动学示意图；
图6(a)为本发明的单级工作点坐标集合示意图；
图6(b)为本发明的单级工作空间密度函数示意图；
图7(a)为本发明的平面扩散过程单级工作空间密度函数示意图；
图7(b)为本发明的平面扩散过程两级工作空间密度函数示意图；
图7(c)为本发明的平面扩散过程三级工作空间密度函数示意图；
图8为本发明的单级桁架单元投影逆运动学；
图9为本发明的实施例计算所得逆运动学。
具体实施方式
下面将参照附图更详细地描述本公开的示例性实施例。虽然附图中显示了本公开的示例性实施例，然而应当理解，可以以各种形式实现本公开而不应被这里阐述的实施例所限制。相反，提供这些实施例是为了能够更透彻地理解本公开，并且能够将本公开的范围完整的传达给本领域的技术人员。需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。下面将参考附图并结合实施例来详细说明本发明。
图1为本发明的实施例提供的空间可变形桁架逆运动学计算方法的流程图。如图1所示，该方法包括以下步骤：
(1)将n个桁架单元级联构成空间可变形桁架，n为偶数，每级桁架单元为六面体，包括底面、四个侧面和顶面。底面和顶面为正方形，各边长度相同，且为常量w。自下而上，将各级桁架单元依次记为第1,…,n级桁架单元。对第1至n/2级桁架单元，前后两个面为主动杆面，左右两个面为被动杆面。对第n/2+1至n级桁架单元，左右两个面为主动杆面，前后两个面为被动杆面。主动杆面中，侧棱和面对角线均为主动调节长度的杆，且每组相对的杆的杆长同步运动。被动杆面中，侧棱与主动杆面的侧棱共用，面对角线为长度被动变化的杆。第k级(k＝1,…,n)桁架单元中，底面和顶面分别记为第k级和第k+1级固定杆面，主动杆面、被动杆面分别记为第k级主动杆面和第k级被动杆面；
(2)建立基坐标系oxyz，原点o位于第1级固定杆面正方形中心，ox轴平行于第1级主动杆面与底面的交线，oy轴平行于第1级被动杆面与底面的交线，oz轴沿底面垂线向上，ox、oy、oz轴构成右手坐标系；
(3)建立各级桁架单元坐标系，对第k级(k＝1,…,n)桁架单元，建立坐标系okxkykzk，原点ok位于第k级固定杆平面正方形中心，xkokyk平面与第k级固定杆面重合，okxk轴平行于第k级主动杆面与固定杆面的交线，okyk轴平行于第k级被动杆面与固定杆面的交线，okzk轴沿第k级平面的垂线向上，okxk、okyk、okzk轴构成右手坐标系。坐标系o1x1y1z1与基坐标系oxyz重合。在坐标系okxkykzk中，落在xkokyk平面第一至四象限内的固定杆面正方形顶点分别记为点pak、pbk、pck、pdk，记第k+1级固定杆面中心为第k级桁架单元的工作点ek，即ek与ok+1重合。第n+1级固定杆面中心记为可变形桁架的工作点e，第n/2+1级固定杆面中心为空间可变形桁架工作点e1；
(4)获取oxyz坐标系中的工作目标点d(xd,yd,zd)，计算第一阶段目标点d1(xd,0,zd/2)，其中，第一阶段为计算第1级到第n/2级桁架单元逆运动学；
(5)计算从第1级到第n/2级的桁架单元的每一级的工作点在其对应桁架单元坐标系的平行于主动杆面的坐标系平面的投影点，这些投影点构成坐标集合，根据坐标集合得到单级工作空间密度函数；定义第u级桁架单元在平面xvovzv中的工作空间密度函数
(6)在步骤(5)中的每一级桁架单元的坐标集合中使得e1在d1处概率密度函数最大的一个工作点作为该级工作点；
(7)由步骤(6)中的工作点得到第1至n/2级桁架单元各杆长度；
(8)开始第二阶段计算，第二阶段是指计算从第n/2+1级到第n级的桁架单元逆运动学；计算从第n/2+1级到第n级的桁架单元的每一级的工作点在其对应桁架单元坐标系的平行于主动杆面的坐标系平面的投影点，这些投影点构成坐标集合，根据坐标集合得到单级工作空间密度函数；在第k级桁架单元中，保持x轴坐标不变，投影在ykokzk平面内，以长度步长rs遍历主动杆杆长，计算第k级桁架单元的工作点ek在平面ykokzk中投影的坐标集合，得到第k级桁架单元在平面ykokzk中的单级工作空间密度函数
(9)在步骤(8)中的每一级桁架单元的坐标集合中使得e在d处概率密度函数最大的一个工作点作为该级工作点；
(10)由步骤(9)中的工作点得到第n/2+1至n级桁架单元主动杆长度。
在本实施例中，首先对整个桁架编号，建立坐标系，在原长状态下，桁架为竖直状态。获取oxyz坐标系中的工作目标点d(xd,yd,zd)，对发明中的n级空间可变形桁架，设置第1至n/2级桁架单元中，平行于xoz平面的两个侧面为同步运动的主动杆面，平行于yoz的两个侧面为被动杆面；第n/2+1至n级桁架单元中，平行于yoz的两个侧面为同步运动的主动杆面，平行于xoz平面的两个侧面为被动杆面。从而可以将计算过程分为两个阶段，首先计算第一阶段目标点d1(xd,0,zd/2)到第1级桁架单元，即1～n/2级桁架单元的逆运动学，而后计算目标点d到第n/2+1级桁架单元，即第n/2+1～n级桁架单元的逆运动学。
在第一阶段，以第n/2级桁架工作点e1到达d1为目标。在该阶段，保持y轴坐标不变，又由于平行于xoz平面的两个侧面为同步运动的主动杆面，则只需计算第1至n/2级桁架单元在xoz平面的投影，首先计算第1至n/2级各级工作点在xoz平面的投影点及工作空间密度函数。然后根据扩散过程，依次寻找各级工作点中使得e1的工作空间概率密度在d1处取值最大的一个工作点作为当级的工作点。从而，得到第1至n/2级的投影工作点。令第n/2级顶面与xoy平面平行，依次求取各级顶点位置，进而求得第1至n/2级的各级主动杆和被动杆杆长。
在第二阶段，以第n级桁架工作点e到达d为目标。在该阶段，保持x轴坐标不变，又由于平行于yoz平面的两个侧面为同步运动的主动杆面，则只需计算第n/2+1至n级桁架单元在yoz平面的投影，首先计算第n/2+1至n级各级工作点在yoz平面的投影点及工作空间密度函数。然后根据扩散过程，依次寻找各级工作点中使得e的工作空间概率密度在d处取值最大的一个工作点作为当级的工作点。从而，得到第1至n/2级的投影工作点。令第n/2级顶面与xoy平面平行，依次求取各级顶点位置，进而求得第n/2+1至n级的各级主动杆和被动杆杆长。
综合第一阶段和第二阶段的结果，可以在已知oxyz坐标系一个末端工作目标点的条件下，求得各级杆长，即空间可变形桁架逆运动学。
如图2所示，为本发明的空间可变形桁架结构示意图。
本发明中，空间可变形桁架由n级桁架单元构成，每级桁架单元为六面体，包括底面、四个侧面和顶面。底面和顶面为正方形，各边长度相同，且为常量w。四个侧面中，一组相对的为主动杆面，另一组相对的为被动杆面。主动杆面中，侧棱和面对角线均为主动调节长度的杆，且每组相对的杆的杆长同步运动。被动杆面中，侧棱与主动杆面的侧棱共用，面对角线为长度被动变化的杆。自下而上，将各级桁架单元依次记为第1,…,n级桁架单元。其中，第k级(k＝1,…,n)桁架单元中，底面和顶面分别记为第k级和第k+1级固定杆面，主动杆面、被动杆面分别记为第k级主动杆面和第k级被动杆面。
如图3所示，为本发明的空间可变形桁架坐标系示意图。
本发明中，建立空间可变形桁架坐标系。建立基坐标系oxyz，原点o位于第1级固定杆面正方形中心，ox轴平行于第1级主动杆面与底面的交线，oy轴平行于第1级被动杆面与底面的交线，oz轴沿底面垂线向上，ox、oy、oz轴构成右手坐标系。建立第k级(k＝1,…,n)桁架单元坐标系okxkykzk，原点ok位于第k级固定杆平面正方形中心，xkokyk平面与第k级固定杆面重合，okxk轴平行于第k级主动杆面与固定杆面的交线，okyk轴平行于第k级被动杆面与固定杆面的交线，okzk轴沿第k级平面的垂线向上，okxk、okyk、okzk轴构成右手坐标系。坐标系o1x1y1z1与基坐标系oxyz重合。在坐标系okxkykzk中，落在xkokyk平面第一至四象限内的固定杆面正方形顶点分别记为点pak、pbk、pck、pdk，记第k+1级固定杆面中心为第k级桁架单元的工作点ek，即esk与ok+1重合。第n+1级固定杆面中心记为可变形桁架的工作点e。
图4为本发明的空间可变形桁架目标点示意图。
获取oxyz坐标系中的工作目标点d(xd,yd,zd)，为桁架末端期望的工作点，计算第一阶段目标点d1(xd,0,zd/2)，为第一阶段末端期望的工作点，其中，第一阶段为计算第1级到第n/2级桁架单元逆运动学。
如图5所示，为本发明的单级投影正运动学示意图。
已知第k级桁架单元在平面xoz上的投影，已知第k级固定杆面投影产生的左、右顶点的坐标分别为[xkl,ykl]t和[xkr,ykr]t，主动杆的投影杆长分别为rk1，rk2和rk3。求取第k+1级固定杆面投影产生的左、右顶点坐标[xk+1,l,yk+1,l]t和[xk+1,r,yk+1,r]t。由几何关系易见
求解该联立方程组。
当xkl＝xkr时，
当xkl≠ykl时，
式中pa、pb、pa、pb、pc为中间变量，表达式为
pb＝2papb-2pbxkl-2ykl
可见取不等号时的解集存在两组解
η1、η2为用于判断增根的中间变量，选取z1、z2中的正值对应的解即为所求解点
采取类似的方法，解方程组
并选择满足要求的解点，可求得
如图6(a)为本发明的单级工作点坐标集合示意图；如图6(b)为本发明的单级工作空间密度函数示意图。
工作空间密度函数是对工作空间精确性的度量，高密度点意味着工作点可以更精确地到达。在坐标系onxnynzn中，考虑第n级桁架单元，通过该级主动杆使得工作点en能够到达的位置称为可达点，可达点数量σ，所有可达点的集合为en的工作空间，记为d。根据d中x的变化范围[xmin,xmax]、y的变化范围[ymin,ymax]和z的变化范围[zmin,zmax]，可以得到覆盖d的立方体v
v＝[xmin,xmax]×[ymin,ymax]×[zmin,zmax]
则
式中，δx、δy、δz表示网格小矩形的覆盖区域，δxe、δye、δze为边界处小网格的覆盖区域，表达式为：
δx＝[xmin+(i-1)s0,xmin+is0)
δy＝[ymin+(j-1)s0,ymin+js0)
δz＝[zmin+(k-1)s0,zmin+ks0)
δxe＝[xmin+(i-1)s0,max{xmin+is0,xmax})
δye＝[ymin+(j-1)s0,max{ymin+js0,ymax})
δze＝[zmin+(k-1)s0,max{zmin+ks0,zmax})
落在网格小立方体(ωi,ωj,ωk)内的可达点数量为σijk，则网格(ωi,ωj,ωk)的工作空间密度gijk为
所有网格的工作空间密度构成了en在坐标系onxnynzn内的工作空间密度函数
当桁架工作点仅在平面运动时，退化为二维情形。即在坐标系xnonyn中，考虑第n级桁架单元，此时d和v为二维区域，沿x、y以步长s0将v划分为ωx×ωy个平面网格。落在网格(ωi,ωj)(i＝1,…,ωx；j＝1,…,ωy)内的可达点数量为σij，则网格(ωi,ωj)的工作空间密度为
所有网格的工作空间密度构成了en在坐标系xnonyn内的工作空间密度函数
如图7(a)为本发明的平面扩散过程单级工作空间密度函数示意图；如图7(b)为本发明的平面扩散过程两级工作空间密度函数示意图；如图7(c)为本发明的平面扩散过程三级工作空间密度函数示意图。
由p个桁架单元串联而成的空间可变形桁架，已知第n个(n＝1,…,p)桁架单元的工作点en在坐标系xnonyn内的工作空间密度函数为gn。由于工作空间密度函数具有马尔科夫性质，因而工作空间的生成可以描述为扩散过程。即p个桁架单元末端工作点ep在坐标系x1o1y1内的工作空间密度函数
其中，*为二维卷积，定义如下：
设有两个二维离散函数
式中，fa、fb为函数f的坐标范围，ga、gb为函数g的坐标范围，f、g的二维离散卷积为
式中，fs与gs分别为f与g的周期化函数，表达式为
其中，
m＝fa+ga-1,n＝fb+gb-1
m×n是二维离散卷积的一个周期。
如图8所示为单级桁架单元投影逆运动学。
已知投影中，第k级固定杆面四个顶点坐标分别为pak、pbk、pck和pdk，第k+1级固定杆面投影产生的四个顶点坐标pak+1、pbk+1、pck+1和pdk+1，则第k级桁架单元所有杆的长度分别为
raak＝||pak-pak+1||2
rbak＝||pbk-pak+1||2
rbbk＝||pbk-pbk+1||2
rcbk＝||pck-pbk+1||2
rcck＝||pck-pck+1||2
rcdk＝||pck-pdk+1||2
rddk＝||pdk-pdk+1||2
rdak＝||pdk-pak+1||2
式中，||·||2为2-范数，即两点的距离。
实施例：
以某空间可变形桁架逆运动学为例，桁架级数n＝6，固定杆长度w＝1m，给定期望末端工作点d(-2,-2,8)。
首先对各级桁架单元编号，建立坐标系。在原长状态下，桁架为竖直状态。设置第1至3级桁架单元中，平行于xoz平面的两个侧面为同步运动的主动杆面，平行于yoz的两个侧面为被动杆面；第4至6级桁架单元中，平行于yoz的两个侧面为同步运动的主动杆面，平行于xoz平面的两个侧面为被动杆面。
设定第一阶段目标点d1(-2,0,4)，第3级桁架单元工作点e1。保持y轴坐标不变，计算第1至n/2级桁架单元在xoz平面的投影。首先计算第1至3级各级桁架单元在xoz平面投影的工作点及工作空间密度函数。然后根据扩散过程，依次寻找各级工作点中使得e1的工作空间概率密度在d1处取值最大的一个工作点作为该级的工作点，求得e1坐标为(-1.9333,0,3.9873)。从而，得到第1至3级的投影工作点。令第3级顶面与xoy平面平行，依次求取各级顶点位置
pa1＝[0.50000.50000]t,pb1＝[-0.50000.50000]t
pc1＝[-0.5000-0.50000]t,pd1＝[0.5000-0.50000]t
pa2＝[-0.64420.50001.6485]t,pb2＝[-1.43490.50001.0362]t
pc2＝[-1.4349-0.50001.0362]t,pd2＝[-0.6442-0.50001.6485]t
pa3＝[-1.28400.50002.8185]t,pb3＝[-2.15930.50002.3348]t
pc3＝[-2.1593-0.50002.3348]t,pd3＝[-1.2840-0.50002.8185]t
pa4＝[-1.43330.50003.9873]t,pb4＝[-2.43330.50003.9873]t
pc4＝[-2.4333-0.50003.9873]t,pd4＝[-1.4333-0.50003.9873]t
进而求得第1至3级的各级主动杆和被动杆杆长。
raa1＝2.0067,rba1＝1.6548,rbb1＝1.3956,rcb1＝1.7169
rcc1＝1.3956,rcd1＝1.6548,rdd1＝2.0067,rda1＝2.2421
raa2＝1.3335,rba2＝1.7886,rbb2＝1.4870,rcb2＝1.7920
rcc2＝1.4870,rcd2＝1.7886,rdd2＝1.3335,rda2＝1.6668
raa3＝1.1784,rba3＝1.8050,rbb3＝1.6751,rcb3＝1.9509
rcc3＝1.6751,rcd3＝1.8050,rdd3＝1.1784,rda3＝1.5455
而后设定第二阶段目标点d(-2,-2,8)，第6级桁架单元工作点e。保持x轴坐标不变，计算第1至n/2级桁架单元在xoz平面的投影。首先计算第4至6级各级桁架单元在xoz平面投影的工作点及工作空间密度函数。然后根据扩散过程，依次寻找各级工作点中使得e的工作空间概率密度在d处取值最大的一个工作点作为该级的工作点，求得e坐标为(-1.9333,-1.9333,7.9747)。从而，得到第4至6级的投影工作点。令第6级顶面与xoy平面平行，依次求取各级顶点位置
pa5＝[-1.4333-0.64425.6358]t,pb5＝[-2.4333-0.64425.6358]t
pc5＝[-2.4333-1.43495.0235]t,pd5＝[-1.4333-1.43495.0235]t
pa6＝[-1.4333-1.28406.8058]t,pb6＝[-2.4333-1.28406.8058]t
pc6＝[-2.4333-2.15936.3221]t,pd6＝[-1.4333-2.15936.3221]t
pa7＝[-1.43330.50002.8185]t,pb7＝[-2.4333-1.43337.9747]t
pc7＝[-2.4333-2.43337.9747]t,pd7＝[-1.4333-2.43337.9747]t
进而求得第4至6级的各级主动杆和被动杆杆长。
raa4＝2.0067,rba4＝2.2421,rbb4＝2.0067,rcb4＝1.6548
rcc4＝1.3956,rcd4＝1.7169,rdd4＝1.3956,rda4＝1.6548
raa5＝1.3335,rba5＝1.6668,rbb5＝1.3335,rcb5＝1.7886
rcc5＝1.4870,rcd5＝1.7920,rdd5＝1.4870,rda5＝1.7886
raa6＝1.1784,rba6＝1.5455,rbb6＝1.1784,rcb6＝1.8050
rcc6＝1.6751,rcd6＝1.9509,rdd6＝1.6751,rda6＝1.8050
图9所示为实施例计算所得逆运动学。
本实施例还提供了一种空间可变形桁架逆运动学计算系统，包括：第一模块，用于将n个桁架单元级联构成空间可变形桁架，其中，n为偶数；第二模块，用于根据第1级桁架单元建立基坐标系oxyz；第三模块，用于建立各级桁架单元坐标系，将第n+1级桁架单元的固定杆面中心记为可变形桁架的工作点e，第n/2+1级桁架单元的固定杆面中心记为空间可变形桁架工作点e1；第四模块，用于获取oxyz坐标系中的工作目标点d(xd,yd,zd)，计算第一阶段目标点d1(xd,0,zd/2)，其中，第一阶段为计算第1级到第n/2级桁架单元逆运动学；第五模块，用于计算从第1级到第n/2级的桁架单元的每一级的工作点在其对应桁架单元坐标系的平行于主动杆面的坐标系平面的投影点，这些投影点构成坐标集合，根据坐标集合得到单级工作空间密度函数；第六模块，用于在第五模块中的每一级桁架单元的坐标集合中使得e1在d1处概率密度函数最大的一个工作点作为该级工作点；第七模块，用于由第六模块中的工作点得到第1至n/2级桁架单元各杆长度；第八模块，用于计算从第n/2+1级到第n级的桁架单元的每一级的工作点在其对应桁架单元坐标系的平行于主动杆面的坐标系平面的投影点，这些投影点构成坐标集合，根据坐标集合得到单级工作空间密度函数；第九模块，用于在第八模块中的每一级桁架单元的坐标集合中使得e在d处概率密度函数最大的一个工作点作为该级工作点；第十模块，用于由第九模块中的工作点得到第n/2+1至n级桁架单元主动杆长度。
本实施例还提供了一个或多个机器可读介质，其上存储有指令，当由一个或多个处理器执行时，使得设备执行航天器喷气控制器参数辅助设计方法。
现有的蛮力搜索方法在解决多级桁架逆运动学问题时运算量大，计算时间长，实施例利用工作空间概率密度函数的扩散过程，将多级桁架逆运动学的求解转化为递推形式，显著减小了运算量，缩短了计算时间；本实施例通过选取使得末端工作点在期望位置工作空间密度函数值最大的一组工作点作为各级桁架单元的工作点，解决了现有的单级桁架逆运动学算法无法解决多级桁架的冗余性问题；本实施例能够推广到同类型的多级重复结构逆运动学的求解中。
以上所述的实施例只是本发明较优选的具体实施方式，本领域的技术人员在本发明技术方案范围内进行的通常变化和替换都应包含在本发明的保护范围内。