Quadrature Error: unsupported operand type(s) for /: 'function' and 'Symbol'

so I am new to Python and using it as a tool for research. So excuse my noob-ness :D

I wanted to plot a function which is defined by an integral (that cannot be solved analytically) and the independent variable is the lower bound of this integral. Before I get to that step, I wanted to check if the numerical integration works at a single point (i.e. lower limit of the integral is 4, say) just to see if everything is in order. In the code below 'Int' is the integrand, which has the functions 'b' and 'f' in it as well, both functions of 'r'.

import numpy as np
import matplotlib.pyplot as plt
import scipy as sp
import sympy as smp
from scipy.integrate import quad

r = smp.symbols('r', real = True)
s = 1
Mu = 0.5
lamb = 0.6
b = lambda r: ((r - s) * lamb + (s**2)*r)/(s*r)
f = lambda r: (1 - lamb/(s*r))**((s**2 + lamb)/lamb)
Int = 1/(r*((1 - b/r)*((r**2)/(Mu**2)*(1/f) - 1))**smp.Rational(1/2))
quad(Int, 4, smp.oo)

On running the above I get this error: "TypeError: unsupported operand type(s) for /: 'function' and 'Symbol'". Any beginner-friendly breakdown of this issue will be helpful, thanks :)

Package scipy and sympy are not really meant to work together. Below is a simple example of numerical integration based on scipy only.

import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate

s = 1
Mu = 0.5
lamb = 0.6

def b(r):
    return ((r - s) * lamb + (s ** 2) * r) / (s * r)

def f(r):
    return (1 - lamb / (s * r)) ** ((s ** 2 + lamb) / lamb)

def integrand(r): 
    return 1. / (r * np.sqrt((1 - b(r) / r) * ((r ** 2) / (Mu ** 2) * (1. / f(r)) - 1.)))

I, sI = integrate.quad(integrand, 4., np.inf)
# (0.12556464889996616, 1.585413063883688e-09)

Then you can vectorize the integration:

@np.vectorize
def solve(a):
    return integrate.quad(integrand, a, np.inf)[0]

And compute the integral for several lower bounds:

a = np.linspace(1., 20., 200)
I = solve(a)

Which renders as follow: