diff --git a/doc/api/figure_api.rst b/doc/api/figure_api.rst
index 5dd3adbfec9f..9d00bad65d93 100644
--- a/doc/api/figure_api.rst
+++ b/doc/api/figure_api.rst
@@ -60,6 +60,8 @@ Annotating
    :nosignatures:
 
    Figure.colorbar
+   Figure.colorbar_bivar
+   Figure.colorbar_multivar
    Figure.legend
    Figure.text
    Figure.suptitle
@@ -254,6 +256,8 @@ Annotating
    :nosignatures:
 
    SubFigure.colorbar
+   SubFigure.colorbar_bivar
+   SubFigure.colorbar_multivar
    SubFigure.legend
    SubFigure.text
    SubFigure.suptitle
diff --git a/doc/api/pyplot_summary.rst b/doc/api/pyplot_summary.rst
index 97d9c576cc86..5885e747bed2 100644
--- a/doc/api/pyplot_summary.rst
+++ b/doc/api/pyplot_summary.rst
@@ -263,6 +263,8 @@ Colormapping
 
    clim
    colorbar
+   colorbar_bivar
+   colorbar_multivar
    gci
    sci
    get_cmap
diff --git a/lib/matplotlib/_constrained_layout.py b/lib/matplotlib/_constrained_layout.py
index ce488d555898..507204f3acc9 100644
--- a/lib/matplotlib/_constrained_layout.py
+++ b/lib/matplotlib/_constrained_layout.py
@@ -398,12 +398,18 @@ def make_layout_margins(layoutgrids, fig, renderer, *, w_pad=0, h_pad=0,
         # make margin for colorbars.  These margins go in the
         # padding margin, versus the margin for Axes decorators.
         for cbax in ax._colorbars:
+            if cbax._colorbar_info["type"] == "MultivarColorbar":
+                #  a matplotlib.colorbar.MultivarColorbar object
+                fig = cbax.axes[0].get_figure(root=False)
+                tightbbox = cbax.get_tightbbox(renderer, for_layout_only=True)
+                cbbbox = tightbbox.transformed(fig.transFigure.inverted())
+            else:
+                cbpos, cbbbox = get_pos_and_bbox(cbax, renderer)
             # note pad is a fraction of the parent width...
             pad = colorbar_get_pad(layoutgrids, cbax)
             # colorbars can be child of more than one subplot spec:
             cbp_rspan, cbp_cspan = get_cb_parent_spans(cbax)
             loc = cbax._colorbar_info['location']
-            cbpos, cbbbox = get_pos_and_bbox(cbax, renderer)
             if loc == 'right':
                 if cbp_cspan.stop == ss.colspan.stop:
                     # only increase if the colorbar is on the right edge
@@ -691,7 +697,7 @@ def reposition_colorbar(layoutgrids, cbax, renderer, *, offset=None, compress=Fa
     Parameters
     ----------
     layoutgrids : dict
-    cbax : `~matplotlib.axes.Axes`
+    cbax : `~matplotlib.colorbar.MultiColorbars` or `~matplotlib.axes.Axes`
         Axes for the colorbar.
     renderer : `~matplotlib.backend_bases.RendererBase` subclass.
         The renderer to use.
@@ -701,10 +707,18 @@ def reposition_colorbar(layoutgrids, cbax, renderer, *, offset=None, compress=Fa
     compress : bool
         Whether we're in compressed layout mode.
     """
+    cb_info = cbax._colorbar_info
+    if cb_info["type"] == "MultivarColorbar":
+        multi_cbar = cbax
+        cbaxes = multi_cbar.axes
+        multi = True
+    else:
+        multi = False
+        cbaxes = [cbax]
 
-    parents = cbax._colorbar_info['parents']
+    parents = cb_info['parents']
     gs = parents[0].get_gridspec()
-    fig = cbax.get_figure(root=False)
+    fig = cbaxes[0].get_figure(root=False)
     trans_fig_to_subfig = fig.transFigure - fig.transSubfigure
 
     cb_rspans, cb_cspans = get_cb_parent_spans(cbax)
@@ -712,11 +726,11 @@ def reposition_colorbar(layoutgrids, cbax, renderer, *, offset=None, compress=Fa
                                                  cols=cb_cspans)
     pb = layoutgrids[gs].get_inner_bbox(rows=cb_rspans, cols=cb_cspans)
 
-    location = cbax._colorbar_info['location']
-    anchor = cbax._colorbar_info['anchor']
-    fraction = cbax._colorbar_info['fraction']
-    aspect = cbax._colorbar_info['aspect']
-    shrink = cbax._colorbar_info['shrink']
+    location = cb_info['location']
+    anchor = cb_info['anchor']
+    fraction = cb_info['fraction']
+    aspect = cb_info['aspect']
+    shrink = cb_info['shrink']
 
     # For colorbars with a single parent in compressed layout,
     # use the actual visual size of the parent axis after apply_aspect()
@@ -736,14 +750,20 @@ def reposition_colorbar(layoutgrids, cbax, renderer, *, offset=None, compress=Fa
             # Keep the pb x-coordinates but use actual y-coordinates
             pb = Bbox.from_extents(pb.x0, actual_pos_fig.y0,
                                    pb.x1, actual_pos_fig.y1)
-        elif location in ('top', 'bottom'):
+        else:  # location in ('top', 'bottom'):
             # For horizontal colorbars, use the actual parent bbox width
             # for colorbar sizing
             # Keep the pb y-coordinates but use actual x-coordinates
             pb = Bbox.from_extents(actual_pos_fig.x0, pb.y0,
                                    actual_pos_fig.x1, pb.y1)
 
-    cbpos, cbbbox = get_pos_and_bbox(cbax, renderer)
+    if multi:
+        tightbbox = martist._get_tightbbox_for_layout_only(multi_cbar, renderer)
+        cbbbox = tightbbox.transformed(fig.transFigure.inverted())
+        cbpos = multi_cbar._get_original_position()
+        cbpos = cbpos.transformed(fig.transSubfigure - fig.transFigure)
+    else:
+        cbpos, cbbbox = get_pos_and_bbox(cbax, renderer)
 
     # Colorbar gets put at extreme edge of outer bbox of the subplotspec
     # It needs to be moved in by: 1) a pad 2) its "margin" 3) by
@@ -754,6 +774,8 @@ def reposition_colorbar(layoutgrids, cbax, renderer, *, offset=None, compress=Fa
         pbcb = pb.shrunk(fraction, shrink).anchored(anchor, pb)
         # The colorbar is at the left side of the parent.  Need
         # to translate to right (or left)
+        if multi:
+            pbcb.x1 = pbcb.x0 + cbbbox.width
         if location == 'right':
             lmargin = cbpos.x0 - cbbbox.x0
             dx = bboxparent.x1 - pbcb.x0 + offset['right']
@@ -768,6 +790,8 @@ def reposition_colorbar(layoutgrids, cbax, renderer, *, offset=None, compress=Fa
             pbcb = pbcb.translated(dx, 0)
     else:  # horizontal axes:
         pbcb = pb.shrunk(shrink, fraction).anchored(anchor, pb)
+        if multi:
+            pbcb.y1 = pbcb.y0 + cbbbox.height
         if location == 'top':
             bmargin = cbpos.y0 - cbbbox.y0
             dy = bboxparent.y1 - pbcb.y0 + offset['top']
@@ -781,14 +805,19 @@ def reposition_colorbar(layoutgrids, cbax, renderer, *, offset=None, compress=Fa
             offset['bottom'] += cbbbox.height + cbpad
             pbcb = pbcb.translated(0, dy)
 
-    pbcb = trans_fig_to_subfig.transform_bbox(pbcb)
-    cbax.set_transform(fig.transSubfigure)
-    cbax._set_position(pbcb)
-    cbax.set_anchor(anchor)
     if location in ['bottom', 'top']:
         aspect = 1 / aspect
-    cbax.set_box_aspect(aspect)
-    cbax.set_aspect('auto')
+    if multi:
+        new_bboxs = multi_cbar._get_tight_packing_inside(pbcb)
+    else:
+        new_bboxs = [pbcb]
+    for cbax, new_bbox in zip(cbaxes, new_bboxs):
+        transformed_bbox = trans_fig_to_subfig.transform_bbox(new_bbox)
+        cbax.set_transform(fig.transSubfigure)
+        cbax._set_position(transformed_bbox)
+        cbax.set_anchor(anchor)
+        cbax.set_box_aspect(aspect)
+        cbax.set_aspect('auto')
     return offset
 
 
diff --git a/lib/matplotlib/backend_bases.py b/lib/matplotlib/backend_bases.py
index 384987e3d036..00a8e585d3b1 100644
--- a/lib/matplotlib/backend_bases.py
+++ b/lib/matplotlib/backend_bases.py
@@ -3275,7 +3275,8 @@ def press_zoom(self, event):
         # to the edge of the Axes bbox in the other dimension. To do that we
         # store the orientation of the colorbar for later.
         parent_ax = axes[0]
-        if hasattr(parent_ax, "_colorbar"):
+        if hasattr(parent_ax, "_colorbar") and hasattr(parent_ax._colorbar,
+                                                       "orientation"):
             cbar = parent_ax._colorbar.orientation
         else:
             cbar = None
diff --git a/lib/matplotlib/colorbar.py b/lib/matplotlib/colorbar.py
index 100807dd7f18..fb72637897ae 100644
--- a/lib/matplotlib/colorbar.py
+++ b/lib/matplotlib/colorbar.py
@@ -15,6 +15,7 @@
 import logging
 
 import numpy as np
+from collections.abc import Sequence
 
 import matplotlib as mpl
 from matplotlib import _api, cbook, collections, colors, contour, ticker
@@ -28,6 +29,7 @@
 
 _log = logging.getLogger(__name__)
 
+
 _docstring.interpd.register(
     _make_axes_kw_doc="""
 location : None or {'left', 'right', 'top', 'bottom'}
@@ -114,7 +116,79 @@
 spacing : {'uniform', 'proportional'}
     For discrete colorbars (`.BoundaryNorm` or contours), 'uniform' gives each
     color the same space; 'proportional' makes the space proportional to the
-    data interval.""")
+    data interval.""",
+    _make_bivar_axes_kw_doc="""
+location : None or {'left', 'right', 'top', 'bottom'}
+    The location, relative to the parent Axes, where the colorbar Axes
+    is created. Also determines the position of the ticks and labels,
+    which will favour being away from the parent axes.
+
+fraction : float, default: 0.15
+    Fraction of original Axes to use for colorbar.
+
+shrink : float, default: 1.0
+    Fraction by which to multiply the size of the colorbar.
+
+pad : float, default: 0.05 if left or right, 0.15 if top or bottom
+    Fraction of original Axes between colorbar and new image Axes.
+
+anchor : (float, float), optional
+    The anchor point of the colorbar Axes.
+    Defaults to (0.0, 0.5) if left or right; (0.5, 1.0) if top or bottom.
+
+panchor : (float, float), or *False*, optional
+    The anchor point of the colorbar parent Axes. If *False*, the parent
+    axes' anchor will be unchanged.
+    Defaults to (1.0, 0.5) if left or right; (0.5, 0.0) if top or bottom.
+""",
+    _bivar_colormap_kw_doc="""
+ticklocations : tuple describing the ticklocation of the y and x axis
+    The first element must be {'auto', 'left', 'right'}
+    The second element must be {'auto', 'top', 'bottom''}.
+    If 'auto', the ticklocations are determined by the *location*.
+
+""",
+    _make_multivar_axes_kw_doc="""
+location : None or {'left', 'right', 'top', 'bottom'}
+    The location, relative to the parent Axes, where the colorbar Axes
+    is created.  It also determines the *orientation* of the colorbar
+    (colorbars on the left and right are vertical, colorbars at the top
+    and bottom are horizontal).  If None, the location will come from the
+    *orientation* if it is set (vertical colorbars on the right, horizontal
+    ones at the bottom), or default to 'right' if *orientation* is unset.
+
+orientation : None or {'vertical', 'horizontal'}
+    The orientation of the colorbars.  It is preferable to set the *location*
+    of the colorbar, as that also determines the *orientation*; passing
+    incompatible values for *location* and *orientation* raises an exception.
+
+fraction : float, default: 0.15
+    Fraction of original Axes to use for colorbars.
+
+shrink : float, default: 1.0
+    Fraction by which to multiply the size.
+
+aspect : float, default: 20
+    Ratio of long to short dimensions.
+
+pad : float, default: 0.05 if vertical, 0.15 if horizontal
+    Fraction of original Axes between colorbars and new image Axes.
+
+anchor : (float, float), optional
+    The anchor point of the colorbars Axes.
+    Defaults to (0.0, 0.5) if vertical; (0.5, 1.0) if horizontal.
+
+panchor : (float, float), or *False*, optional
+    The anchor point of the colorbars parent Axes. If *False*, the parent
+    axes' anchor will be unchanged.
+    Defaults to (1.0, 0.5) if vertical; (0.5, 0.0) if horizontal.
+
+major_pad : float
+    Spacing between colorbars along the long axis
+
+major_pad : float
+    Spacing between colorbars along the short axis
+""")
 
 
 def _set_ticks_on_axis_warn(*args, **kwargs):
@@ -531,7 +605,6 @@ def update_normal(self, mappable=None):
         if self.mappable.norm != self.norm:
             self.norm = self.mappable.norm
             self._reset_locator_formatter_scale()
-
         self._draw_all()
         if isinstance(self.mappable, contour.ContourSet):
             CS = self.mappable
@@ -1063,6 +1136,8 @@ def remove(self):
 
         try:
             ax = self.mappable.axes
+            if ax is None:
+                return
         except AttributeError:
             return
         try:
@@ -1357,9 +1432,575 @@ def drag_pan(self, button, key, x, y):
 ColorbarBase = Colorbar  # Backcompat API
 
 
+class BivarColorbar:
+    r"""
+    Draw a bivariate colorbar in an existing Axes.
+
+    Typically, bivariate colorbars are created using `.Figure.colorbar_bivar`
+    and associated with `.ColorizingArtist`\s (such as an
+    `.AxesImage` generated via `~.axes.Axes.imshow`).
+
+    Unlike `Colorbar`, `BivarColorbar` does not support
+    customizing the ticks, and ticks must be customized on the axes instead.
+    """
+
+    n_rasterize = 256  # rasterize solids if number of colors >= n_rasterize
+
+    def __init__(
+        self, ax, mappable,
+        *,
+        alpha=None,
+        location=None,
+        ticklocations=('auto', 'auto'),
+        aspect=1.0,
+                ):
+        """
+        Parameters
+        ----------
+        ax : `~matplotlib.axes.Axes`
+            The `~.axes.Axes` instance in which the colorbar is drawn.
+
+        mappable : `.ColorizingArtist`
+            The mappable whose colormap and norm will be used.
+
+        alpha : float
+            The colorbars transparency between 0 (transparent) and 1 (opaque).
+
+        location : None or {'left', 'right', 'top', 'bottom'}
+            Set the bivariate colorbars's location
+
+        Other Parameters
+        ----------------
+        ticklocations : tuple describing the ticklocation of the y and x axis
+            The first element must be {'auto', 'left', 'right'}
+            The second element must be {'auto', 'top', 'bottom''}.
+            If 'auto', the ticklocations are determined by the *location*.
+        """
+
+        self.ax = ax
+        self.location = location
+        if isinstance(mappable, mpl.colorizer.Colorizer):
+            mappable = mcolorizer.ColorizingArtist(mappable)
+
+        self.mappable = mappable
+        self.aspect = aspect
+        self.colorizer = mappable.colorizer
+
+        mappable.colorbar = self
+        mappable.colorbar_cid = mappable.callbacks.connect(
+            'changed', self.update_normals)
+
+        ticklocations = list(ticklocations)
+        if len(ticklocations) != 2:
+            raise ValueError("ticklocations must be a tuple of length 2")
+        _api.check_in_list(['auto', 'left', 'right'],
+                           ticklocation=ticklocations[0])
+        _api.check_in_list(['auto', 'top', 'bottom'],
+                           ticklocation=ticklocations[1])
+
+        location_ticklocs = _get_bivar_ticklocations_from_location(location)
+        for i in range(2):
+            if ticklocations[i] == 'auto':
+                ticklocations[i] = location_ticklocs[i]
+        self.ticklocations = ticklocations
+        self.ax.yaxis.set(label_position=self.ticklocations[0],
+                          ticks_position=self.ticklocations[0])
+        self.ax.xaxis.set(label_position=self.ticklocations[1],
+                          ticks_position=self.ticklocations[1])
+
+        self._image = None
+        self.alpha = None
+        # Call set_alpha to handle array-like alphas properly
+        self.set_alpha(alpha)
+
+        self.update_normals()  # also calls _draw_all()
+        self.ax._colorbar = self
+        self._interactive_funcs = ["_get_view", "_set_view",
+                                   "_set_view_from_bbox", "drag_pan"]
+        for x in self._interactive_funcs:
+            setattr(self.ax, x, getattr(self, x))
+        self.ax.cla = self._cbar_cla
+
+    @property
+    def aspect(self):
+        return self._aspect
+
+    @aspect.setter
+    def aspect(self, aspect):
+        aspect = float(aspect)
+        self._aspect = aspect
+        self.ax.set_box_aspect(aspect)
+        if hasattr(self.ax, "_colorbar_info"):
+            self.ax._colorbar_info["aspect"] = aspect
+
+    def _draw_all(self):
+        """
+        Calculate any free parameters based on the current cmap and norm,
+        and do all the drawing.
+        """
+
+        # transform from 0-1 to vmin-vmax:
+        if self.mappable.get_array() is not None:
+            self.mappable.autoscale_None()
+        if not self.colorizer.norm.scaled():
+            # If we still aren't scaled after autoscaling, use 0, 1 as default
+            self._set_view([0, 1, 0, 1])
+        n, m = self.colorizer.cmap.N, self.colorizer.cmap.M
+        x = self.colorizer.norm.norms[1].inverse(np.linspace(0, 1, m + 1))
+        y = self.colorizer.norm.norms[0].inverse(np.linspace(0, 1, n + 1))
+        X, Y = np.meshgrid(x, y)
+
+        if self.alpha is None:
+            lut = self.colorizer.cmap.lut
+        else:
+            lut = np.copy(self.colorizer.cmap.lut)
+            lut[:, :, 3] *= self.alpha
+
+        if n * m > self.n_rasterize:
+            rasterized = True
+        else:
+            rasterized = False
+
+        if self._image is not None:
+            self._image.remove()
+        self._image = self.ax.pcolormesh(
+                X, Y, lut,
+                alpha=self.alpha,
+                rasterized=rasterized,
+                edgecolors='none', shading='flat')
+        # Apply norm scaling (supports LogNorm etc.)
+        if getattr(self.colorizer.norm.norms[0], '_scale', None):
+            # use the norm's scale (if it exists and is not None):
+            self.ax.set_yscale(self.colorizer.norm.norms[0]._scale)
+        else:
+            # fallback for custom norms, or NoNorm()
+            self.ax.set_yscale(
+                'function',
+                functions=(
+                    self.colorizer.norm.norms[0],
+                    self.colorizer.norm.norms[0].inverse
+                )
+            )
+
+        if getattr(self.colorizer.norm.norms[1], '_scale', None):
+            # use the norm's scale (if it exists and is not None):
+            self.ax.set_xscale(self.colorizer.norm.norms[1]._scale)
+        else:
+            # fallback for custom norms, or NoNorm()
+            self.ax.set_xscale(
+                'function',
+                functions=(
+                    self.colorizer.norm.norms[1],
+                    self.colorizer.norm.norms[1].inverse
+                )
+            )
+        # Manually set limits (image is in Axes coordinates)
+        extent = [
+            self.colorizer.norm.norms[1].vmin,
+            self.colorizer.norm.norms[1].vmax,
+            self.colorizer.norm.norms[0].vmin,
+            self.colorizer.norm.norms[0].vmax,
+        ]
+        self.ax.set_ylim(extent[2:4])
+        self.ax.set_xlim(extent[0:2])
+
+    def set_xlabel(self, label):
+        self.ax.set_xlabel(label)
+
+    def set_ylabel(self, label):
+        self.ax.set_ylabel(label)
+
+    @property
+    def xaxis(self):
+        return self.ax.xaxis
+
+    @property
+    def yaxis(self):
+        return self.ax.yaxis
+
+    def update_normals(self, mappable=None):
+        self.set_alpha(self.mappable.get_alpha())
+        self._draw_all()
+
+    def set_alpha(self, alpha):
+        """
+        Set the transparency between 0 (transparent) and 1 (opaque).
+
+        If an array is provided, *alpha* will be set to None to use the
+        transparency values associated with the colormap.
+        """
+        self.alpha = None if isinstance(alpha, np.ndarray) else alpha
+
+    def remove(self):
+        """
+        Remove this colorbar from the figure.
+
+        If the colorbar was created with ``use_gridspec=True`` the previous
+        gridspec is restored.
+        """
+        if hasattr(self.ax, '_colorbar_info'):
+            parents = self.ax._colorbar_info['parents']
+            for a in parents:
+                if self.ax in a._colorbars:
+                    a._colorbars.remove(self.ax)
+        self.ax.remove()
+        self.mappable.callbacks.disconnect(self.mappable.colorbar_cid)
+        self.mappable.colorbar = None
+        self.mappable.colorbar_cid = None
+
+        ax = self.mappable.axes
+        if ax is None:
+            return
+        try:
+            subplotspec = self.ax.get_subplotspec().get_gridspec()._subplot_spec
+        except AttributeError:  # use_gridspec was False
+            pos = ax.get_position(original=True)
+            ax._set_position(pos)
+        else:  # use_gridspec was True
+            ax.set_subplotspec(subplotspec)
+
+    def _get_view(self):
+        ynorm, xnorm = self.colorizer.norm.norms
+        return ynorm.vmin, ynorm.vmax, xnorm.vmin, xnorm.vmax
+
+    def _set_view(self, view):
+        ynorm, xnorm = self.colorizer.norm.norms
+        if (view[0] != ynorm.vmin
+           or view[1] != ynorm.vmax
+           or view[2] != xnorm.vmin
+           or view[3] != xnorm.vmax):
+            with self.colorizer.norm.callbacks.blocked():
+                ynorm.vmin, ynorm.vmax, xnorm.vmin, xnorm.vmax = view
+            self.colorizer.norm._changed()
+
+    def _set_view_from_bbox(self, bbox, direction='in',
+                            mode=None, twinx=False, twiny=False):
+        new_xbound, new_ybound = self.ax._prepare_view_from_bbox(
+            bbox, direction=direction, mode=mode, twinx=twinx, twiny=twiny)
+        ynorm, xnorm = self.colorizer.norm.norms
+        ynorm.vmin, ynorm.vmax = new_ybound
+        xnorm.vmin, xnorm.vmax = new_xbound
+
+    def drag_pan(self, button, key, x, y):
+        points = self.ax._get_pan_points(button, key, x, y)
+        if points is not None:
+            ynorm, xnorm = self.colorizer.norm.norms
+            xnorm.vmin, xnorm.vmax = points[:, 0]
+            ynorm.vmin, ynorm.vmax = points[:, 1]
+
+    def _cbar_cla(self):
+        """Function to clear the interactive colorbar state."""
+        for x in self._interactive_funcs:
+            delattr(self.ax, x)
+        # We now restore the old cla() back and can call it directly
+        del self.ax.cla
+        self.ax.cla()
+
+
+class MultivarColorbar(Sequence):
+    r"""
+    Draw a multivariate colorbar in existing Axes.
+
+    Typically, multivariate colorbars are created using `.Figure.colorbar_multivar`
+    and associated with `.ColorizingArtist`\s (such as an
+    `.AxesImage` generated via `~.axes.Axes.imshow`).
+
+    MultivarColorbar is iterable, and the constituent Colorbar objects can be accessed
+    by index.
+    """
+
+    def __init__(self, axes, mappable=None, **kwargs):
+        """
+        Parameters
+        ----------
+        axes : list of `~matplotlib.axes.Axes`
+            The `~.axes.Axes` instances in which the colorbars are drawn.
+
+        mappable : `.ColorizingArtist`
+            The mappable whose colormap and norm will be used.
+
+        alpha : float
+            The colorbars transparency between 0 (transparent) and 1 (opaque).
+
+        location : None or {'left', 'right', 'top', 'bottom'}
+            Set the multivariate colorbars's location
+
+        Other Parameters
+        ----------------
+        orientation : None or {'vertical', 'horizontal'}
+            If None, use the value determined by *location*. If both
+            *orientation* and *location* are None then defaults to 'vertical'.
+
+        ticklocation : {'auto', 'left', 'right', 'top', 'bottom'}
+            The location of the colorbar ticks. The *ticklocation* must match
+            *orientation*. For example, a horizontal colorbar can only have ticks
+            at the top or the bottom. If 'auto', the ticks will be the same as
+            *location*, so a colorbar to the left will have ticks to the left. If
+            *location* is None, the ticks will be at the bottom for a horizontal
+            colorbar and at the right for a vertical.
+        """
+
+        if isinstance(mappable, mpl.colorizer.Colorizer):
+            mappable = mcolorizer.ColorizingArtist(mappable)
+
+        self.mappable = mappable
+        self.colorizer = mappable.colorizer
+        cmap = self.colorizer.cmap
+        norm = self.colorizer.norm
+        n = cmap.n_variates
+
+        self._colorbars = [Colorbar(axes[i],
+                                    norm=norm.norms[i],
+                                    cmap=cmap[i],
+                                    **kwargs,
+                                    ) for i in range(n)]
+
+        mappable.colorbar = self
+        mappable.colorbar_cid = mappable.callbacks.connect(
+            'changed', self.update_normals)
+
+        self.axes = axes
+
+    def _set_colorbar_info(self, colorbar_info):
+        self._colorbar_info = colorbar_info
+        if colorbar_info is not None:
+            parents = colorbar_info['parents']
+            for a in parents:
+                a._colorbars.append(self)
+
+    def update_normals(self, mappable=None):
+        [c.update_normal() for c in self._colorbars]
+
+    def remove(self):
+        if hasattr(self, '_colorbar_info'):
+            parents = self._colorbar_info['parents']
+            for a in parents:
+                if self in a._colorbars:
+                    a._colorbars.remove(self)
+
+        for ax in self.axes:
+            ax.remove()
+
+        self.mappable.callbacks.disconnect(self.mappable.colorbar_cid)
+        self.mappable.colorbar = None
+        self.mappable.colorbar_cid = None
+
+        try:
+            ax = self.mappable.axes
+            if ax is None:
+                return
+        except AttributeError:
+            return
+        try:
+            subplotspec = self.ax.get_subplotspec().get_gridspec()._subplot_spec
+        except AttributeError:  # use_gridspec was False
+            pos = ax.get_position(original=True)
+            ax._set_position(pos)
+        else:  # use_gridspec was True
+            ax.set_subplotspec(subplotspec)
+
+    def __getitem__(self, index):
+        return self._colorbars[index]
+
+    def __len__(self):
+        return len(self._colorbars)
+
+    def get_tightbbox(self, renderer=None, for_layout_only=False):
+        if for_layout_only and hasattr(self, '_colorbar_info'):
+            # figure out the maximum size of the tight boxes
+            # then multiply that up to the correct size
+            bounds = self.axes[0].get_tightbbox(renderer=renderer,
+                                                for_layout_only=True).bounds
+            x0, y0, width, height = bounds
+            for ax in self.axes[1:]:
+                bb = ax.get_tightbbox(renderer=renderer,
+                                      for_layout_only=True).bounds
+                x0 = min(x0, bb[0])
+                y0 = min(y0, bb[1])
+                width = max(width, bb[2])
+                height = max(height, bb[3])
+            n_major = self._colorbar_info["n_major"]
+            n_minor = self._colorbar_info["n_minor"]
+
+            maj_p = self._colorbar_info["major_pad"]
+            min_p = self._colorbar_info["minor_pad"]
+
+            if self._colorbar_info["orientation"] == 'vertical':
+                # height *= n_major
+                width *= n_minor + min_p * (n_minor - 1)
+                height *= n_major + maj_p * (n_major - 1)
+            else:
+                width *= n_major + maj_p * (n_major - 1)
+                height *= n_minor + min_p * (n_minor - 1)
+
+            bbox = mtransforms.Bbox.from_bounds(x0, y0, width, height)
+        else:
+            # calculate the minimum size of the bbox that
+            # fits the current distribution of colormaps
+            # colorbars in the required grid
+            bbox = self.axes[0].get_tightbbox(renderer=renderer,
+                                              for_layout_only=for_layout_only)
+            for ax in self.axes[1:]:
+                bb = ax.get_tightbbox(renderer=renderer,
+                                      for_layout_only=for_layout_only)
+                if bb.x0 < bbox.x0:
+                    bbox.x0 = bb.x0
+                if bb.y0 < bbox.y0:
+                    bbox.y0 = bb.y0
+                if bb.x1 > bbox.x1:
+                    bbox.x1 = bb.x1
+                if bb.y1 > bbox.y1:
+                    bbox.y1 = bb.y1
+        return bbox
+
+    def _get_tight_packing_inside(self, bbox):
+        """
+        Positions the colorbars in a grid contained in bbox
+
+        This function relates to get_tightbbox(for_layout_only=True)
+        It requires that self._colorbar_info exists
+
+        It will position colorbars so that the bboxes
+        of each cmponent spans the assigned bbox
+
+        The procedure for vertical colorbars is as follows:
+        1. Calculate the height of each colorbar, based on n_minor and the padding
+        2. From the height calculate the width using the aspect
+        3. With both the width and height known, position the colorbars equispaced.
+
+        Horizontal colorbars follow the same procedure but the height and width
+        are swapped.
+        """
+        x_b, y_b, width_b, height_b = bbox.bounds  # b for box
+        if hasattr(self, '_colorbar_info'):
+            n_major = self._colorbar_info["n_major"]
+            n_minor = self._colorbar_info["n_minor"]
+            aspect = self._colorbar_info["aspect"]
+            major_pad = self._colorbar_info["major_pad"]
+            # pad_minor = 0.6
+            if self._colorbar_info["orientation"] == 'vertical':
+                if n_major > 1:
+                    bar_height = height_b * (1 - major_pad) / n_major
+                    y_spacing = (height_b - bar_height * n_major)/(n_major-1)
+                else:
+                    bar_height = height_b
+                    y_spacing = 0
+                bar_width = bar_height / aspect
+                if n_minor > 1:
+                    x_spacing = (width_b - bar_width * n_minor)/(n_minor)
+                else:
+                    x_spacing = 0
+            else:
+                if n_major > 1:
+                    bar_width = width_b * (1 - major_pad) / n_major
+                    x_spacing = (width_b - bar_width * n_major)/(n_major-1)
+                else:
+                    bar_width = width_b
+                    x_spacing = 0
+                bar_height = bar_width / aspect
+                if n_minor > 1:
+                    y_spacing = (height_b - bar_height * n_minor)/(n_minor)
+                else:
+                    y_spacing = 0
+
+            x_step = bar_width + x_spacing
+            y_step = bar_height + y_spacing
+            bboxes = []
+            if self._colorbar_info["orientation"] == 'vertical':
+                for i in range(n_minor):
+                    for j in range(n_major):
+                        bboxes.append(mtransforms.Bbox([[x_b + i * x_step,
+                                                         y_b + j * y_step],
+                                                        [x_b + i * x_step + bar_width,
+                                                         y_b + j * y_step + bar_height],
+                                                        ]))
+            else:
+                for i in range(n_minor):
+                    for j in range(n_major):
+                        bboxes.append(mtransforms.Bbox([[x_b + j * x_step,
+                                                         y_b + i * y_step],
+                                                        [x_b + j * x_step + bar_width,
+                                                         y_b + i * y_step + bar_height],
+                                                        ]))
+            return bboxes
+        else:
+            raise ValueError("_set_tight_packing cannot be called "
+                             "unless the MultivarColorbar was created "
+                             "by fig.multicolorbar(mappable).")
+
+    def _get_original_position(self):
+        # comparable to axes.get_position(original=True)
+        bbox = self.axes[0].get_position(original=True)
+        for ax in self.axes[1:]:
+            bb = ax.get_position(original=True)
+            if bb.x0 < bbox.x0:
+                bbox.x0 = bb.x0
+            if bb.y0 < bbox.y0:
+                bbox.y0 = bb.y0
+            if bb.x1 > bbox.x1:
+                bbox.x1 = bb.x1
+            if bb.y1 > bbox.y1:
+                bbox.y1 = bb.y1
+        return bbox
+
+    @staticmethod
+    def _subdivide_bbox(bbox,
+                        n_major,
+                        n_minor,
+                        orientation,
+                        major_pad=0.1,
+                        minor_pad=0.6,
+                        ):
+        major_width = (1-major_pad) / n_major
+        if n_major > 1:
+            major_space = major_pad / (n_major - 1)
+        else:
+            major_space = 0
+        major_split = np.empty(2 * (n_major - 1))
+        major_split[0::2] = major_width
+        major_split[1::2] = major_space
+        major_split = np.cumsum(major_split)
+
+        minor_width = (1-minor_pad) / n_minor
+        if n_minor > 1:
+            minor_space = minor_pad / (n_minor - 1)
+        else:
+            minor_space = 0
+        minor_split = np.empty(2 * (n_minor - 1))
+        minor_split[0::2] = minor_width
+        minor_split[1::2] = minor_space
+        minor_split = np.cumsum(minor_split)
+
+        # make the colorbar bboxes
+        sub_bboxs = []
+        v = orientation == "vertical"
+        bboxs = bbox.splitx(*minor_split) if v else bbox.splity(*minor_split)[::-1]
+        for i, bboxs_i in enumerate(bboxs):
+            if i % 2 == 1:
+                continue
+            bboxs_j = bboxs_i.splity(*major_split
+                                     )[::-1] if v else bboxs_i.splitx(*major_split)
+            for j, sub_bbox in enumerate(bboxs_j):
+                if j % 2 == 1:
+                    continue
+                sub_bboxs.append(sub_bbox)
+        return sub_bboxs
+
+
 def _normalize_location_orientation(location, orientation):
     if location is None:
         location = _get_ticklocation_from_orientation(orientation)
+    loc_settings = _normalize_location(location)
+    loc_settings["orientation"] = _get_orientation_from_location(location)
+    if orientation is not None and orientation != loc_settings["orientation"]:
+        # Allow the user to pass both if they are consistent.
+        raise TypeError("location and orientation are mutually exclusive")
+    return loc_settings
+
+
+def _normalize_location(location):
+    if location is None:
+        location = 'right'
     loc_settings = _api.getitem_checked({
         "left":   {"location": "left", "anchor": (1.0, 0.5),
                    "panchor": (0.0, 0.5), "pad": 0.10},
@@ -1370,10 +2011,6 @@ def _normalize_location_orientation(location, orientation):
         "bottom": {"location": "bottom", "anchor": (0.5, 1.0),
                    "panchor": (0.5, 0.0), "pad": 0.15},
     }, location=location)
-    loc_settings["orientation"] = _get_orientation_from_location(location)
-    if orientation is not None and orientation != loc_settings["orientation"]:
-        # Allow the user to pass both if they are consistent.
-        raise TypeError("location and orientation are mutually exclusive")
     return loc_settings
 
 
@@ -1389,39 +2026,18 @@ def _get_ticklocation_from_orientation(orientation):
         orientation=orientation)
 
 
-@_docstring.interpd
-def make_axes(parents, location=None, orientation=None, fraction=0.15,
-              shrink=1.0, aspect=20, **kwargs):
-    """
-    Create an `~.axes.Axes` suitable for a colorbar.
+def _get_bivar_ticklocations_from_location(location):
+    loc_0 = _api.getitem_checked(
+        {None: "right", "left": "left", "right": "right",
+         "top": "left", "bottom": "left"}, location=location)
+    loc_1 = _api.getitem_checked(
+        {None: "bottom", "left": "bottom", "right": "bottom",
+         "top": "top", "bottom": "bottom"}, location=location)
+    return (loc_0, loc_1)
 
-    The Axes is placed in the figure of the *parents* Axes, by resizing and
-    repositioning *parents*.
 
-    Parameters
-    ----------
-    parents : `~matplotlib.axes.Axes` or iterable or `numpy.ndarray` of `~.axes.Axes`
-        The Axes to use as parents for placing the colorbar.
-    %(_make_axes_kw_doc)s
-
-    Returns
-    -------
-    cax : `~matplotlib.axes.Axes`
-        The child Axes.
-    kwargs : dict
-        The reduced keyword dictionary to be passed when creating the colorbar
-        instance.
-    """
-    loc_settings = _normalize_location_orientation(location, orientation)
-    # put appropriate values into the kwargs dict for passing back to
-    # the Colorbar class
-    kwargs['orientation'] = loc_settings['orientation']
-    location = kwargs['ticklocation'] = loc_settings['location']
-
-    anchor = kwargs.pop('anchor', loc_settings['anchor'])
-    panchor = kwargs.pop('panchor', loc_settings['panchor'])
-    aspect0 = aspect
-    # turn parents into a list if it is not already.  Note we cannot
+def _normalize_parents(parents):
+    # Turn parents into a list if it is not already. Note we cannot
     # use .flatten or .ravel as these copy the references rather than
     # reuse them, leading to a memory leak
     if isinstance(parents, np.ndarray):
@@ -1433,18 +2049,20 @@ def make_axes(parents, location=None, orientation=None, fraction=0.15,
 
     fig = parents[0].get_figure()
 
-    pad0 = 0.05 if fig.get_constrained_layout() else loc_settings['pad']
-    pad = kwargs.pop('pad', pad0)
-
     if not all(fig is ax.get_figure() for ax in parents):
         raise ValueError('Unable to create a colorbar Axes as not all '
                          'parents share the same figure.')
+    return parents, fig
+
 
+def _get_bbox_shrink_parents(parents, location, fraction,
+                             pad, shrink, anchor, panchor):
+    """Shrink parents and get the bbox for the colorbar."""
     # take a bounding box around all of the given Axes
-    parents_bbox = mtransforms.Bbox.union(
+    pb = mtransforms.Bbox.union(
         [ax.get_position(original=True).frozen() for ax in parents])
 
-    pb = parents_bbox
+    # calculate the new bounding boxes
     if location in ('left', 'right'):
         if location == 'left':
             pbcb, _, pb1 = pb.splitx(fraction, fraction + pad)
@@ -1458,64 +2076,164 @@ def make_axes(parents, location=None, orientation=None, fraction=0.15,
             pb1, _, pbcb = pb.splity(1 - fraction - pad, 1 - fraction)
         pbcb = pbcb.shrunk(shrink, 1.0).anchored(anchor, pbcb)
 
-        # define the aspect ratio in terms of y's per x rather than x's per y
-        aspect = 1.0 / aspect
-
     # define a transform which takes us from old axes coordinates to
     # new axes coordinates
-    shrinking_trans = mtransforms.BboxTransform(parents_bbox, pb1)
+    shrinking_trans = mtransforms.BboxTransform(pb, pb1)
 
-    # transform each of the Axes in parents using the new transform
     for ax in parents:
         new_posn = shrinking_trans.transform(ax.get_position(original=True))
         new_posn = mtransforms.Bbox(new_posn)
         ax._set_position(new_posn)
         if panchor is not False:
             ax.set_anchor(panchor)
+    return pbcb
 
-    cax = fig.add_axes(pbcb, label="<colorbar>")
-    for a in parents:
-        a._colorbars.append(cax)  # tell the parent it has a colorbar
-    cax._colorbar_info = dict(
+
+def _make_axes_get_pbcb(loc_settings, fraction, shrink, aspect, kwargs, parents):
+
+    location = loc_settings['location']
+    kwargs['location'] = location
+
+    anchor = kwargs.pop('anchor', loc_settings['anchor'])
+    panchor = kwargs.pop('panchor', loc_settings['panchor'])
+    parents, fig = _normalize_parents(parents)
+
+    pad0 = 0.05 if fig.get_constrained_layout() else loc_settings['pad']
+    pad = kwargs.pop('pad', pad0)
+
+    # shrink the parents and get the bbox for the colorbar
+    pbcb = _get_bbox_shrink_parents(parents, location, fraction,
+                                    pad, shrink, anchor, panchor)
+    colorbar_info = dict(
         parents=parents,
         location=location,
         shrink=shrink,
         anchor=anchor,
         panchor=panchor,
         fraction=fraction,
-        aspect=aspect0,
+        aspect=aspect,
         pad=pad)
-    # and we need to set the aspect ratio by hand...
-    cax.set_anchor(anchor)
+
+    return fig, pbcb, colorbar_info
+
+
+def _make_axes_helper(loc_settings, fraction, shrink, aspect, kwargs, parents):
+    """
+    Help function for `make_axes` and `make_bivar_axes`.
+
+    `make_axes` and `make_bivar_axes` are identical
+    except for the fact that `make_axes` also deals with:
+    1. the aspect
+    2. orientation.
+
+    `make_bivar_axes` on the other hand has no concept of orientation
+    and the aspect is handled by the `BivarColormap` instance, not during
+    creation of the axes.
+    """
+    fig, pbcb, colorbar_info = _make_axes_get_pbcb(loc_settings, fraction, shrink,
+                                                   aspect, kwargs, parents)
+    cax = fig.add_axes(pbcb, label="<colorbar>")
+    for a in colorbar_info["parents"]:
+        a._colorbars.append(cax)  # tell the parent it has a colorbar
+
+    cax._colorbar_info = colorbar_info
+    cax.set_anchor(colorbar_info["anchor"])
+    return cax
+
+
+@_docstring.interpd
+def make_axes(parents, location=None, orientation=None, fraction=0.15,
+              shrink=1.0, aspect=20, **kwargs):
+    """
+    Create an `~.axes.Axes` suitable for a colorbar.
+
+    The Axes is placed in the figure of the *parents* Axes, by resizing and
+    repositioning *parents*.
+
+    Parameters
+    ----------
+    parents : `~matplotlib.axes.Axes` or iterable or `numpy.ndarray` of `~.axes.Axes`
+        The Axes to use as parents for placing the colorbar.
+    %(_make_axes_kw_doc)s
+
+    Returns
+    -------
+    cax : `~matplotlib.axes.Axes`
+        The child Axes.
+    kwargs : dict
+        The reduced keyword dictionary to be passed when creating the colorbar
+        instance.
+    """
+    loc_settings = _normalize_location_orientation(location, orientation)
+    # put appropriate values into the kwargs dict for passing back to
+    # the Colorbar class
+    kwargs['orientation'] = loc_settings['orientation']
+    kwargs['ticklocation'] = loc_settings['location']
+
+    cax = _make_axes_helper(loc_settings, fraction, shrink, aspect, kwargs, parents)
+    cax._colorbar_info["type"] = 'Colorbar'
+    if loc_settings["location"] in ('top', 'bottom'):
+        aspect = 1.0 / aspect
     cax.set_box_aspect(aspect)
     cax.set_aspect('auto')
+    # and we need to set the aspect ratio by hand...
 
     return cax, kwargs
 
 
 @_docstring.interpd
-def make_axes_gridspec(parent, *, location=None, orientation=None,
-                       fraction=0.15, shrink=1.0, aspect=20, **kwargs):
+def make_bivar_axes(parents, location=None, fraction=0.15,
+                    shrink=1.0, aspect=1.0, **kwargs):
     """
-    Create an `~.axes.Axes` suitable for a colorbar.
+    Create an `~.axes.Axes` suitable for a bivariate colorbar.
 
-    The Axes is placed in the figure of the *parent* Axes, by resizing and
-    repositioning *parent*.
+    The Axes is placed in the figure of the *parents* Axes, by resizing and
+    repositioning *parents*.
 
-    This function is similar to `.make_axes` and mostly compatible with it.
-    Primary differences are
+    Parameters
+    ----------
+    parents : `~matplotlib.axes.Axes` or iterable or `numpy.ndarray` of `~.axes.Axes`
+        The Axes to use as parents for placing the colorbar.
+    %(_make_bivar_axes_kw_doc)s
 
-    - `.make_axes_gridspec` requires the *parent* to have a subplotspec.
-    - `.make_axes` positions the Axes in figure coordinates;
-      `.make_axes_gridspec` positions it using a subplotspec.
-    - `.make_axes` updates the position of the parent.  `.make_axes_gridspec`
-      replaces the parent gridspec with a new one.
+    Returns
+    -------
+    cax : `~matplotlib.axes.Axes`
+        The child Axes.
+    kwargs : dict
+        The reduced keyword dictionary to be passed when creating the colorbar
+        instance.
+    """
+    loc_settings = _normalize_location(location)
+    cax = _make_axes_helper(loc_settings, fraction, shrink, aspect, kwargs, parents)
+    cax._colorbar_info["type"] = 'BivarColorbar'
+    # need to add aspect to kwargs so it propagates to the BivarColorbar
+    kwargs["aspect"] = aspect
+
+    return cax, kwargs
+
+
+@_docstring.interpd
+def make_multivar_axes(parents, n_variates, n_major, location=None, orientation=None,
+                       fraction=0.15, shrink=1.0, aspect=20, major_pad=0.2,
+                       minor_pad=0.6, **kwargs):
+    """
+    Create an `~.axes.Axes` suitable for a mulitvariate colorbar.
+
+    The Axes is placed in the figure of the *parents* Axes, by resizing and
+    repositioning *parents*.
 
     Parameters
     ----------
-    parent : `~matplotlib.axes.Axes`
-        The Axes to use as parent for placing the colorbar.
-    %(_make_axes_kw_doc)s
+    parents : `~matplotlib.axes.Axes` or iterable or `numpy.ndarray` of `~.axes.Axes`
+        The Axes to use as parents for placing the colorbar.
+
+    n_variates : int
+        The number of colorbars to be made
+
+    n_major : int
+        Number of colorbars along the long axis of the colorbars
+    %(_make_multivar_axes_kw_doc)s
 
     Returns
     -------
@@ -1525,12 +2243,88 @@ def make_axes_gridspec(parent, *, location=None, orientation=None,
         The reduced keyword dictionary to be passed when creating the colorbar
         instance.
     """
+    if n_major == -1:
+        n_major = n_variates
+        n_minor = 1
+    else:
+        if n_major == 0:
+            raise ValueError("n_major cannot be zero")
+        n_minor = n_variates // n_major
+        if n_major * n_minor < n_variates:
+            n_minor += 1
+
+    aspect = aspect / n_major
+    fraction = fraction * n_minor
 
     loc_settings = _normalize_location_orientation(location, orientation)
-    kwargs['orientation'] = loc_settings['orientation']
-    location = kwargs['ticklocation'] = loc_settings['location']
 
-    aspect0 = aspect
+    # put appropriate values into the kwargs dict for passing back to
+    # the Colorbar class
+    orientation = loc_settings['orientation']
+    kwargs['orientation'] = orientation
+    kwargs['ticklocation'] = loc_settings['location']
+
+    # get the shape of the grid of new colorbars
+
+    if n_minor > 1:
+        location = loc_settings["location"]
+        if location == 'left':
+            loc_settings["anchor"] = (1, 0.5)
+        elif location == 'right':
+            loc_settings["anchor"] = (0, 0.5)
+        elif location == 'top':
+            loc_settings["anchor"] = (0.5, 0)
+        else:
+            loc_settings["anchor"] = (0.5, 1)
+
+    fig, pbcb, colorbar_info = _make_axes_get_pbcb(loc_settings, fraction, shrink,
+                                                   aspect, kwargs, parents)
+    colorbar_info["type"] = 'MultivarColorbar'
+
+    # split pbcb into the required parts
+    sub_bboxes = MultivarColorbar._subdivide_bbox(pbcb,
+                                                  n_major,
+                                                  n_minor,
+                                                  orientation,
+                                                  major_pad=major_pad,
+                                                  minor_pad=minor_pad,
+                                                  )[:n_variates]
+    caxes = [fig.add_axes(sub_bbox, label="<colorbar>")
+             for sub_bbox in sub_bboxes]
+
+    # adjust aspect
+    if loc_settings["location"] in ('top', 'bottom'):
+        aspect = 1.0 / aspect
+
+    for cax in caxes:
+        cax.set_anchor(colorbar_info["anchor"])
+        cax.set_box_aspect(aspect)
+        cax.set_aspect('auto')
+
+    colorbar_info["n_major"] = n_major
+    colorbar_info["n_minor"] = n_minor
+    colorbar_info["major_pad"] = major_pad
+    colorbar_info["minor_pad"] = minor_pad
+    colorbar_info["orientation"] = orientation
+    return caxes, kwargs, colorbar_info
+
+
+def _make_axes_gridspec_helper(loc_settings, fraction, shrink, aspect, kwargs, parent):
+    """
+    Help function for `make_axes_gridspec` and `make_bivar_axes_gridspec`.
+
+    `make_axes_gridspec` and `make_bivar_axes_gridspec` are identical
+    except for the fact that `make_axes_gridspec` also deals with:
+    1. the aspect
+    2. orientation.
+
+    `make_bivar_axes_gridspec` on the other hand has no concept of orientation
+    and the aspect is handled by the `BivarColormap` instance, not during creation
+    of the axes.
+    """
+    location = loc_settings['location']
+    kwargs['location'] = location
+
     anchor = kwargs.pop('anchor', loc_settings['anchor'])
     panchor = kwargs.pop('panchor', loc_settings['panchor'])
     pad = kwargs.pop('pad', loc_settings["pad"])
@@ -1560,7 +2354,6 @@ def make_axes_gridspec(parent, *, location=None, orientation=None,
             gs.set_height_ratios([1 - fraction - pad, fraction])
             ss_main = gs[0, :]
             ss_cb = gs[1, 1]
-        aspect = 1 / aspect
 
     parent.set_subplotspec(ss_main)
     if panchor is not False:
@@ -1570,8 +2363,7 @@ def make_axes_gridspec(parent, *, location=None, orientation=None,
     cax = fig.add_subplot(ss_cb, label="<colorbar>")
     parent._colorbars.append(cax)  # tell the parent it has a colorbar
     cax.set_anchor(anchor)
-    cax.set_box_aspect(aspect)
-    cax.set_aspect('auto')
+
     cax._colorbar_info = dict(
         location=location,
         parents=[parent],
@@ -1579,7 +2371,98 @@ def make_axes_gridspec(parent, *, location=None, orientation=None,
         anchor=anchor,
         panchor=panchor,
         fraction=fraction,
-        aspect=aspect0,
+        aspect=aspect,
         pad=pad)
 
+    return cax
+
+
+@_docstring.interpd
+def make_axes_gridspec(parent, *, location=None, orientation=None,
+                       fraction=0.15, shrink=1.0, aspect=20, **kwargs):
+    """
+    Create an `~.axes.Axes` suitable for a colorbar.
+
+    The Axes is placed in the figure of the *parent* Axes, by resizing and
+    repositioning *parent*.
+
+    This function is similar to `.make_axes` and mostly compatible with it.
+    Primary differences are
+
+    - `.make_axes_gridspec` requires the *parent* to have a subplotspec.
+    - `.make_axes` positions the Axes in figure coordinates;
+      `.make_axes_gridspec` positions it using a subplotspec.
+    - `.make_axes` updates the position of the parent.  `.make_axes_gridspec`
+      replaces the parent gridspec with a new one.
+
+    Parameters
+    ----------
+    parent : `~matplotlib.axes.Axes`
+        The Axes to use as parent for placing the colorbar.
+    %(_make_axes_kw_doc)s
+
+    Returns
+    -------
+    cax : `~matplotlib.axes.Axes`
+        The child Axes.
+    kwargs : dict
+        The reduced keyword dictionary to be passed when creating the colorbar
+        instance.
+    """
+    loc_settings = _normalize_location_orientation(location, orientation)
+    kwargs['orientation'] = loc_settings['orientation']
+    kwargs['ticklocation'] = loc_settings['location']
+
+    cax = _make_axes_gridspec_helper(loc_settings, fraction, shrink,
+                                     aspect, kwargs, parent)
+    cax._colorbar_info["type"] = 'Colorbar'
+    if loc_settings["location"] in ('top', 'bottom'):
+        aspect = 1 / aspect
+    cax.set_box_aspect(aspect)
+    cax.set_aspect('auto')
+
+    return cax, kwargs
+
+
+@_docstring.interpd
+def make_bivar_axes_gridspec(parent, *, location=None,
+                             fraction=0.15, shrink=1.0, aspect=1.0, **kwargs):
+    """
+    Create an `~.axes.Axes` suitable for a bivariate colorbar.
+
+    The Axes is placed in the figure of the *parent* Axes, by resizing and
+    repositioning *parent*.
+
+    This function is similar to `.make_axes` and mostly compatible with it.
+    Primary differences are
+
+    - `.make_axes_gridspec` requires the *parent* to have a subplotspec.
+    - `.make_axes` positions the Axes in figure coordinates;
+      `.make_axes_gridspec` positions it using a subplotspec.
+    - `.make_axes` updates the position of the parent.  `.make_axes_gridspec`
+      replaces the parent gridspec with a new one.
+
+    Parameters
+    ----------
+    parent : `~matplotlib.axes.Axes`
+        The Axes to use as parent for placing the colorbar.
+    %(_make_bivar_axes_kw_doc)s
+
+    Returns
+    -------
+    cax : `~matplotlib.axes.Axes`
+        The child Axes.
+    kwargs : dict
+        The reduced keyword dictionary to be passed when creating the colorbar
+        instance.
+    """
+    loc_settings = _normalize_location(location)
+
+    cax = _make_axes_gridspec_helper(loc_settings, fraction, shrink,
+                                     aspect, kwargs, parent)
+    cax._colorbar_info["type"] = 'BivarColorbar'
+
+    # need to add aspect to kwargs so it propagates to the BivarColorbar
+    kwargs["aspect"] = aspect
+
     return cax, kwargs
diff --git a/lib/matplotlib/colorbar.pyi b/lib/matplotlib/colorbar.pyi
index d1401b238793..b64d326ef5b9 100644
--- a/lib/matplotlib/colorbar.pyi
+++ b/lib/matplotlib/colorbar.pyi
@@ -1,5 +1,5 @@
 import matplotlib.spines as mspines
-from matplotlib import cm, collections, colors, contour, colorizer
+from matplotlib import cm, collections, colors, contour, colorizer as mcolorizer
 from matplotlib.axes import Axes
 from matplotlib.axis import Axis
 from matplotlib.backend_bases import RendererBase
@@ -21,7 +21,7 @@ class _ColorbarSpine(mspines.Spine):
 
 class Colorbar:
     n_rasterize: int
-    mappable: cm.ScalarMappable | colorizer.ColorizingArtist
+    mappable: cm.ScalarMappable | mcolorizer.ColorizingArtist
     ax: Axes
     alpha: float | None
     cmap: colors.Colormap
@@ -43,7 +43,7 @@ class Colorbar:
     def __init__(
         self,
         ax: Axes,
-        mappable: cm.ScalarMappable | colorizer.ColorizingArtist | None = ...,
+        mappable: cm.ScalarMappable | mcolorizer.ColorizingArtist | None = ...,
         *,
         cmap: str | colors.Colormap | None = ...,
         norm: colors.Normalize | None = ...,
@@ -117,6 +117,59 @@ class Colorbar:
 
 ColorbarBase = Colorbar
 
+class BivarColorbar:
+    n_rasterize: int
+    mappable: mcolorizer.ColorizingArtist
+    ax: Axes
+    alpha: float | None
+    colorizer: mcolorizer.Colorizer
+    ticklocations: tuple[Literal["auto", "left", "right"], Literal["auto", "top", "bottom"]]
+    def __init__(
+        self,
+        ax: Axes,
+        mappable: mcolorizer.ColorizingArtist | mcolorizer.Colorizer,
+        *,
+        alpha: float | None = ...,
+        location: Literal["left", "right", "top", "bottom"] | None = ...,
+        ticklocations: tuple[Literal["auto", "left", "right"], Literal["auto", "top", "bottom"]] = ...,
+        aspect: float = ...,
+    ) -> None: ...
+    @property
+    def aspect(self) -> float: ...
+    @aspect.setter
+    def aspect(self, aspect: float) -> None: ...
+    def set_xlabel(self, label: str) -> None: ...
+    def set_ylabel(self, label: str) -> None: ...
+    @property
+    def xaxis(self) -> Axis: ...
+    @property
+    def yaxis(self) -> Axis: ...
+    def update_normals(self, mappable: mcolorizer.ColorizingArtist | None = ...) -> None: ...
+    def set_alpha(self, alpha: float | None) -> None: ...
+    def remove(self) -> None: ...
+    def drag_pan(self, button: Any, key: Any, x: float, y: float) -> None: ...
+
+class MultivarColorbar(Sequence[Colorbar]):
+    mappable: mcolorizer.ColorizingArtist
+    colorizer: mcolorizer.Colorizer
+    axes: Sequence[Axes]
+    _colorbars: list[Colorbar]
+    def __init__(
+        self,
+        axes: Sequence[Axes],
+        mappable: mcolorizer.ColorizingArtist | mcolorizer.Colorizer | None = ...,
+        **kwargs: Any,
+    ) -> None: ...
+
+    def update_normals(self, mappable: mcolorizer.ColorizingArtist | None = ...) -> None: ...
+    def remove(self) -> None: ...
+    @overload
+    def __getitem__(self, index: int, /) -> Colorbar: ...
+    @overload
+    def __getitem__(self, index: slice[int | None, int | None, int | None], /) -> Sequence[Colorbar]: ...
+    def __len__(self) -> int: ...
+    def get_tightbbox(self, renderer: RendererBase | None = ..., for_layout_only: bool = ...) -> Bbox: ...
+
 def make_axes(
     parents: Axes | list[Axes] | np.ndarray,
     location: Literal["left", "right", "top", "bottom"] | None = ...,
@@ -126,6 +179,27 @@ def make_axes(
     aspect: float = ...,
     **kwargs
 ) -> tuple[Axes, dict[str, Any]]: ...
+def make_bivar_axes(
+    parents: Axes | list[Axes] | np.ndarray,
+    location: Literal["left", "right", "top", "bottom"] | None = ...,
+    fraction: float = ...,
+    shrink: float = ...,
+    aspect: float = ...,
+    **kwargs
+) -> tuple[Axes, dict[str, Any]]: ...
+def make_multivar_axes(
+    parents: Axes | list[Axes] | np.ndarray,
+    n_variates: int,
+    n_major: int,
+    location: Literal["left", "right", "top", "bottom"] | None = ...,
+    orientation: Literal["vertical", "horizontal"] | None = ...,
+    fraction: float = ...,
+    shrink: float = ...,
+    aspect: float = ...,
+    major_pad: float = ...,
+    minor_pad: float = ...,
+    **kwargs
+) -> tuple[Axes, dict[str, Any], dict[str, Any]]: ...
 def make_axes_gridspec(
     parent: Axes,
     *,
@@ -136,3 +210,12 @@ def make_axes_gridspec(
     aspect: float = ...,
     **kwargs
 ) -> tuple[Axes, dict[str, Any]]: ...
+def make_bivar_axes_gridspec(
+    parent: Axes | list[Axes] | np.ndarray,
+    *,
+    location: Literal["left", "right", "top", "bottom"] | None = ...,
+    fraction: float = ...,
+    shrink: float = ...,
+    aspect: float = ...,
+    **kwargs
+) -> tuple[Axes, dict[str, Any]]: ...
diff --git a/lib/matplotlib/colorizer.py b/lib/matplotlib/colorizer.py
index 99cbcf157db4..e2fb7e3b31d0 100644
--- a/lib/matplotlib/colorizer.py
+++ b/lib/matplotlib/colorizer.py
@@ -805,7 +805,8 @@ def _ensure_norm(norm, n_components=1):
             _api.check_isinstance((colors.MultiNorm, None, tuple), norm=norm)
         if norm is None:
             norm = colors.MultiNorm(['linear']*n_components)
-        else:  # iterable, i.e. multiple strings or Normalize objects
+        elif not isinstance(norm, colors.MultiNorm):
+            # iterable, i.e. multiple strings or Normalize objects
             norm = colors.MultiNorm(norm)
         if isinstance(norm, colors.MultiNorm) and norm.n_components == n_components:
             return norm
diff --git a/lib/matplotlib/colors.py b/lib/matplotlib/colors.py
index 51b440a3a99d..5d9561b95d62 100644
--- a/lib/matplotlib/colors.py
+++ b/lib/matplotlib/colors.py
@@ -3573,11 +3573,16 @@ def autoscale_None(self, A):
             - If structured array, must have `n_components` fields. Each field
               is used for the limits of one constituent norm.
         """
+        changed = False
         with self.callbacks.blocked():
             A = self._iterable_components_in_data(A, self.n_components)
             for n, a in zip(self.norms, A):
+                vmin, vmax = n.vmin, n.vmax
                 n.autoscale_None(a)
-        self._changed()
+                if vmin != n.vmin or vmax != n.vmax:
+                    changed = True
+        if changed:
+            self._changed()
 
     def scaled(self):
         """Return whether both *vmin* and *vmax* are set on all constituent norms."""
diff --git a/lib/matplotlib/figure.py b/lib/matplotlib/figure.py
index 9920f6d908b3..d3aef0c545ff 100644
--- a/lib/matplotlib/figure.py
+++ b/lib/matplotlib/figure.py
@@ -46,6 +46,8 @@
 import matplotlib._api as _api
 import matplotlib.cbook as cbook
 import matplotlib.colorbar as cbar
+import matplotlib.colors as mcolors
+import matplotlib.colorizer as mcolorizer
 import matplotlib.image as mimage
 
 from matplotlib.axes import Axes
@@ -1353,6 +1355,13 @@ def colorbar(
         therefore, this workaround is not used by default (see issue #1188).
 
         """
+        if isinstance(mappable.cmap, mcolors.BivarColormap):
+            raise ValueError("`Figure.colorbar` can only be together with a "
+                             "scalar colormap, please use `Figure.colorbar_bivar` "
+                             "when working with a bivariate colormap")
+        if isinstance(mappable.cmap, mcolors.MultivarColormap):
+            raise ValueError("colorbar can only be together with a"
+                             "scalar colormap")
 
         if ax is None:
             ax = getattr(mappable, "axes", None)
@@ -1399,6 +1408,188 @@ def colorbar(
         cax.get_figure(root=False).stale = True
         return cb
 
+    @_docstring.interpd
+    def colorbar_bivar(
+            self, mappable, *, cax=None, ax=None, use_gridspec=True, **kwargs):
+        """
+        Add a bivariate colorbar to a plot.
+
+        Parameters
+        ----------
+        mappable
+            The `matplotlib.colorizer.ColorizingArtist` (i.e., `.AxesImage`
+            etc.) described by this bivariate colorbar.
+            This argument is mandatory for the `.Figure.colorbar_bivar` method
+            but optional for the`.pyplot.colorbar_bivar` function, which sets
+            the default to the current image.
+
+        cax : `~matplotlib.axes.Axes`, optional
+            Axes into which the colorbar will be drawn.  If `None`, then a new
+            Axes is created and the space for it will be stolen from the Axes(s)
+            specified in *ax*.
+
+        ax : `~matplotlib.axes.Axes` or iterable or `numpy.ndarray` of Axes, optional
+            The one or more parent Axes from which space for a new colorbar Axes
+            will be stolen. This parameter is only used if *cax* is not set.
+
+            Defaults to the Axes that contains the mappable used to create the
+            colorbar.
+
+        use_gridspec : bool, optional
+            If *cax* is ``None``, a new *cax* is created as an instance of
+            Axes.  If *ax* is positioned with a subplotspec and *use_gridspec*
+            is ``True``, then *cax* is also positioned with a subplotspec.
+
+        Returns
+        -------
+        bivariate_colorbar : `~matplotlib.colorbar.BivarColorbar`
+
+        Other Parameters
+        ----------------
+        %(_make_bivar_axes_kw_doc)s
+        %(_bivar_colormap_kw_doc)s
+
+        """
+
+        if isinstance(mappable, mpl.colorizer.Colorizer):
+            mappable = mcolorizer.ColorizingArtist(mappable)
+        if not isinstance(mappable.colorizer.cmap, mcolors.BivarColormap):
+            raise ValueError("A bivariate colorbar can only be used together with a "
+                             f"bivariate colormap, not {type(mappable.colorizer.cmap)}")
+        if ax is None:
+            ax = getattr(mappable, "axes", None)
+
+        if cax is None:
+            if ax is None:
+                raise ValueError(
+                    'Unable to determine Axes to steal space for Colorbar. '
+                    'Either provide the *cax* argument to use as the Axes for '
+                    'the Colorbar, provide the *ax* argument to steal space '
+                    'from it, or add *mappable* to an Axes.')
+            fig = (  # Figure of first Axes; logic copied from make_axes.
+                [*ax.flat] if isinstance(ax, np.ndarray)
+                else [*ax] if np.iterable(ax)
+                else [ax])[0].get_figure(root=False)
+            current_ax = fig.gca()
+            if (fig.get_layout_engine() is not None and
+                    not fig.get_layout_engine().colorbar_gridspec):
+                use_gridspec = False
+            if (use_gridspec
+                    and isinstance(ax, mpl.axes._base._AxesBase)
+                    and ax.get_subplotspec()):
+                cax, kwargs = cbar.make_bivar_axes_gridspec(ax, **kwargs)
+            else:
+                cax, kwargs = cbar.make_bivar_axes(ax, **kwargs)
+            # make_axes calls add_{axes,subplot} which changes gca; undo that.
+            fig.sca(current_ax)
+            cax.grid(visible=False, which='both', axis='both')
+
+        if (hasattr(mappable, "get_figure") and
+                (mappable_host_fig := mappable.get_figure(root=True)) is not None):
+            # Warn in case of mismatch
+            if mappable_host_fig is not self._root_figure:
+                _api.warn_external(
+                        f'Adding colorbar to a different Figure '
+                        f'{repr(mappable_host_fig)} than '
+                        f'{repr(self._root_figure)} which '
+                        f'fig.colorbar is called on.')
+        NON_COLORBAR_KEYS = [  # remove kws that cannot be passed to Colorbar
+            'fraction', 'pad', 'shrink', 'anchor', 'panchor']
+        cb = cbar.BivarColorbar(cax, mappable, **{
+              k: v for k, v in kwargs.items() if k not in NON_COLORBAR_KEYS})
+        cax.get_figure(root=False).stale = True
+        return cb
+
+    @_docstring.interpd
+    def colorbar_multivar(
+            self, mappable, *, caxes=None, ax=None,
+            n_major=-1, **kwargs):
+        """
+        Add a bivariate colorbar to a plot.
+
+        Parameters
+        ----------
+        mappable
+            The `matplotlib.colorizer.ColorizingArtist` (i.e., `.AxesImage`
+            etc.) described by this multivariate colorbar.
+            This argument is mandatory for the `.Figure.colorbar_multivar` method
+            but optional for the`.pyplot.colorbar_multivar` function, which sets
+            the default to the current image.
+
+        caxes : `~matplotlib.axes.Axes`, optional
+            Axes into which the colorbar will be drawn.  If `None`, then new
+            Axes are created and the space for it will be stolen from the Axes(s)
+            specified in *ax*.
+
+        ax : `~matplotlib.axes.Axes` or iterable or `numpy.ndarray` of Axes, optional
+            The one or more parent Axes from which space for a new colorbar Axes
+            will be stolen. This parameter is only used if *cax* is not set.
+
+            Defaults to the Axes that contains the mappable used to create the
+            colorbar.
+
+        Returns
+        -------
+        multivariate_colorbar : `~matplotlib.colorbar.MultivarColorbar`
+
+        Other Parameters
+        ----------------
+        %(_make_multivar_axes_kw_doc)s
+
+        """
+
+        if isinstance(mappable, mpl.colorizer.Colorizer):
+            mappable = mcolorizer.ColorizingArtist(mappable)
+        if not isinstance(mappable.colorizer.cmap, mcolors.MultivarColormap):
+            raise ValueError("A multivariate colorbar can only be used together "
+                             "with a multivariate colormap, not "
+                             f"{type(mappable.colorizer.cmap)}")
+
+        n_variates = mappable.colorizer.cmap.n_variates
+
+        if ax is None:
+            ax = getattr(mappable, "axes", None)
+
+        cbar_info = None
+        if caxes is None:
+            if ax is None:
+                raise ValueError(
+                    'Unable to determine Axes to steal space for Colorbar. '
+                    'Either provide the *cax* argument to use as the Axes for '
+                    'the Colorbar, provide the *ax* argument to steal space '
+                    'from it, or add *mappable* to an Axes.')
+            fig = (  # Figure of first Axes; logic copied from make_axes.
+                [*ax.flat] if isinstance(ax, np.ndarray)
+                else [*ax] if np.iterable(ax)
+                else [ax])[0].get_figure(root=False)
+            current_ax = fig.gca()
+            caxes, kwargs, cbar_info = cbar.make_multivar_axes(ax, n_variates,
+                                                               n_major, **kwargs)
+            # make_axes calls add_{axes,subplot} which changes gca; undo that.
+            fig.sca(current_ax)
+            for cax in caxes:
+                cax.grid(visible=False, which='both', axis='both')
+
+        if (hasattr(mappable, "get_figure") and
+                (mappable_host_fig := mappable.get_figure(root=True)) is not None):
+            # Warn in case of mismatch
+            if mappable_host_fig is not self._root_figure:
+                _api.warn_external(
+                        f'Adding colorbar to a different Figure '
+                        f'{repr(mappable_host_fig)} than '
+                        f'{repr(self._root_figure)} which '
+                        f'fig.colorbar is called on.')
+        NON_COLORBAR_KEYS = [  # remove kws that cannot be passed to Colorbar
+            'fraction', 'pad', 'shrink', 'anchor', 'panchor']
+
+        cb = cbar.MultivarColorbar(caxes, mappable, **{
+              k: v for k, v in kwargs.items() if k not in NON_COLORBAR_KEYS})
+        cb._set_colorbar_info(cbar_info)
+
+        for cax in caxes:
+            cax.get_figure(root=False).stale = True
+        return cb
+
     def subplots_adjust(self, left=None, bottom=None, right=None, top=None,
                         wspace=None, hspace=None):
         """
diff --git a/lib/matplotlib/figure.pyi b/lib/matplotlib/figure.pyi
index cf17f4694dbd..148c2420b2b8 100644
--- a/lib/matplotlib/figure.pyi
+++ b/lib/matplotlib/figure.pyi
@@ -14,7 +14,7 @@ from matplotlib.backend_bases import (
     RendererBase,
 )
 from matplotlib.colors import Colormap, Normalize
-from matplotlib.colorbar import Colorbar
+from matplotlib.colorbar import Colorbar, BivarColorbar, MultivarColorbar
 from matplotlib.colorizer import ColorizingArtist, Colorizer
 from matplotlib.cm import ScalarMappable
 from matplotlib.gridspec import GridSpec, SubplotSpec, SubplotParams as SubplotParams
@@ -185,6 +185,24 @@ class FigureBase(Artist):
         use_gridspec: bool = ...,
         **kwargs
     ) -> Colorbar: ...
+    def colorbar_bivar(
+        self,
+        mappable: ColorizingArtist,
+        *,
+        cax: Axes | None = ...,
+        ax: Axes | Iterable[Axes] | None = ...,
+        use_gridspec: bool = ...,
+        **kwargs
+    ) -> BivarColorbar: ...
+    def colorbar_multivar(
+        self,
+        mappable: ColorizingArtist,
+        *,
+        caxes: Iterable[Axes] | None = ...,
+        ax: Axes | Iterable[Axes] | None = ...,
+        n_major: int = ...,
+        **kwargs
+    ) -> MultivarColorbar: ...
     def subplots_adjust(
         self,
         left: float | None = ...,
diff --git a/lib/matplotlib/pyplot.py b/lib/matplotlib/pyplot.py
index 8315056c81a2..2ac20f6d55e2 100644
--- a/lib/matplotlib/pyplot.py
+++ b/lib/matplotlib/pyplot.py
@@ -126,7 +126,11 @@
         EventCollection,
         QuadMesh,
     )
-    from matplotlib.colorbar import Colorbar
+    from matplotlib.colorbar import (
+        Colorbar,
+        BivarColorbar,
+        MultivarColorbar,
+    )
     from matplotlib.container import (
         BarContainer,
         ErrorbarContainer,
@@ -2684,6 +2688,42 @@ def colorbar(
     return ret
 
 
+@_copy_docstring_and_deprecators(Figure.colorbar_bivar)
+def colorbar_bivar(
+    mappable: ColorizingArtist | None = None,
+    cax: matplotlib.axes.Axes | None = None,
+    ax: matplotlib.axes.Axes | Iterable[matplotlib.axes.Axes] | None = None,
+    **kwargs
+) -> BivarColorbar:
+    if mappable is None:
+        mappable = gci()
+        if mappable is None:
+            raise RuntimeError('No mappable was found to use for colorbar '
+                               'creation. First define a mappable such as '
+                               'an image (with imshow) or a contour set ('
+                               'with contourf).')
+    ret = gcf().colorbar_bivar(mappable, cax=cax, ax=ax, **kwargs)
+    return ret
+
+
+@_copy_docstring_and_deprecators(Figure.colorbar_multivar)
+def colorbar_multivar(
+    mappable: ColorizingArtist | None = None,
+    caxes: Iterable[matplotlib.axes.Axes] | None = None,
+    ax: matplotlib.axes.Axes | Iterable[matplotlib.axes.Axes] | None = None,
+    **kwargs
+) -> MultivarColorbar:
+    if mappable is None:
+        mappable = gci()
+        if mappable is None:
+            raise RuntimeError('No mappable was found to use for colorbar '
+                               'creation. First define a mappable such as '
+                               'an image (with imshow) or a contour set ('
+                               'with contourf).')
+    ret = gcf().colorbar_multivar(mappable, caxes=caxes, ax=ax, **kwargs)
+    return ret
+
+
 def clim(vmin: float | None = None, vmax: float | None = None) -> None:
     """
     Set the color limits of the current image.
diff --git a/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_locationing.png b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_locationing.png
new file mode 100644
index 000000000000..a81d0fb1d4ee
Binary files /dev/null and b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_locationing.png differ
diff --git a/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_locationing_constrained.png b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_locationing_constrained.png
new file mode 100644
index 000000000000..d412a50f4fcc
Binary files /dev/null and b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_locationing_constrained.png differ
diff --git a/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_locationing_gridspec.png b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_locationing_gridspec.png
new file mode 100644
index 000000000000..5abdeceae8a9
Binary files /dev/null and b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_locationing_gridspec.png differ
diff --git a/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_not_rasterized.png b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_not_rasterized.png
new file mode 100644
index 000000000000..55436d3f07ce
Binary files /dev/null and b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_not_rasterized.png differ
diff --git a/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_sharing.png b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_sharing.png
new file mode 100644
index 000000000000..243a5a944351
Binary files /dev/null and b/lib/matplotlib/tests/baseline_images/test_colorbar/bivar_cbar_sharing.png differ
diff --git a/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_locationing.png b/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_locationing.png
new file mode 100644
index 000000000000..7600a8ab99d0
Binary files /dev/null and b/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_locationing.png differ
diff --git a/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_locationing_constrained.png b/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_locationing_constrained.png
new file mode 100644
index 000000000000..57b300920c0f
Binary files /dev/null and b/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_locationing_constrained.png differ
diff --git a/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_n_major.png b/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_n_major.png
new file mode 100644
index 000000000000..8c9569bd79b8
Binary files /dev/null and b/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_n_major.png differ
diff --git a/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_sharing.png b/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_sharing.png
new file mode 100644
index 000000000000..6423085fb9e0
Binary files /dev/null and b/lib/matplotlib/tests/baseline_images/test_colorbar/multivar_cbar_sharing.png differ
diff --git a/lib/matplotlib/tests/test_colorbar.py b/lib/matplotlib/tests/test_colorbar.py
index 33174a6df45f..74d897c4e124 100644
--- a/lib/matplotlib/tests/test_colorbar.py
+++ b/lib/matplotlib/tests/test_colorbar.py
@@ -1254,3 +1254,510 @@ def test_colorbar_format_string_and_old():
     plt.imshow([[0, 1]])
     cb = plt.colorbar(format="{x}%")
     assert isinstance(cb._formatter, StrMethodFormatter)
+
+
+@pytest.mark.parametrize('use_gridspec', [True, False])
+@image_comparison(['bivar_cbar_locationing.png',
+                   ], style='mpl20',
+                  remove_text=True, savefig_kwarg={'dpi': 40}, tol=0.05)
+def test_colorbar_bivar_location(use_gridspec):
+    data = (np.arange(12).reshape((3, 4)) % 4,
+            np.arange(12).reshape((3, 4)) % 3)
+    # -------------------
+    locations = ['left', 'right', 'top', 'bottom']
+    fig, axes = plt.subplots(2, 2)
+    for i, ax in enumerate(axes.ravel()):
+        mim = ax.imshow(data, cmap='BiOrangeBlue')
+        fig.colorbar_bivar(mim, location=locations[i], use_gridspec=use_gridspec)
+
+
+@image_comparison(['bivar_cbar_locationing_constrained.png',
+                   ], style='mpl20',
+                  remove_text=True, savefig_kwarg={'dpi': 40}, tol=0.05)
+def test_colorbar_bivar_location_constrained():
+    data = (np.arange(12).reshape((3, 4)) % 4,
+            np.arange(12).reshape((3, 4)) % 3)
+    # -------------------
+    locations = ['left', 'right', 'top', 'bottom']
+    fig, axes = plt.subplots(2, 2, constrained_layout='constrained')
+    for i, ax in enumerate(axes.ravel()):
+        mim = ax.imshow(data, cmap='BiOrangeBlue')
+        fig.colorbar_bivar(mim, location=locations[i])
+
+
+@image_comparison(['multivar_cbar_locationing.png',
+                   ], style='mpl20',
+                  remove_text=True, savefig_kwarg={'dpi': 40}, tol=0.05)
+def test_colorbar_multivar_location():
+    data = (np.arange(12).reshape((3, 4)) % 4,
+            np.arange(12).reshape((3, 4)) % 3)
+    # -------------------
+    locations = ['left', 'right', 'top', 'bottom']
+    fig, axes = plt.subplots(2, 2)
+    for i, ax in enumerate(axes.ravel()):
+        mim = ax.imshow(data, cmap='2VarAddA')
+        fig.colorbar_multivar(mim, location=locations[i])
+
+
+@image_comparison(['multivar_cbar_locationing_constrained.png',
+                   ], style='mpl20',
+                  remove_text=True, savefig_kwarg={'dpi': 40}, tol=0.05)
+def test_colorbar_multivar_location_constrained():
+    data = (np.arange(12).reshape((3, 4)) % 4,
+            np.arange(12).reshape((3, 4)) % 3)
+    # -------------------
+    locations = ['left', 'right', 'top', 'bottom']
+    fig, axes = plt.subplots(2, 2, constrained_layout='constrained')
+    for i, ax in enumerate(axes.ravel()):
+        mim = ax.imshow(data, cmap='2VarAddA')
+        fig.colorbar_multivar(mim, location=locations[i])
+
+
+@pytest.mark.parametrize('use_gridspec', [True, False])
+@image_comparison(['bivar_cbar_sharing.png',
+                   ], style='mpl20',
+                  remove_text=True, savefig_kwarg={'dpi': 40}, tol=0.05)
+def test_colorbar_bivar_sharing(use_gridspec):
+    data = (np.arange(12).reshape((3, 4)) % 4,
+            np.arange(12).reshape((3, 4)) % 3)
+    # -------------------
+    plt.figure()
+    ax1 = plt.subplot(211, anchor='NE', aspect='equal')
+    plt.imshow(data, cmap='BiOrangeBlue')
+    ax2 = plt.subplot(223)
+    plt.imshow(data, cmap='BiOrangeBlue')
+    ax3 = plt.subplot(224)
+    plt.imshow(data, cmap='BiOrangeBlue')
+
+    plt.colorbar_bivar(ax=[ax2, ax3, ax1], location='right', pad=0.0, shrink=0.5,
+                       panchor=False, use_gridspec=use_gridspec)
+    plt.colorbar_bivar(ax=[ax2, ax3, ax1], location='left', shrink=0.5,
+                       panchor=False, use_gridspec=use_gridspec)
+    plt.colorbar_bivar(ax=[ax1], location='bottom', panchor=False,
+                       anchor=(0.8, 0.5), shrink=0.6, use_gridspec=use_gridspec)
+
+
+@image_comparison(['multivar_cbar_sharing.png',
+                   ], style='mpl20',
+                  remove_text=True, savefig_kwarg={'dpi': 40}, tol=0.05)
+def test_colorbar_multivar_sharing():
+    data = (np.arange(12).reshape((3, 4)) % 4,
+            np.arange(12).reshape((3, 4)) % 3,
+            np.arange(12).reshape((3, 4)) % 5)
+    # -------------------
+    plt.figure()
+    ax1 = plt.subplot(211, anchor='NE', aspect='equal')
+    plt.imshow(data, cmap='3VarAddA')
+    ax2 = plt.subplot(223)
+    plt.imshow(data, cmap='3VarAddA')
+    ax3 = plt.subplot(224)
+    plt.imshow(data, cmap='3VarAddA')
+
+    plt.colorbar_multivar(ax=[ax2, ax3, ax1], location='right')
+    plt.colorbar_multivar(ax=[ax2, ax3, ax1], location='left')
+    plt.colorbar_multivar(ax=[ax1], location='bottom')
+
+
+@pytest.mark.parametrize('constrained', [False, True],
+                         ids=['standard', 'constrained'])
+def test_bivar_cbar_single_ax_panchor_east(constrained):
+    fig, ax = plt.subplots(constrained_layout=constrained)
+    ax.set_anchor('N')
+    assert ax.get_anchor() == 'N'
+    mp = ax.imshow([[[0, 1], [1, 2]], [[2, 1], [0, 1]]], cmap='BiOrangeBlue')
+    fig.colorbar_bivar(mp, panchor='E')
+    assert ax.get_anchor() == 'E'
+
+
+def test_bivar_cbar_set_xylabel():
+    fig, axes = plt.subplots(1, 2)
+    mp = axes[0].imshow([[[0, 1], [1, 2]], [[2, 1], [0, 1]]], cmap='BiOrangeBlue')
+    cb = fig.colorbar_bivar(mp, cax=axes[1])
+    cb.set_ylabel('y')
+    cb.set_xlabel('x')
+    assert axes[1].get_ylabel() == 'y'
+    assert axes[1].get_xlabel() == 'x'
+    assert cb.xaxis is axes[1].xaxis
+    assert cb.yaxis is axes[1].yaxis
+
+
+def test_bivar_cbar_log_no_scale():
+    fig, axes = plt.subplots(1, 2)
+    mp = axes[0].imshow([[[100, 1], [10, 1]], [[0.5, 0], [0.3, 1]]],
+                        cmap='BiOrangeBlue',
+                        norm=['log', mcolors.NoNorm()],
+                        )
+    assert axes[1].get_yscale() == 'linear'
+    assert axes[1].get_xscale() == 'linear'
+    fig.colorbar_bivar(mp, cax=axes[1])
+    assert axes[1].get_yscale() == 'log'
+    assert axes[1].get_xscale() == 'function'
+
+
+def test_bivar_cbar_change_vmin_vmax():
+    fig, axes = plt.subplots(1, 2)
+    mp = axes[0].imshow([[[100, 1], [10, 1]], [[0.5, 0], [0.3, 1]]],
+                        cmap='BiOrangeBlue',
+                        )
+    fig.colorbar_bivar(mp, cax=axes[1])
+    assert np.all(axes[1].get_ylim() == np.array([1, 100]))
+    assert np.all(axes[1].get_xlim() == np.array([0, 1]))
+    mp.colorizer.norm.vmin = [-1, -2]
+    mp.colorizer.norm.vmax = [3, 5]
+    assert np.all(axes[1].get_ylim() == np.array([-1, 3]))
+    assert np.all(axes[1].get_xlim() == np.array([-2, 5]))
+
+
+def test_bivar_cbar_change_norms():
+    fig, axes = plt.subplots(1, 2)
+    mp = axes[0].imshow([[[100, 1], [10, 1]], [[0.5, 0.2], [0.3, 1]]],
+                        cmap='BiOrangeBlue',
+                        )
+    fig.colorbar_bivar(mp, cax=axes[1])
+    assert axes[1].get_yscale() == 'linear'
+    assert axes[1].get_xscale() == 'linear'
+    mp.colorizer.norm = ['log', mcolors.NoNorm()]
+    assert axes[1].get_yscale() == 'log'
+    assert axes[1].get_xscale() == 'function'
+
+
+def test_bivar_cbar_anchor():
+    # right
+    fig, ax = plt.subplots(figsize=(6, 2))
+    mp = ax.imshow([[[0, 1], [1, 2]], [[2, 1], [0, 1]]], cmap='BiOrangeBlue')
+    anchor = (1, 0.3)
+    shrink = 0.3
+    cbar = fig.colorbar_bivar(mp, anchor=anchor, shrink=shrink)
+
+    x0, y0, x1, y1 = ax.get_position().extents
+    cx0, cy0, cx1, cy1 = cbar.ax.get_position().extents
+    p0 = (y1 - y0) * anchor[1] + y0
+    np.testing.assert_allclose(
+            [cy1, cy0],
+            [y1 * shrink + (1 - shrink) * p0, p0 * (1 - shrink) + y0 * shrink])
+
+    # left
+    fig, ax = plt.subplots(figsize=(6, 2))
+    mp = ax.imshow([[[0, 1], [1, 2]], [[2, 1], [0, 1]]], cmap='BiOrangeBlue')
+    anchor = (1, 0.7)
+    shrink = 0.3
+    cbar = fig.colorbar_bivar(mp, anchor=anchor, shrink=shrink, location='left')
+
+    x0, y0, x1, y1 = ax.get_position().extents
+    cx0, cy0, cx1, cy1 = cbar.ax.get_position().extents
+    p0 = (y1 - y0) * anchor[1] + y0
+    np.testing.assert_allclose(
+            [cy1, cy0],
+            [y1 * shrink + (1 - shrink) * p0, p0 * (1 - shrink) + y0 * shrink])
+
+    # top
+    fig, ax = plt.subplots(figsize=(2, 6))
+    mp = ax.imshow([[[0, 1], [1, 2]], [[2, 1], [0, 1]]], cmap='BiOrangeBlue')
+    anchor = (0.3, 1)
+    shrink = 0.3
+    cbar = fig.colorbar_bivar(mp, anchor=anchor, shrink=shrink, location='top')
+
+    x0, y0, x1, y1 = ax.get_position().extents
+    cx0, cy0, cx1, cy1 = cbar.ax.get_position().extents
+    p0 = (x1 - x0) * anchor[0] + x0
+    np.testing.assert_allclose(
+            [cx1, cx0],
+            [x1 * shrink + (1 - shrink) * p0, p0 * (1 - shrink) + x0 * shrink])
+
+    # bottom
+    fig, ax = plt.subplots(figsize=(2, 6))
+    mp = ax.imshow([[[0, 1], [1, 2]], [[2, 1], [0, 1]]], cmap='BiOrangeBlue')
+    anchor = (0.3, 1)
+    shrink = 0.3
+    cbar = fig.colorbar_bivar(mp, anchor=anchor, shrink=shrink, location='bottom')
+
+    x0, y0, x1, y1 = ax.get_position().extents
+    cx0, cy0, cx1, cy1 = cbar.ax.get_position().extents
+    p0 = (x1 - x0) * anchor[0] + x0
+    np.testing.assert_allclose(
+            [cx1, cx0],
+            [x1 * shrink + (1 - shrink) * p0, p0 * (1 - shrink) + x0 * shrink])
+
+
+def test_bivar_cbar_aspect():
+    fig, ax = plt.subplots(1, 1)
+    mp = ax.imshow([[[100, 1], [10, 1]], [[0.5, 0.2], [0.3, 1]]],
+                   cmap='BiOrangeBlue',
+                   )
+    cbar = fig.colorbar_bivar(mp, aspect=0.3)
+    assert cbar.ax.get_box_aspect() == 0.3
+    assert cbar.ax._colorbar_info["aspect"] == 0.3
+    assert cbar.aspect == 0.3
+    cbar.aspect = 4
+    assert cbar.ax.get_box_aspect() == 4
+    assert cbar.ax._colorbar_info["aspect"] == 4
+    assert cbar.aspect == 4
+
+
+def test_bivar_cbar_alpha():
+    fig, ax = plt.subplots(1, 1)
+    mp = ax.imshow([[[100, 1], [10, 1]], [[0.5, 0.2], [0.3, 1]]],
+                   cmap='BiOrangeBlue',
+                   alpha=0.5,
+                   )
+    cbar = fig.colorbar_bivar(mp)
+    assert cbar._image._alpha == 0.5
+
+
+# If we decide in the future to disallow calling colorbar() on the "wrong" figure,
+# just delete this test.
+def test_bivar_cbar_wrong_figure():
+    fig0, ax0 = plt.subplots()
+    fig1, ax1 = plt.subplots()
+    mp = ax0.imshow([[[100, 1], [10, 1]], [[0.5, 0.2], [0.3, 1]]],
+                    cmap='BiOrangeBlue',
+                    )
+    with pytest.warns(UserWarning, match="different Figure"):
+        fig1.colorbar_bivar(mp)
+
+
+@pytest.mark.parametrize('use_gridspec', [True, False])
+@pytest.mark.parametrize('nested_gridspecs', [True, False])
+def test_bivar_cbar_remove_from_figure(nested_gridspecs, use_gridspec):
+    """Test `remove` with the specified ``use_gridspec`` setting."""
+    fig = plt.figure()
+    if nested_gridspecs:
+        gs = fig.add_gridspec(2, 2)[1, 1].subgridspec(2, 2)
+        ax = fig.add_subplot(gs[1, 1])
+    else:
+        ax = fig.add_subplot()
+    mp = ax.pcolormesh([[[100, 1], [10, 1]], [[0.5, 0.2], [0.3, 1]]],
+                       cmap='BiOrangeBlue',
+                       )
+    pre_position = ax.get_position()
+    cb = fig.colorbar_bivar(mp, use_gridspec=use_gridspec)
+    fig.subplots_adjust()
+    cb.remove()
+    fig.subplots_adjust()
+    post_position = ax.get_position()
+    assert (pre_position.get_points() == post_position.get_points()).all()
+
+
+def test_multivar_cbar_remove_from_figure():
+    """Test `remove` with the specified ``use_gridspec`` setting."""
+    fig = plt.figure()
+    ax = fig.add_subplot()
+    mp = ax.pcolormesh([[[100, 1], [10, 1]], [[0.5, 0.2], [0.3, 1]]],
+                       cmap='2VarAddA',
+                       )
+    pre_position = ax.get_position()
+    cb = fig.colorbar_multivar(mp)
+    fig.subplots_adjust()
+    cb.remove()
+    fig.subplots_adjust()
+    post_position = ax.get_position()
+    assert (pre_position.get_points() == post_position.get_points()).all()
+
+
+def test_bivar_cbar_remove_with_no_mappable():
+    fig, ax = plt.subplots()
+    norm = mpl.colors.MultiNorm(['linear', 'linear'])
+    ca = mpl.colorizer.Colorizer('BiOrangeBlue', norm)
+    cb = mpl.colorbar.BivarColorbar(ax, ca)
+    cb.remove()
+
+
+def test_multivar_cbar_from_colorizer():
+    fig, ax = plt.subplots()
+    norm = mpl.colors.MultiNorm(['linear', 'linear'])
+    ca = mpl.colorizer.Colorizer('2VarAddA', norm)
+    cb = fig.colorbar_multivar(ca, ax=ax)
+
+
+def test_multivar_cbar_from_colorizer_cax():
+    fig, axes = plt.subplots(1, 3)
+    norm = mpl.colors.MultiNorm(['linear', 'linear', 'linear'])
+    ca = mpl.colorizer.Colorizer('3VarAddA', norm)
+    cb = fig.colorbar_multivar(ca, caxes=axes)
+
+
+def test_bivar_cbar_ticklocations():
+    norm = mpl.colors.MultiNorm(['linear', 'linear'])
+    ca = mpl.colorizer.Colorizer('BiOrangeBlue', norm)
+    fig, ax = plt.subplots()
+    with pytest.raises(ValueError, match='ticklocations must be a tuple of'):
+        cbar = fig.colorbar_bivar(ca, cax=ax, ticklocations=['left'])
+    cbar = fig.colorbar_bivar(ca, cax=ax, ticklocations=['left', 'top'])
+    cbar.ax.yaxis.get_label_position() == 'left'
+    cbar.ax.xaxis.get_label_position() == 'top'
+
+
+def test_wrong_kind_colorbar():
+    fig, ax = plt.subplots(1, 1)
+    mp = ax.imshow([[[100, 1], [10, 1]], [[0.5, 0.2], [0.3, 1]]],
+                   cmap='BiOrangeBlue',
+                   )
+    with pytest.raises(ValueError, match='can only be together with a scalar colormap'):
+        fig.colorbar(mp, ax=ax)
+
+    fig, ax = plt.subplots(1, 1)
+    mp = ax.imshow([[100, 1], [10, 1]])
+    with pytest.raises(ValueError, match='bivariate colorbar can only be used '):
+        fig.colorbar_bivar(mp, ax=ax)
+
+
+def test_colorbar_bivar_set_get_view():
+    # maybe not the best way to test this, this is normally used
+    # interacitively ._get_view() and ._set_view() are normally
+    # used in interactive mode
+    fig, axes = plt.subplots(1, 1)
+    mp = axes.imshow([[[100, 1], [10, 1]], [[0.5, 0.2], [0.3, 1]]],
+                     cmap='BiOrangeBlue',
+                     )
+    cb = fig.colorbar_bivar(mp, ax=axes)
+    # limits are x: (0.2, 1), y: (1, 100)
+    view = cb._get_view()
+    view = np.array(view) + 1
+    cb._set_view(view)
+    # limits are x: (1.2, 2), y: (2, 101)
+    assert mp.colorizer.norm.norms[1].vmin == 1.2
+    assert mp.colorizer.norm.norms[1].vmax == 2
+    assert mp.colorizer.norm.norms[0].vmin == 2
+    assert mp.colorizer.norm.norms[0].vmax == 101
+
+
+def test_colorbar_bivar_no_ax():
+    norm = mpl.colors.MultiNorm(['linear', 'linear'])
+    ca = mpl.colorizer.Colorizer('BiOrangeBlue', norm)
+    fig, ax = plt.subplots()
+    with pytest.raises(ValueError, match='Unable to determine Axes'):
+        fig.colorbar_bivar(ca)
+
+
+def test_colorbar_bivar_custom_norm():
+    """
+    This tests setting the correct scale for norms that do not have a
+    ._scale property.
+
+    In the normal case, colorbar_bivar uses the ._scale property
+    of each norm to set the transform on each axis.
+    """
+    class CustomHalfNorm(mcolors.Normalize):
+        def __init__(self):
+            super().__init__()
+
+        @property
+        def vmin(self):
+            return 0
+
+        @vmin.setter
+        def vmin(self, val):
+            ...
+
+        @property
+        def vmax(self):
+            return 1
+
+        @vmax.setter
+        def vmax(self, val):
+            ...
+
+        @property
+        def clip(self):
+            return False
+
+        @clip.setter
+        def clip(self, val):
+            ...
+
+        def __call__(self, value, clip=None):
+            return value / 2
+
+        def inverse(self, value):
+            return 2 * value
+
+        def autoscale(self, A):
+            pass
+
+        def autoscale_None(self, A):
+            pass
+
+        def scaled(self):
+            return True
+
+        @property
+        def n_components(self):
+            return 1
+
+    norm = mpl.colors.MultiNorm([CustomHalfNorm(), CustomHalfNorm()])
+    ca = mpl.colorizer.Colorizer('BiOrangeBlue', norm)
+    fig, ax = plt.subplots()
+    fig.colorbar_bivar(ca, cax=ax)
+    assert ax.get_yscale() == 'function'
+    assert ax.get_xscale() == 'function'
+
+
+def test_colorbar_bivar_not_via_fig():
+    # This test makes a colorbar without calling
+    # fig.colorbar_bivar() and without a ColorizingArtist
+    norm = mpl.colors.MultiNorm(['linear', 'linear'])
+    ca = mpl.colorizer.Colorizer('BiOrangeBlue', norm)
+    fig, ax = plt.subplots()
+    cb = mpl.colorbar.BivarColorbar(ax, ca)
+    assert cb.colorizer is ca
+
+
+def test_remove_colorbar_with_no_mappable():
+    fig, ax = plt.subplots()
+    cb = mpl.colorbar.Colorbar(ax)
+    cb.remove()
+
+
+@image_comparison(['bivar_cbar_not_rasterized.png',
+                   ], style='mpl20',
+                  remove_text=True, savefig_kwarg={'dpi': 40}, tol=0.05)
+def test_bivar_cbar_not_rasterized():
+    cm = mpl.bivar_colormaps['BiOrangeBlue'].resampled((5, 3))
+    fig, axes = plt.subplots(1, 1)
+    mp = axes.imshow([[[100, 1], [10, 1]], [[0.5, 0.2], [0.3, 1]]],
+                     cmap=cm,
+                     interpolation='nearest'
+                     )
+    fig.colorbar_bivar(mp)
+
+
+@image_comparison(['multivar_cbar_n_major.png',
+                   ], style='mpl20',
+                  remove_text=True, savefig_kwarg={'dpi': 40}, tol=0.05)
+def test_colorbar_multivar_n_major():
+    data = (np.arange(12).reshape((3, 4)) % 4,
+            np.arange(12).reshape((3, 4)) % 3,
+            np.arange(12).reshape((3, 4)) % 5)
+    # -------------------
+    fig, axes = plt.subplots(4, 3, figsize=(8, 8), constrained_layout='constrained')
+    locations = ['left', 'right', 'top', 'bottom']
+    for i, axs in enumerate(axes):
+        for j, ax in enumerate(axs):
+            mim = ax.imshow(data, cmap='3VarAddA')
+            fig.colorbar_multivar(mim, n_major=j + 1, location=locations[i])
+
+    with pytest.raises(ValueError, match="cannot be zero"):
+        fig.colorbar_multivar(mim, n_major=0)
+
+
+def test_cbar_wrong_figures():
+    fig0, ax0 = plt.subplots()
+    fig1, ax1 = plt.subplots()
+    im0 = ax0.imshow([[0, 1], [2, 3]])
+    im1 = ax1.imshow([[0, 1], [2, 3]])
+    with pytest.raises(ValueError, match="not all parents share"):
+        fig0.colorbar(im0, ax=[im0, im1])
+
+
+def test_multivar_cbar_set_label_limits():
+    data = (np.arange(12).reshape((3, 4)) % 4,
+            np.arange(12).reshape((3, 4)) % 3,
+            np.arange(12).reshape((3, 4)) % 5)
+    # -------------------
+    fig, ax = plt.subplots(1, 1)
+    mim = ax.imshow(data, cmap='3VarAddA')
+    cbs = fig.colorbar_multivar(mim)
+    cbs[0].set_label('A')
+    assert len(cbs) == 3
+    mim.norm.vmin = (-1, -1, -1)
+    mim.norm.vmax = (1, 2, 3)
diff --git a/lib/mpl_toolkits/mplot3d/axes3d.py b/lib/mpl_toolkits/mplot3d/axes3d.py
index 9cc1e8eaefbe..e5572361a78c 100644
--- a/lib/mpl_toolkits/mplot3d/axes3d.py
+++ b/lib/mpl_toolkits/mplot3d/axes3d.py
@@ -813,7 +813,7 @@ def set_ybound(self, lower=None, upper=None, view_margin=None):
             is not modified.
         view_margin : float or None
             The margin to apply to the bounds. If *None*, the margin is handled
-            by `.set_ylim`.
+            by `.Axes3D.set_ylim`.
 
         See Also
         --------
@@ -2066,7 +2066,7 @@ def _get_w_centers_ranges(self):
 
     def set_zlabel(self, zlabel, fontdict=None, labelpad=None, **kwargs):
         """
-        Set zlabel.  See doc for `.set_ylabel` for description.
+        Set zlabel.  See doc for `.Axes.set_ylabel` for description.
         """
         if labelpad is not None:
             self.zaxis.labelpad = labelpad