Ben G
Ben G

Reputation: 26789

Solving a simple matrix differential equation

I want to use Sympy to solve the differential matrix equation: du/dt = [[0,1],[1,0]]u with initial value u(0) = [[4],[2]].

The answer is

enter image description here enter image description here

So the complete final answer is: 3e^t[[1],[1]] + e^{-t}[[1],[-1]]

How could I solve this with SymPy?

Upvotes: 0

Views: 181

Answers (3)

Some more modifications

import sympy as sm

def matrix_differential_eq(A, x0):
    A = sm.Matrix(A)
    x_sym = sm.symbols("x1:" + str(A.cols+1), cls=sm.Function)
    print(x_sym)
    t = sm.symbols("t")

    xx = sm.Matrix([i(t) for i in x_sym])

    z = A@xx
    eq_all = [sm.Eq(xx[i, 0].diff(t), z[i]) for i in range(len(xx))]

    funct = [*xx]

    u = [sm.Function(f'x{i+1}') for i in range(len(x_sym))]
    iccs = {u[i](0): val for i, val in enumerate([*x0])}

    return sm.dsolve(eq_all, funct, ics = iccs)

A = [
    [0, 1],
    [1, 0]
    ]
x0 = [4, 2]

print(matrix_differential_eq(A, x0))

Upvotes: 0

Ben G
Ben G

Reputation: 26789

Here's my attempt for a generic function to solve these kind of equations, based on help from @Joao

from sympy import symbols, Eq, Function, init_printing
from sympy.solvers.ode.systems import dsolve_system

# Solves du/dt = Au for matrix A, with initial values vector for t=0.
def matrix_differential_eq(A, initial_vals):
    u_n = symbols("u_:" + str(A.cols), cls=Function)
    t = symbols("t")
    u = Matrix([u_n[i](t) for i in range(0, len(u_n))])

    mul = A @ u

    eqs = []
    for i in range (0,len(mul)):
        eqs.append(Eq(u[i].diff(t), mul[i]))

    ics = {}
    for i, val in enumerate(initial_vals):
        ics[u_n[i](0)] = val

    solutions = dsolve_system(eqs, ics=ics)
    return solutions[0]

A = Matrix([[0,1],[1,0]])
initial_vals = Matrix([4,2])
display(matrix_differential_eq(A, initial_vals))

Upvotes: 1

João Areias
João Areias

Reputation: 1470

This can be achieved with dsolve_system:

Code:

from sympy import symbols, Eq, Function, init_printing
from sympy.solvers.ode.systems import dsolve_system

init_printing()

u1, u2 = symbols("u1 u2", cls=Function)
t = symbols("t")


eqs = [
    Eq(u1(t).diff(t), u2(t)),
    Eq(u2(t).diff(t), u1(t))
]

funcs = [u1(t), u2(t)]
solutions = dsolve_system(eqs, ics={u1(0): 4, u2(0): 2})
solutions[0]

Output:

enter image description here

Output in latex:

\left[ u_{1}{\left(t \right)} = 3 e^{t} + e^{- t}, \  u_{2}{\left(t \right)} = 3 e^{t} - e^{- t}\right]

Edit For a general matrix, you can write this:

import numpy as np

from sympy import symbols, Eq, Function, init_printing
from sympy.solvers.ode.systems import dsolve_system

init_printing()

# Replace with your matrix
A = np.array([
    [0, 1],
    [1, 0]
])
initial_condition = np.array([4, 2])

# Defines the symbols to be used
u = [Function(f'u{i}') for i in range(A.shape[0])]
t = symbols("t")

# Write the expression for the system of differential equations
eqs = []
for i, row in enumerate(A):
    rhs = 0
    for j, coeff in enumerate(row):
        rhs += coeff * u[j](t)

    eqs.append(Eq(u[i](t).diff(t), rhs))

# Solves the system of differential equations
ics = {u[i](0): val for i, val in enumerate(initial_condition)}
solutions = dsolve_system(eqs, ics=ics)

solutions[0]

Upvotes: 2

Related Questions