docs: def bm3_dPdV(v, v0, k0, k0p, precision=1.e-5):
      """
      calculate dP/dV for numerical calculation of bulk modulus
      according to test this differs from analytical result by 1.e-5
  
      :param v: volume
      :param v0: volume at reference conditions
      :param k0: bulk modulus at reference conditions
      :param k0p: pressure derivative of bulk modulus at different conditions
      :param precision: precision for numerical calculation (default = 1.e-5*v0)
      :return: dP/dV
      """
      def f_scalar(v, v0, k0, k0p, precision=precision):
          return derivative(bm3_p, v, args=(v0, k0, k0p), dx=v0 * precision)
      f_v = np.vectorize(f_scalar, excluded=[1, 2, 3, 4])
      return f_v(v, v0, k0, k0p, precision=precision)