system of ODE simpy dsolve

Question

How can I solve multiple ODEs? sympy.dsolve returns the same integrations constants, so I cannot solve it.

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
from sympy.interactive import printing
printing.init_printing(use_latex=True)
from sympy import *

x,l,F,k_0,q_0,n,a,C1,C2,C3,C4,Φ1,Φ2,D1,D2=symbols('x l F k_0 q_0 n a C1 C2 C3 C4 Φ_1 Φ_2 D1 D2')
n=5

q=0
k1=k_0
k2=n*k_0

eq1=Function('eq1')(x)
eq2=Function('eq2')(x)

w1=dsolve(k1.diff(x)*eq1.diff(x)+k1*eq1.diff(x,2)+q).rhs
w2=dsolve(k2.diff(x)*eq2.diff(x)+k2*eq2.diff(x,2)+q).rhs

display(w1)
display(w2)

The code returns:

𝐶1+𝐶2*𝑥

𝐶1+𝐶2*𝑥

While I would want something like:

𝐶1+𝐶2*𝑥

𝐶3+𝐶4*𝑥

Answer 1

I've found this as an answer, but maybe there is something simpler

char='D'
diff=k2.diff(x)*eq2.diff(x)+k2*eq2.diff(x,2)+q

def ODE(diff,char):
    func=dsolve(diff).rhs
    for i in range(1,len(func.free_symbols)):
        old=symbols('C{}'.format(i))
        new=symbols(char + '{}'.format(i))
        func=func.subs(old,new)
    return func

w2=ODE(diff,char)
display(w2)