code: def _expm(A, use_exact_onenorm):
      # Core of expm, separated to allow testing exact and approximate
      # algorithms.
  
      # Avoid indiscriminate asarray() to allow sparse or other strange arrays.
      if isinstance(A, (list, tuple, np.matrix)):
          A = np.asarray(A)
      if len(A.shape) != 2 or A.shape[0] != A.shape[1]:
          raise ValueError("expected a square matrix")
  
      # Trivial case
      if A.shape == (1, 1):
          out = [[np.exp(A[0, 0])]]
  
          # Avoid indiscriminate casting to ndarray to
          # allow for sparse or other strange arrays
          if isspmatrix(A):
              return A.__class__(out)
  
          return np.array(out)
  
      # Ensure input is of float type, to avoid integer overflows etc.
      if (
          isinstance(A, np.ndarray) or isspmatrix(A)
      ) and not np.issubdtype(A.dtype, np.inexact):
          A = A.astype(float)
  
      # Detect upper triangularity.
      structure = UPPER_TRIANGULAR if _is_upper_triangular(A) else None
  
      if use_exact_onenorm == "auto":
          # Hardcode a matrix order threshold for exact vs. estimated one-norms.
          use_exact_onenorm = A.shape[0] < 200
  
      # Track functions of A to help compute the matrix exponential.
      h = _ExpmPadeHelper(
          A, structure=structure, use_exact_onenorm=use_exact_onenorm
      )
  
      # Try Pade order 3.
      eta_1 = max(h.d4_loose, h.d6_loose)
      if eta_1 < 1.495585217958292e-002 and _ell(h.A, 3) == 0:
          U, V = h.pade3()
          return _solve_P_Q(U, V, structure=structure)
  
      # Try Pade order 5.
      eta_2 = max(h.d4_tight, h.d6_loose)
      if eta_2 < 2.539398330063230e-001 and _ell(h.A, 5) == 0:
          U, V = h.pade5()
          return _solve_P_Q(U, V, structure=structure)
  
      # Try Pade orders 7 and 9.
      eta_3 = max(h.d6_tight, h.d8_loose)
      if eta_3 < 9.504178996162932e-001 and _ell(h.A, 7) == 0:
          U, V = h.pade7()
          return _solve_P_Q(U, V, structure=structure)
      if eta_3 < 2.097847961257068e000 and _ell(h.A, 9) == 0:
          U, V = h.pade9()
          return _solve_P_Q(U, V, structure=structure)
  
      # Use Pade order 13.
      eta_4 = max(h.d8_loose, h.d10_loose)
      eta_5 = min(eta_3, eta_4)
      theta_13 = 4.25
  
      # Choose smallest s>=0 such that 2**(-s) eta_5 <= theta_13
      if eta_5 == 0:
          # Nilpotent special case
          s = 0
      else:
          s = max(int(np.ceil(np.log2(eta_5 / theta_13))), 0)
      s = s + _ell(2**-s * h.A, 13)
      U, V = h.pade13_scaled(s)
      X = _solve_P_Q(U, V, structure=structure)
      if structure == UPPER_TRIANGULAR:
          # Invoke Code Fragment 2.1.
          X = _fragment_2_1(X, h.A, s)
      else:
          # X = r_13(A)^(2^s) by repeated squaring.
          for i in range(s):
              X = X.dot(X)
      return X
  
imports: from __future__ import division, print_function, absolute_import
  import math
  import numpy
  import scipy.special
  from scipy.linalg.basic import solve, solve_triangular
  from scipy.sparse.base import isspmatrix
  from scipy.sparse.construct import eye
  from scipy.sparse.linalg import spsolve
  import scipy.sparse
  import scipy.sparse.linalg
  from scipy.sparse.linalg.interface import LinearOperator
label_desc: Use image node_id 10 for calling a global function with example usage: _expm(A, use_exact_onenorm) and returns: X, np, _solve_P_Q, _solve_P_Q, _solve_P_Q, _solve_P_Q, A