Macaulay2 » Documentation
Packages » Oscillators :: findRealSolutions
next | previous | forward | backward | up | index | toc

findRealSolutions -- find real solutions, at least one per component for well-conditioned systems

Description

In this package the main use is to zero in on equilibrium points for oscillator systems associated to graphs.

We use this in the following example.

i1 : G = graph({0,1,2,3,4}, {{0,1},{1,2},{2,3},{0,3},{0,4}})

o1 = Graph{0 => {1, 3, 4}}
           1 => {0, 2}
           2 => {1, 3}
           3 => {0, 2}
           4 => {0}

o1 : Graph
i2 : S = oscRing(G, CoefficientRing => QQ, Reduced => true)

o2 = S

o2 : PolynomialRing
i3 : I = oscSystem(G,S)

                                                                            
o3 = ideal (- y  - y  - y , x y  - x y  + y , - x y  + x y  + x y  - x y , -
               1    3    4   2 1    1 2    1     2 1    1 2    3 2    2 3   
     ------------------------------------------------------------------------
                            2    2       2    2       2    2       2    2
     x y  + x y  + y , y , x  + y  - 1, x  + y  - 1, x  + y  - 1, x  + y  -
      3 2    2 3    3   4   1    1       2    2       3    3       4    4
     ------------------------------------------------------------------------
     1)

o3 : Ideal of S
i4 : dim I

o4 = 1
i5 : Jac = oscJacobian(G,S)

o5 = | -x_1-x_3-x_4 x_1                0                           
     | x_1          -x_1x_2-y_1y_2-x_1 x_1x_2+y_1y_2               
     | 0            x_1x_2+y_1y_2      -x_1x_2-x_2x_3-y_1y_2-y_2y_3
     | x_3          0                  x_2x_3+y_2y_3               
     | x_4          0                  0                           
     ------------------------------------------------------------------------
     x_3                x_4  |
     0                  0    |
     x_2x_3+y_2y_3      0    |
     -x_2x_3-y_2y_3-x_3 0    |
     0                  -x_4 |

             5      5
o5 : Matrix S  <-- S
i6 : assert(det Jac == 0)
i7 : assert(Jac - transpose Jac == 0)
i8 : realsols = findRealSolutions I
warning: some solutions are not regular: {14, 17, 21, 25, 26, 27, 28, 29, 30, 34, 43, 44, 45, 46}

o8 = {{1, 1, 1, -1, 0, 0, 0, 0}, {1, 1, 1, 1, 0, 0, 0, 0}, {-1, 1, -1, -1, 0,
     ------------------------------------------------------------------------
     0, 0, 0}, {-1, 1, -1, 1, 0, 0, 0, 0}, {0, -1, 0, -1, -1, 0, 1, 0}, {0,
     ------------------------------------------------------------------------
     -1, 0, -1, 1, 0, -1, 0}, {0, -1, 0, -1, 1, 0, -1, 0}, {0, -1, 0, 1, 1,
     ------------------------------------------------------------------------
     0, -1, 0}, {0, -1, 0, 1, -1, 0, 1, 0}, {0, -1, 0, 1, 1, 0, -1, 0},
     ------------------------------------------------------------------------
     {.203311, -.91733, -.203311, -1, -.979114, -.398129, .979114, 0},
     ------------------------------------------------------------------------
     {.203311, -.91733, -.203311, 1, -.979114, -.398129, .979114, 0},
     ------------------------------------------------------------------------
     {-.203311, -.91733, .203311, -1, -.979114, .398129, .979114, 0},
     ------------------------------------------------------------------------
     {-.203311, -.91733, .203311, 1, -.979114, .398129, .979114, 0}, {0, -1,
     ------------------------------------------------------------------------
     0, -1, -1, 0, 1, 0}, {0, -1, 0, 1, -1, 0, 1, 0}}

o8 : List
i9 : netList getAngles(4, realsols, Radians=>false)

     +-------+-------+-------+---+
