docs: def pdf(self, d, n=None):
          r'''Computes the probability density function of a
          continuous particle size distribution at a specified particle diameter,
          an optionally in a specified basis. The evaluation function varies with
          the distribution chosen. The interconversion between distribution 
          orders is performed using the following formula [1]_:
              
          .. math::
              q_s(d) = \frac{x^{(s-r)} q_r(d) dd}
              { \int_0^\infty d^{(s-r)} q_r(d) dd}
          
          Parameters
          ----------
          d : float
              Particle size diameter, [m]
          n : int, optional
              None (for the `order` specified when the distribution was created),
              0 (number), 1 (length), 2 (area), 3 (volume/mass),
              or any integer, [-]
  
          Returns
          -------
          pdf : float
              The probability density function at the specified diameter and
              order, [-]
              
          Notes
          -----
          The pdf order conversions are typically available analytically after 
          some work. They have been verified numerically. See the various
          functions with names ending with 'basis_integral' for the formulations.
          The distributions normally do not have analytical limits for diameters
          of 0 or infinity, but large values suffice to capture the area of the
          integral.
          
          Examples
          --------
          >>> psd = PSDLognormal(s=0.5, d_characteristic=5E-6, order=3)
          >>> psd.pdf(1e-5)
          30522.765209509154
          >>> psd.pdf(1e-5, n=3)
          30522.765209509154
          >>> psd.pdf(1e-5, n=0)
          1238.661379483343
          
          References
          ----------
          .. [1] Masuda, Hiroaki, Ko Higashitani, and Hideto Yoshida. Powder 
             Technology: Fundamentals of Particles, Powder Beds, and Particle 
             Generation. CRC Press, 2006.
          '''
          ans = self._pdf(d=d)
          if n is not None and n != self.order:
              power = n - self.order
              numerator = d**power*ans
              denominator = self._pdf_basis_integral_definite(d_min=0.0, d_max=self.d_excessive, n=power)
              ans = numerator/denominator
          # Handle splines which might go below zero
          ans = max(ans, 0.0)
          if self.truncated:
              if d < self.d_min or d > self.d_max:
                  return 0.0
              ans = (ans)/(self._cdf_d_max - self._cdf_d_min)
          return ans