docs: def alpha_s(scale, f, alphasMZ=0.1185, loop=3):
      """3-loop computation of alpha_s for f flavours
      with initial condition alpha_s(MZ) = 0.1185"""
      if scale == MZ and f == 5:
          return alphasMZ  # nothing to do
      _sane(scale, f)
      crd = rundec.CRunDec()
      if f == 5:
          return_value = crd.AlphasExact(alphasMZ, MZ, scale, f, loop)
      elif f == 6:
          crd.nfMmu.Mth = 170
          crd.nfMmu.muth = 170
          crd.nfMmu.nf = 6
          return_value = crd.AlL2AlH(alphasMZ, MZ, crd.nfMmu, scale, loop)
      elif f == 4:
          crd.nfMmu.Mth = 4.8
          crd.nfMmu.muth = 4.8
          crd.nfMmu.nf = 5
          return_value = crd.AlH2AlL(alphasMZ, MZ, crd.nfMmu, scale, loop)
      elif f == 3:
          crd.nfMmu.Mth = 4.8
          crd.nfMmu.muth = 4.8
          crd.nfMmu.nf = 5
          mc = 1.3
          asmc =  crd.AlH2AlL(alphasMZ, MZ, crd.nfMmu, mc, loop)
          crd.nfMmu.Mth = mc
          crd.nfMmu.muth = mc
          crd.nfMmu.nf = 4
          return_value = crd.AlH2AlL(asmc, mc, crd.nfMmu, scale, loop)
      else:
          raise ValueError("Invalid input: f={}, scale={}".format(f, scale))
      if return_value == 0:
          raise ValueError("Return value is 0, probably `scale={}` is too small.".format(scale))
      else:
          return return_value