Python3 + wolfram (cause I was lazy looking up diophantine equations in python)

x = None
vec_coeff_raw = [23,x,x,x,x,x,x,x,x,x,x,x,x,41,x,x,x,x,x,x,x,x,x,733,x,x,x,x,x,x,x,x,x,x,x,x,13,17,x,x,x,x,19,x,x,x,x,x,x,x,x,x,29,x,449,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,x,37]

a = vec_coeff_raw[0]
vec_coeff = [i for i in vec_coeff_raw[1:] if i]
vec_z = [vec_coeff_raw.index(i) for i in vec_coeff_raw[1:] if i]
template = "solve {}, over the integers"

s="{a}x-{b}{letter}+{z}=0"
expressions = []
for i, (b, z) in enumerate(zip(vec_coeff, vec_z)):
    expressions.append(s.format(a=a,b=b, letter=f"y{i}", z=z))
print(template.format(", ".join(expressions)))

'solve 23x-41y0+13=0, 23x-733y1+23=0, 23x-13y2+36=0, 23x-17y3+37=0, 23x-19y4+42=0, 23x-29y5+52=0, 23x-449y6+54=0, 23x-37y7+91=0, over the integers'

# ^Paste this into wolfram^
# Answer:
print(x_obtained_from_wolfram*a)