o9 = |0      |0      |0      |180|
     +-------+-------+-------+---+
     |0      |0      |0      |0  |
     +-------+-------+-------+---+
     |180    |0      |180    |180|
     +-------+-------+-------+---+
     |180    |0      |180    |0  |
     +-------+-------+-------+---+
     |270    |180    |90     |180|
     +-------+-------+-------+---+
     |90     |180    |270    |180|
     +-------+-------+-------+---+
     |90     |180    |270    |180|
     +-------+-------+-------+---+
     |90     |180    |270    |0  |
     +-------+-------+-------+---+
     |270    |180    |90     |0  |
     +-------+-------+-------+---+
     |90     |180    |270    |0  |
     +-------+-------+-------+---+
     |281.731|203.461|101.731|180|
     +-------+-------+-------+---+
     |281.731|203.461|101.731|0  |
     +-------+-------+-------+---+
     |258.269|156.539|78.2694|180|
     +-------+-------+-------+---+
     |258.269|156.539|78.2694|0  |
     +-------+-------+-------+---+
     |270    |180    |90     |180|
     +-------+-------+-------+---+
     |270    |180    |90     |0  |
     +-------+-------+-------+---+
i10 : C = decompose I

                                                         2    2            
o10 = {ideal (y , y , y  + y , x  - 1, x  + 1, x  - x , x  + y  - 1), ideal
               4   2   1    3   4       2       1    3   3    3            
      -----------------------------------------------------------------------
                                                  2    2                     
      (y , y , y  + y , x  - 1, x  + 1, x  + x , x  + y  - 1), ideal (y , y ,
        4   2   1    3   4       2       1    3   3    3               4   2 
      -----------------------------------------------------------------------
                                         2    2                              
      y  + y , x  + 1, x  + 1, x  - x , x  + y  - 1), ideal (y , y , y  + y ,
       1    3   4       2       1    3   3    3               4   2   1    3 
      -----------------------------------------------------------------------
                                2    2                                       
      x  + 1, x  + 1, x  + x , x  + y  - 1), ideal (y , y  + y , x  - 1, x  +
       4       2       1    3   3    3               4   1    3   4       1  
      -----------------------------------------------------------------------
            2                                           2                 
      x , 2y  + x  - 1, 2x y  - y , x y  - x y  - y , 2x  - x  - 1, x x  +
       3    3    2        3 3    2   3 2    2 3    3    3    2       2 3  
      -----------------------------------------------------------------------
                  2    2                                              2     
      y y  - x , x  + y  - 1), ideal (y , y  + y , x  + 1, x  + x , 2y  + x 
       2 3    3   2    2               4   1    3   4       1    3    3    2
      -----------------------------------------------------------------------
                                           2                              2  
      - 1, 2x y  - y , x y  - x y  - y , 2x  - x  - 1, x x  + y y  - x , x  +
             3 3    2   3 2    2 3    3    3    2       2 3    2 3    3   2  
      -----------------------------------------------------------------------
       2
      y  - 1), ideal (y , y , y , y , x  - 1, x  - 1, x  - 1, x  - 1), ideal
       2               4   3   2   1   4       3       2       1            
      -----------------------------------------------------------------------
      (y , y , y , y , x  - 1, x  + 1, x  - 1, x  + 1), ideal (y , y , y ,
        4   3   2   1   4       3       2       1               4   3   2 
      -----------------------------------------------------------------------
      y , x  + 1, x  - 1, x  - 1, x  - 1), ideal (y , y , y , y , x  + 1, x 
       1   4       3       2       1               4   3   2   1   4       3
      -----------------------------------------------------------------------
      + 1, x  - 1, x  + 1)}
            2       1

o10 : List
i11 : for i in C list Jac % i

