code: def mpf_ci_si(x, prec, rnd=round_fast, which=2):
      """
      Calculation of Ci(x), Si(x) for real x.
  
      which = 0 -- returns (Ci(x), -)
      which = 1 -- returns (Si(x), -)
      which = 2 -- returns (Ci(x), Si(x))
  
      Note: if x < 0, Ci(x) needs an additional imaginary term, pi*i.
      """
      wp = prec + 20
      sign, man, exp, bc = x
      ci, si = None, None
      if not man:
          if x == fzero:
              return (fninf, fzero)
          if x == fnan:
              return (x, x)
          ci = fzero
          if which != 0:
              if x == finf:
                  si = mpf_shift(mpf_pi(prec, rnd), -1)
              if x == fninf:
                  si = mpf_neg(
                      mpf_shift(mpf_pi(prec, negative_rnd[rnd]), -1)
                  )
          return (ci, si)
      # For small x: Ci(x) ~ euler + log(x), Si(x) ~ x
      mag = exp + bc
      if mag < -wp:
          if which != 0:
              si = mpf_perturb(x, 1 - sign, prec, rnd)
          if which != 1:
              y = mpf_euler(wp)
              xabs = mpf_abs(x)
              ci = mpf_add(y, mpf_log(xabs, wp), prec, rnd)
          return ci, si
      # For huge x: Ci(x) ~ sin(x)/x, Si(x) ~ pi/2
      elif mag > wp:
          if which != 0:
              if sign:
                  si = mpf_neg(mpf_pi(prec, negative_rnd[rnd]))
              else:
                  si = mpf_pi(prec, rnd)
              si = mpf_shift(si, -1)
          if which != 1:
              ci = mpf_div(mpf_sin(x, wp), x, prec, rnd)
          return ci, si
      else:
          wp += abs(mag)
      # Use an asymptotic series? The smallest value of n!/x^n
      # occurs for n ~ x, where the magnitude is ~ exp(-x).
      asymptotic = mag - 1 > math.log(wp, 2)
      # Case 1: convergent series near 0
      if not asymptotic:
          if which != 0:
              si = mpf_pos(mpf_ci_si_taylor(x, wp, 1), prec, rnd)
          if which != 1:
              ci = mpf_ci_si_taylor(x, wp, 0)
              ci = mpf_add(ci, mpf_euler(wp), wp)
              ci = mpf_add(ci, mpf_log(mpf_abs(x), wp), prec, rnd)
          return ci, si
      x = mpf_abs(x)
      # Case 2: asymptotic series for x >> 1
      xf = to_fixed(x, wp)
      xr = (MPZ_ONE << (2 * wp)) // xf  # 1/x
      s1 = MPZ_ONE << wp
      s2 = xr
      t = xr
      k = 2
      while t:
          t = -t
          t = (t * xr * k) >> wp
          k += 1
          s1 += t
          t = (t * xr * k) >> wp
          k += 1
          s2 += t
      s1 = from_man_exp(s1, -wp)
      s2 = from_man_exp(s2, -wp)
      s1 = mpf_div(s1, x, wp)
      s2 = mpf_div(s2, x, wp)
      cos, sin = mpf_cos_sin(x, wp)
      # Ci(x) = sin(x)*s1-cos(x)*s2
      # Si(x) = pi/2-cos(x)*s1-sin(x)*s2
      if which != 0:
          si = mpf_add(mpf_mul(cos, s1), mpf_mul(sin, s2), wp)
          si = mpf_sub(mpf_shift(mpf_pi(wp), -1), si, wp)
          if sign:
              si = mpf_neg(si)
          si = mpf_pos(si, prec, rnd)
      if which != 1:
          ci = mpf_sub(mpf_mul(sin, s1), mpf_mul(cos, s2), prec, rnd)
      return ci, si
  
imports: import operator
  import math
  from .backend import MPZ_ZERO, MPZ_ONE, BACKEND, xrange, exec_
  from .libintmath import gcd
  from .libmpf import ComplexResult, round_fast, round_nearest, negative_rnd, bitcount, to_fixed, from_man_exp, from_int, to_int, from_rational, fzero, fone, fnone, ftwo, finf, fninf, fnan, mpf_sign, mpf_add, mpf_abs, mpf_pos, mpf_cmp, mpf_lt, mpf_le, mpf_gt, mpf_min_max, mpf_perturb, mpf_neg, mpf_shift, mpf_sub, mpf_mul, mpf_div, sqrt_fixed, mpf_sqrt, mpf_rdiv_int, mpf_pow_int, to_rational
  from .libelefun import mpf_pi, mpf_exp, mpf_log, pi_fixed, mpf_cos_sin, mpf_cos, mpf_sin, mpf_sqrt, agm_fixed
  from .libmpc import mpc_one, mpc_sub, mpc_mul_mpf, mpc_mul, mpc_neg, complex_int_pow, mpc_div, mpc_add_mpf, mpc_sub_mpf, mpc_log, mpc_add, mpc_pos, mpc_shift, mpc_is_infnan, mpc_zero, mpc_sqrt, mpc_abs, mpc_mpf_div, mpc_square, mpc_exp
  from .libintmath import ifac
  from .gammazeta import mpf_gamma_int, mpf_euler, euler_fixed
label_desc: Use image node_id 16 for calling a global function with example usage: mpf_ci_si(x, prec, rnd, which) and returns: ci,  si, ci,  si, ci,  si, ci,  si, fninf,  fzero, x,  x, ci,  si