Ask Your Question

Revision history [back]

Here is a way to build a linear system as in the question.

First, define a function a to create symbolic variables $a_{i,j}$.

Here we represent $a_{i, j}$ as a_i_j, and we use m for the minus sign, so m1 stands for $-1$.

def a(t):
    return SR.var(('a' + '_{}'*len(t)).format(*[str(i).replace('-','m') for i in t]))

Test it:

sage: a((1, 1)), a((12, 437)), a((-1, 2)), a((-1, 2, 8)), a((2, 4, -1, 3))
(a_1_1, a_12_437, a_m1_2, a_m1_2_8, a_2_4_m1_3)

Then, we need a function to build the linear system in the question.

def thesystem(p, q):
    s = []
    for i in range((p-1)//2 + 1):
        for j in range(q-2):
            s.append(a((i, j)) - 3*a((i-1, j)) - 4*a((i, j-1)) + 10*a((i-1, j-1)) == 0)
    for j in range(q-1):
        s.append(a((-1, j)) + a(((p-1)//2 - 1, j)) == 0)
    for i in range((p-1)//2 + 1):
        s.append(a((i, -1)) - a((i, q-2)) == 0)
    return s

Test it with $(p, q) = (5, 4)$:

sage: s_5_4 = thesystem(5, 4)
sage: s_5_4
[a_0_0 - 4*a_0_m1 - 3*a_m1_0 + 10*a_m1_m1 == 0,
 -4*a_0_0 + a_0_1 + 10*a_m1_0 - 3*a_m1_1 == 0,
 -3*a_0_0 + 10*a_0_m1 + a_1_0 - 4*a_1_m1 == 0,
 10*a_0_0 - 3*a_0_1 - 4*a_1_0 + a_1_1 == 0,
 -3*a_1_0 + 10*a_1_m1 + a_2_0 - 4*a_2_m1 == 0,
 10*a_1_0 - 3*a_1_1 - 4*a_2_0 + a_2_1 == 0,
 a_1_0 + a_m1_0 == 0,
 a_1_1 + a_m1_1 == 0,
 a_1_2 + a_m1_2 == 0,
 -a_0_2 + a_0_m1 == 0,
 -a_1_2 + a_1_m1 == 0,
 -a_2_2 + a_2_m1 == 0]

Now a function to extract the variables from the linear system:

def sysvar(s):
    v = set()
    for e in s:
        v.update(e.variables())
    return sorted(v, key=str)

Try it:

sage: v_5_4 = sysvar(s_5_4)
sage: v_5_4
[a_0_0,
 a_0_1,
 a_0_2,
 a_0_m1,
 a_1_0,
 a_1_1,
 a_1_2,
 a_1_m1,
 a_2_0,
 a_2_1,
 a_2_2,
 a_2_m1,
 a_m1_0,
 a_m1_1,
 a_m1_2,
 a_m1_m1]

Now we can solve the linear system:

sage: solve(s_5_4, v_5_4)
[[a_0_0 == r8,
 a_0_1 == -1/5*r5 + 17/5*r8,
 a_0_2 == -1/10*r5 + 2/5*r7 + 3/10*r8,
 a_0_m1 == -1/10*r5 + 2/5*r7 + 3/10*r8,
 a_1_0 == r5,
 a_1_1 == 17/5*r5 + 1/5*r8,
 a_1_2 == r7,
 a_1_m1 == r7,
 a_2_0 == -1/20*r5 + 1/4*r6 - 3/20*r8,
 a_2_1 == r6,
 a_2_2 == -61/80*r5 + 1/16*r6 + 5/2*r7 - 3/80*r8,
 a_2_m1 == -61/80*r5 + 1/16*r6 + 5/2*r7 - 3/80*r8,
 a_m1_0 == -r5,
 a_m1_1 == -17/5*r5 - 1/5*r8,
 a_m1_2 == -r7,
 a_m1_m1 == -17/50*r5 + 4/25*r7 + 1/50*r8]]

The solution involves some real parameters rj whose numbering may vary depending on the history of your Sage session.

Alternatively, we could build a matrix for the linear system.

For this, we define a function, inspired by this answer to this question, but taking advantage of the fact that there are no constant terms here.

def matrix_from_system(s, v):
    return matrix([[(e.lhs() - e.rhs()).coefficient(a) for a in v] for e in s])

The following matrix is obtained by applying this function to our example:

sage: matrix_from_system(s_5_4, v_5_4)
[ 1  0  0 -4  0  0  0  0  0  0  0  0 -3  0  0 10]
[-4  1  0  0  0  0  0  0  0  0  0  0 10 -3  0  0]
[-3  0  0 10  1  0  0 -4  0  0  0  0  0  0  0  0]
[10 -3  0  0 -4  1  0  0  0  0  0  0  0  0  0  0]
[ 0  0  0  0 -3  0  0 10  1  0  0 -4  0  0  0  0]
[ 0  0  0  0 10 -3  0  0 -4  1  0  0  0  0  0  0]
[ 0  0  0  0  1  0  0  0  0  0  0  0  1  0  0  0]
[ 0  0  0  0  0  1  0  0  0  0  0  0  0  1  0  0]
[ 0  0  0  0  0  0  1  0  0  0  0  0  0  0  1  0]
[ 0  0 -1  1  0  0  0  0  0  0  0  0  0  0  0  0]
[ 0  0  0  0  0  0 -1  1  0  0  0  0  0  0  0  0]
[ 0  0  0  0  0  0  0  0  0  0 -1  1  0  0  0  0]

and we can then use linear algebra.

For more on linear algebra, the approach using matrices and vectors, and the approach of solving linear systems, see