二、Galerkin法的应用

为了求出{a}，引入无量纲时间τ=ωt，则（5-10）式可以写成：

式中[q₀  q₁   … q_N  ]^T的各分量为联轴器恢复力（，，，τ）的傅立叶级数各简谐分量的系数分量。[f₁₀  f₁₁    …  f_1N  ]^T，…[f_jo  f_j1    …  f_jN ]^T，…[f_n0  f_n1   … f_nN  ]^T分别为弹性恢复力和阻尼恢复力的F₁=k_（n+1）1+c_（n+1）1…，F_j=k_（n+1）j+ c_（n+1）…，F_n=k_（n+1）n+c_（n+1）n的傅立叶级数各简谐分量的系数向量。

由以上推导可知，求解带有非线性联轴器轴系的稳态响应，可以归结为联合求解2N+1个非线性代数方程组（5-15）和求解n个微分方程组（5-4)而重点又在于由（5-15）式求得2N+1个未知量［a₀  a₁   … a_N  ］^T。由(5-15）式可知，它不同于一般的非线性方程组，在一般的非线性方程组中，未知量是显含的，而（5-15）式的未知量{a}是隐含的，而且{q｝是{a｝的非线性函数。不能通过将yb的傅立叶展开式代人恢复力Q的表达式中来直接求得，这是因为Q表达式中的振幅不能表示成{a}的显函数的缘故。此要想由（5-15）式解得｛a}，还需要解决两个问题：一是如何求（，，，τ），二是用什么来解隐含{a}的方程组（5-15)。下面就来解决这些问题。

三、ILM法

将非线性方程组简写成：

由于向量函数{r（{a}）}中非线性项{q（{a}）}的复杂性，可行的求解方法是采用数值法。

对于式（5-17）的数值法求解问题，我们将其转化为与之等效的目标函数：

的极小值问题。即非线性方程组（5-17）的解{a}可以通过求此极小值问题的最小二乘解得到。在此采用一种改进的LM算法（ILM）来求解这一问题。ILM法的基本思想简述如下：

当LM法用于（5-18）问题时，迭代公式为：

式中μ_k是LM参数，I为单位矩阵，{p^（k）（μ_k）}为搜索逼近（5-17）式解的校正向量，它依赖于参数μ_k值，控制搜索方向和长度。

LM算法的优点是，在计算（5-20）式时，可以不考虑[J（{a^（k）}）]是否满秩，μ_k可以起到（{a}）下降参数的作用，{p^（k）（μ_k）}是（{a}）在{a^（k）}处的下降方向，而且总存在这样的μ≥0，使（{a^（k）}+ p^（k））<{a^（k）}）（对所有k），随着μ值的增大，{p^（k）（μ_k）}的方向越来越接近（{a}的负梯度方向，此时，收敛速度变慢，但放宽了对初始近似{a^（0）}的要求。缺点是选择μ满足（{a^（k+1）}）<（{a^（k）}）时，在一个迭代步骤中需要多次求解方程组（5-20），大大增加了计算量，因此有必要对LM法进行改进，以使选择μk时即不增加太多的计算量，又能保证满足下降性质，即有一定的收敛速度，又能改善（5-20）方程组的病态性。改进后的迭代公式为：

式中为矩阵^（k）的元素。，分别为矩阵L^（k），^（k）的元素。公式（5-23）和（5-27）和（5-28）中不含参数μ_k，因此改变μ_k时，只需重新计算(5-24),而不需要重新进行三角分解。显然，改进后的算法减少了由于选择合适的μ_k而带来的计算量。这种算法与LM法的不同之处在于用正定矩阵L_k 代替了单位矩阵I，它不仅调整了矩阵^（k）的主对角元素，而且对整个矩阵进行了调整，它比LM法有更好的收敛速度。

为了得到合适的μ_k值，使得搜索向量{p^（k）（μ_k）}更快更好地满足（{a^（k）}+ p^（k））<({a^（k）}），加速收敛，改善LM法对初始值{a^（0）}要求过高的缺点，本文采用一种选择合适μk值的新算法，对LM法进行了改进，选择计算步骤如下：

①.取定初始近似{a^（0）}和初始值，为了调整μ_k值，引人调整系数v,可以取得小一些，取=10^-2，v=5。

②.按（5-21) 式首先算得［J（{a^（0）}）]，然后按（5-26）式算得，再根据（5-27）、（5-28）式算得矩阵L(0）和^（0）中的各元素，按(5-25）、(5-23）和（5-24）算得{P^（0）（）}；按（5-17）算得{r（{a^（0）}）}；按（5-18）算（{a^（0）}）。并计算{a^（1）}={a^（0）}+{P^（0）（）}和（{a^（1）}）

③.(l）如果（{a^（1）}）≤（{a^（0）}），则以/v代替重新计算{P^（0）（/v）}和新的{^（1）}，（{^（1）}）。（a）如果有（{^（1）}）≤（{a^（0）}），则取定μ₀=/v；（b）如果（{^（1）}）>（{a^（0）}），则取定μ₀=，（2）（{a^（1）}）>（{a^（0）}），则需增大μ值，取=v^m，以代替，逐次对m=1，2，…时的代入（5-24）进行计算，最后取定使不等式（{a^（0）}）+{ P^（0）（）}）<（{a^（0）}）成立的最小正整数对应的=作为μ₀以及{a^（1）}={a^（0）}+{P（）}。

④.假定已求得第k次近似{ a^（k）}和相应的=μ_k-1，算出[J（{ a^（k）}）]，^（k），L^（k），^（k），r（{ a^（k）}）}，（{ a^（k）}），然后算得{这p（）}，并计算{ a^（k+1）}={ a^（k）}+{P（）}和（{ a^（k+1）}）

⑤.（1）如果（{ a^（k+1）}）≤（{ a^（k）}），则减少以代替重新计算{p（）}和对应的{ ^（k+1）}及（{ ^（k+1）}）（a）如果（{ ^（k+1）}）≤（{ a^（k）}），则取定μ_k=；（b）如果（{ ^（k+1）}）>（{ a^（k）}），则取定μ_k=；（2）如果（{a^（k+1）}）>（{ a^（k）}），则增大值，可取=v^m,以代替,逐次时m=1，2，…时的代入（5-24）式进行计算，最后取定使不定式({ a^（k）}+{p（）}<（{ a^（k）}）成立的最小正整数所对应的作为以及{a^(k+1)}={a^(k)}+{p()}。

⑥.目标函数（{a}）的k次迭代和k+1次迭代如果满足绝对误差和相对误差精度，即可将{ a^（k+1）}作为近似解。如果不满足，则以{ a^（k+1）}代替{ a^（k）}，μ_k代替μ_k-1去执行步骤⑤，直到满足要求为止。

至此，（5-17）式的求解方法已经建立起来，日（5-17）式求解{a}还需要求出（，，，τ），{q}和（f_j）（j=1，2，…，n）。另外，在用（5-21）式求[J（{a}）]时遇到r_i（{a}）对隐含向量{a}求偏导数的麻烦，因此在计算[J（{a}）]时，需要做些处理。由（5-15）式可知，它是在频域内的2N+1个非线性代数方程组，而非线性迟滞恢复力（，，，τ）是时域函数，因此在求解（5-15）时需要进行频域与时域相互变换，并做一些数学上的处理，下面逐一叙述。

上一页

下一页