docs: def contraction_conical_Crane(Di1, Di2, l=None, angle=None):
      r'''Returns loss coefficient for a conical pipe contraction
      as shown in Crane TP 410M [1]_ between 0 and 180 degrees.
  
      If :math:`\theta < 45^{\circ}`:
          
      .. math::
          K_2 = {0.8 \sin \frac{\theta}{2}(1 - \beta^2)}
  
      otherwise:
          
      .. math::
          K_2 = {0.5\sqrt{\sin \frac{\theta}{2}} (1 - \beta^2)}
          
      .. math::
          \beta = d_2/d_1
  
      Parameters
      ----------
      Di1 : float
          Inside pipe diameter of the larger, upstream, pipe, [m]    
      Di2 : float
          Inside pipe diameter of the smaller, downstream, pipe, [m]    
      l : float
          Length of the contraction, optional [m]
      angle : float
          Angle of contraction, optional [degrees]
  
      Returns
      -------
      K : float
          Loss coefficient in terms of the following (smaller) pipe [-]
  
      Notes
      -----
      Cheap and has substantial impact on pressure drop. Note that the 
      nomenclature in [1]_ is somewhat different - the smaller pipe is called 1,
      and the larger pipe is called 2; and so the beta ratio is reversed, and the
      fourth power of beta used in their equation is not necessary.
  
      Examples
      --------
      >>> contraction_conical_Crane(Di1=0.0779, Di2=0.0525, l=0)
      0.2729017979998056
  
      References
      ----------
      .. [1] Crane Co. Flow of Fluids Through Valves, Fittings, and Pipe. Crane,
         2009.
      '''
      if l is None and angle is None:
          raise Exception('One of `l` or `angle` must be specified')
      beta = Di2/Di1
      beta2 = beta*beta
      if angle is not None:
          angle = radians(angle)
          #l = (Di1 - Di2)/(2.0*tan(0.5*angle)) # L is not needed in this calculation
      elif l is not None:
          try:
              angle = 2.0*atan((Di1-Di2)/(2.0*l))
          except ZeroDivisionError:
              angle = pi
      if angle < 0.25*pi:
          # Formula 1
          K2 = 0.8*sin(0.5*angle)*(1.0 - beta2)
      else:
          # Formula 2
          K2 = 0.5*(sin(0.5*angle)**0.5*(1.0 - beta2))
      return K2