docs: def vertical_modes(self):
          """ Calculate standard vertical modes. Simply
              the eigenvectors of the stretching matrix S """
  
          evals,evecs = np.linalg.eig(-self.S)
  
          asort = evals.argsort()
  
          # deformation wavenumbers and radii
          self.kdi2 = evals[asort]
          self.radii = np.zeros_like(self.kdi2)
          self.radii[0] = self.g*self.H/np.abs(self.f) # barotropic def. radius
          self.radii[1:] = 1./np.sqrt(self.kdi2[1:])
  
          # eigenstructure
          self.pmodes = evecs[:,asort]
  
          # normalize to have unit L2-norm
          Ai = (self.H / (self.Hi[:,np.newaxis]*(self.pmodes**2)).sum(axis=0))**0.5
          self.pmodes = Ai[np.newaxis,:]*self.pmodes