【问题标题】:Solving system of quadratic equations求解二次方程组
【发布时间】:2011-04-11 17:57:46
【问题描述】:

任何人都可以找到解决以下系统的方法吗?我尝试了Reduce,但评估需要一段时间,所以我不确定它是否有效

terms = {{g^2, g h, h^2, -o^2, -o p, -p^2}, {g^2, g k, 
    k^2, -o^2, -o q, -q^2}, {g^2, g m, m^2, -o^2, -o r, -r^2}, {g^2, 
    g n, n^2, -o^2, -o s, -s^2}, {h^2, h k, 
    k^2, -p^2, -p q, -q^2}, {h^2, h m, m^2, -p^2, -p r, -r^2}, {h^2, 
    h n, n^2, -p^2, -p s, -s^2}, {k^2, k m, 
    m^2, -q^2, -q r, -r^2}, {k^2, k n, n^2, -q^2, -q s, -s^2}, {m^2, 
    m n, n^2, -r^2, -r s, -s^2}};
vars = Variables@Flatten@terms;
coefs = Array[c, Dimensions[terms]];
eqs = MapThread[#1.#2 == 0 &, {terms, coefs}];
Reduce[eqs, vars, Reals]

【问题讨论】:

  • 所有 c[,] 都不相关?
  • 全零是一种解决方案,可能是无趣的。对于解决方案 vars,-vars 也是一种解决方案。因此,您可以节省制定问题的工作,以便仅找到具有非负分量的解决方案。如果您只是在寻找通用解决方案,Solve 将比 Reduce 更有效。否则,您的问题将是一个棘手的问题,特别是如果您打算将解决方案作为 60 个参数的函数并找出对它们的限制以使解决方案成为现实。
  • c 是与问题相关的浮点值,如果将Real 更改为Complex,问题是否变得容易处理?
  • 可能无法用符号系数处理。可能会看:forums.wolfram.com/mathgroup/archive/2011/Mar/msg00410.html

标签: wolfram-mathematica


【解决方案1】:

您可以从优化的角度解决问题,构建 r.h.s. 的平方和。你的方程式。

定义你的矩阵:

mat[{g_, h_, k_, m_, n_, o_, p_, q_, r_, s_}] := {{g^2, g h, 
    h^2, -o^2, -o p, -p^2}, {g^2, g k, k^2, -o^2, -o q, -q^2}, {g^2, 
    g m, m^2, -o^2, -o r, -r^2}, {g^2, g n, 
    n^2, -o^2, -o s, -s^2}, {h^2, h k, k^2, -p^2, -p q, -q^2}, {h^2, 
    h m, m^2, -p^2, -p r, -r^2}, {h^2, h n, 
    n^2, -p^2, -p s, -s^2}, {k^2, k m, m^2, -q^2, -q r, -r^2}, {k^2, 
    k n, n^2, -q^2, -q s, -s^2}, {m^2, m n, n^2, -r^2, -r s, -s^2}};

现在将求解代数方程的代码定义为约束优化:

Clear[SolveAlgebraic];
SolveAlgebraic[
  coefs_ /; Dimensions[coefs] == {10, 6} && MatrixQ[coefs, NumberQ], 
  opts : OptionsPattern[NMinimize]] := 
 Module[{g, h, k, m, n, o, p, q, r, s, eqs, vars, val, sol},
  eqs = MapThread[#1.#2 &, {mat[
      vars = {g, h, k, m, n, o, p, q, r, s}], coefs}];
  {val, sol} = 
   NMinimize[{Total[eqs^2], 
     vars.vars > 1 && Apply[And, Thread[vars >= 0]]}, vars, opts];
  {val, vars /. sol}
  ]

现在定义一个函数,用给定的解构造一组 c[,]:

CoefficientWithSolution[sol_ /; Length[sol] == 10] := 
 Block[{cc, 
   v}, ((Array[
       cc, {10, 6}]) /. (First[
        Quiet@Solve[(MapThread[
              Dot, {mat[Array[v, 10]], Array[cc, {10, 6}]}] == 0 // 
            Thread), Array[cc, {10, 6}] // Flatten]] /. 
       Thread[Array[v, 10] -> (sol)]) /. _cc :> 1)]

生成矩阵:

In[188]:= coefs = 
 CoefficientWithSolution[{1, 2, 3, 4, 5, 6, 7, 8, 9, 10}]

Out[188]= {{1, 1, 1, 1, 1, -(71/49)}, {1, 1, 1, 1, 1, -(71/64)}, {1, 
  1, 1, 1, 1, -(23/27)}, {1, 1, 1, 1, 1, -(13/20)}, {1, 1, 1, 1, 
  1, -(43/32)}, {1, 1, 1, 1, 1, -(28/27)}, {1, 1, 1, 1, 
  1, -(4/5)}, {1, 1, 1, 1, 1, -(11/9)}, {1, 1, 1, 1, 1, -(19/20)}, {1,
   1, 1, 1, 1, -(11/10)}}

求解具有更高工作精度的方程,并强制转换为机器数:

In[196]:= SolveAlgebraic[coefs, WorkingPrecision -> 30] // N

Out[196]= {1.41177*10^-28, {0.052633, 0.105266, 0.157899, 0.210532, 
  0.263165, 0.315798, 0.368431, 0.421064, 0.473697, 0.52633}}

验证是否找到了预期的解决方案:

In[197]:= Rest[Last[%]]/First[Last[%]]

Out[197]= {2., 3., 4., 5., 6., 7., 8., 9., 10.}

希望这会有所帮助。

【讨论】:

    猜你喜欢
    • 1970-01-01
    • 2013-03-02
    • 1970-01-01
    • 1970-01-01
    • 1970-01-01
    • 2017-04-24
    • 2018-08-05
    • 2017-08-06
    • 2016-08-06
    相关资源
    最近更新 更多