How to solve this system of equation in sage?

edit

Considering the notation of your first (general) example, simply copy the following code (as a script) and adjust $k$ as desired.

k = 3
r = var('r')
c = list(var('c%d' % i) for i in [1..k])
l = list(var('l%d' % i) for i in [1..k])
x = list(var('x%d' % i) for i in [1..k])
eqn_list = list(sum(x[i]*l[j]^i for i in range(k))==c[j]+r^(j+1) for j in range(k))
sol = solve(eqn_list,x)

Output of sol:

sage: sol
[[x1 == ((l1^2*l2 - l1*l2^2)*r^3 + (c3*l2 - c2*l3)*l1^2 - (l1^2*l3 - l1*l3^2)*r^2 + (l2^2*l3 - l2*l3^2)*c1 - (c3*l2^2 - c2*l3^2)*l1 + (l2^2*l3 - l2*l3^2)*r)/(l1^2*(l2 - l3) + l2^2*l3 - l2*l3^2 - (l2^2 - l3^2)*l1), x2 == -((l1^2 - l2^2)*r^3 - (c2 - c3)*l1^2 - c3*l2^2 + c2*l3^2 - (l1^2 - l3^2)*r^2 + (l2^2 - l3^2)*c1 + (l2^2 - l3^2)*r)/(l1^2*(l2 - l3) + l2^2*l3 - l2*l3^2 - (l2^2 - l3^2)*l1), x3 == ((l1 - l2)*r^3 - (l1 - l3)*r^2 - (c2 - c3)*l1 + c1*(l2 - l3) - c3*l2 + c2*l3 + (l2 - l3)*r)/(l1^2*(l2 - l3) + l2^2*l3 - l2*l3^2 - (l2^2 - l3^2)*l1)]]

An alternative way would be to go via a matrix (which is notably faster).

k = 3
r = var('r')
c = list(var('c%d' % i) for i in [1..k])
l = list(var('l%d' % i) for i in [1..k])
A = matrix(list(l[j]^i for i in range(k)) for j in range(k))
v = vector(c[j]+r^(j+1) for j in range(k))
sol = A\v

Output of sol:

sage: sol
((r^3 - c1 + c3 - (r^2 - c1 + c2 - r)*(l1 - l3)/(l1 - l2) - r)*(l1^2 - (l1^2 - l2^2)*l1/(l1 - l2))/(l1^2 - l3^2 - (l1^2 - l2^2)*(l1 - l3)/(l1 - l2)) + c1 + (r^2 - c1 + c2 - r)*l1/(l1 - l2) + r, (r^3 - c1 + c3 - (r^2 - c1 + c2 - r)*(l1 - l3)/(l1 - l2) - r)*(l1^2 - l2^2)/((l1^2 - l3^2 - (l1^2 - l2^2)*(l1 - l3)/(l1 - l2))*(l1 - l2)) - (r^2 - c1 + c2 - r)/(l1 - l2), -(r^3 - c1 + c3 - (r^2 - c1 + c2 - r)*(l1 - l3)/(l1 - l2) - r)/(l1^2 - l3^2 - (l1^2 - l2^2)*(l1 - l3)/(l1 - l2)))

This is indeed the same solution.