o11 = {| -2x_3-1 x_3  0    x_3  1  |, | -1   -x_3 0    x_3  1  |, | -2x_3+1
       | x_3     0    -x_3 0    0  |  | -x_3 0    x_3  0    0  |  | x_3    
       | 0       -x_3 2x_3 -x_3 0  |  | 0    x_3  0    -x_3 0  |  | 0      
       | x_3     0    -x_3 0    0  |  | x_3  0    -x_3 0    0  |  | x_3    
       | 1       0    0    0    -1 |  | 1    0    0    0    -1 |  | -1     
      -----------------------------------------------------------------------
      x_3  0    x_3  -1 |, | 1    -x_3 0    x_3  -1 |, | -1   -x_3 0    x_3  
      0    -x_3 0    0  |  | -x_3 0    x_3  0    0  |  | -x_3 2x_3 -x_3 0    
      -x_3 2x_3 -x_3 0  |  | 0    x_3  0    -x_3 0  |  | 0    -x_3 0    x_3  
      0    -x_3 0    0  |  | x_3  0    -x_3 0    0  |  | x_3  0    x_3  -2x_3
      0    0    0    1  |  | -1   0    0    0    1  |  | 1    0    0    0    
      -----------------------------------------------------------------------
      1  |, | 1    -x_3 0    x_3   -1 |, | -3 1  0  1  1  |, | 1  -1 0  -1 1 
      0  |  | -x_3 2x_3 -x_3 0     0  |  | 1  -2 1  0  0  |  | -1 2  -1 0  0 
      0  |  | 0    -x_3 0    x_3   0  |  | 0  1  -2 1  0  |  | 0  -1 2  -1 0 
      0  |  | x_3  0    x_3  -2x_3 0  |  | 1  0  1  -2 0  |  | -1 0  -1 2  0 
      -1 |  | -1   0    0    0     1  |  | 1  0  0  0  -1 |  | 1  0  0  0  -1
      -----------------------------------------------------------------------
      |, | -1 1  0  1  -1 |, | 3  -1 0  -1 -1 |}
      |  | 1  -2 1  0  0  |  | -1 2  -1 0  0  |
      |  | 0  1  -2 1  0  |  | 0  -1 2  -1 0  |
      |  | 1  0  1  -2 0  |  | -1 0  -1 2  0  |
      |  | -1 0  0  0  1  |  | -1 0  0  0  1  |

o11 : List
i12 : M = Jac % C_1

o12 = | -1   -x_3 0    x_3  1  |
      | -x_3 0    x_3  0    0  |
      | 0    x_3  0    -x_3 0  |
      | x_3  0    -x_3 0    0  |
      | 1    0    0    0    -1 |

              5      5
o12 : Matrix S  <-- S
i13 : M1 = matrix{{-1,0,0,0,1},{0,0,0,0,0},{0,0,0,0,0},{0,0,0,0,0},{1,0,0,0,-1}}

o13 = | -1 0 0 0 1  |
      | 0  0 0 0 0  |
      | 0  0 0 0 0  |
      | 0  0 0 0 0  |
      | 1  0 0 0 -1 |

               5       5
o13 : Matrix ZZ  <-- ZZ
i14 : use S

o14 = S

o14 : PolynomialRing
i15 : eigenvalues lift(sub(M - M1, x_3=>1), QQ)

o15 = {-2          }
      {-5.07407e-17}
      {2           }
      {0           }
      {0           }

o15 : VerticalList
i16 : eigenvalues lift(sub(M, x_3=>1), QQ)

o16 = {1.8662     }
      {-2.65544   }
      {1.03516e-16}
      {-1.21076   }
      {2.96672e-17}

o16 : VerticalList
i17 : eigenvalues lift(sub(M, x_3=>1/100000), QQ)

o17 = {-2          }
      {.0000173205 }
      {-.0000173205}
      {1.04779e-16 }
      {-2.24504e-17}

o17 : VerticalList
i18 : eigenvalues lift(sub(M, x_3=>100), QQ)

o18 = {199.752                   }
      {-200.252                  }
      {-1.49997                  }
      {-3.9177e-17+6.40979e-15*ii}
      {-3.9177e-17-6.40979e-15*ii}

o18 : VerticalList

Caveat

If the system is positive dimensional, then of course it won't find all roots. There might be real roots on a component, but none are found. If a component (even a point) is singular, then it might have problems, depending on the situation. In the latter case, there is a warning message emitted (about non-regular solutions).

See also

Ways to use findRealSolutions:

  • findRealSolutions(Graph)
  • findRealSolutions(Ideal)

For the programmer

The object findRealSolutions is a method function.


The source of this document is in /build/reproducible-path/macaulay2-1.25.06+ds/M2/Macaulay2/packages/Oscillators/Documentation.m2:373:0.