Cleanup: use doxy sections for mathutils types

Also minor improvements & corrections to comments.

10 files changed, 1782 insertions, 1171 deletions

diff --git a/source/blender/python/mathutils/mathutils_Color.c b/source/blender/python/mathutils/mathutils_Color.c
index 0fcde229907..24744d0cb6e 100644
--- a/source/blender/python/mathutils/mathutils_Color.c
+++ b/source/blender/python/mathutils/mathutils_Color.c
@@ -22,8 +22,40 @@
 
 #define COLOR_SIZE 3
 
-/* ----------------------------------mathutils.Color() ------------------- */
-/* makes a new color for you to play with */
+/* -------------------------------------------------------------------- */
+/** \name Utilities
+ * \{ */
+
+/**
+ * \note #BaseMath_ReadCallback must be called beforehand.
+ */
+static PyObject *Color_to_tuple_ex(ColorObject *self, int ndigits)
+{
+  PyObject *ret;
+  int i;
+
+  ret = PyTuple_New(COLOR_SIZE);
+
+  if (ndigits >= 0) {
+    for (i = 0; i < COLOR_SIZE; i++) {
+      PyTuple_SET_ITEM(ret, i, PyFloat_FromDouble(double_round((double)self->col[i], ndigits)));
+    }
+  }
+  else {
+    for (i = 0; i < COLOR_SIZE; i++) {
+      PyTuple_SET_ITEM(ret, i, PyFloat_FromDouble(self->col[i]));
+    }
+  }
+
+  return ret;
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: `__new__` / `mathutils.Color()`
+ * \{ */
+
 static PyObject *Color_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
 {
   float col[3] = {0.0f, 0.0f, 0.0f};
@@ -54,29 +86,11 @@ static PyObject *Color_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
   return Color_CreatePyObject(col, type);
 }
 
-/* -----------------------------METHODS---------------------------- */
-
-/* NOTE: BaseMath_ReadCallback must be called beforehand. */
-static PyObject *Color_ToTupleExt(ColorObject *self, int ndigits)
-{
-  PyObject *ret;
-  int i;
-
-  ret = PyTuple_New(COLOR_SIZE);
-
-  if (ndigits >= 0) {
-    for (i = 0; i < COLOR_SIZE; i++) {
-      PyTuple_SET_ITEM(ret, i, PyFloat_FromDouble(double_round((double)self->col[i], ndigits)));
-    }
-  }
-  else {
-    for (i = 0; i < COLOR_SIZE; i++) {
-      PyTuple_SET_ITEM(ret, i, PyFloat_FromDouble(self->col[i]));
-    }
-  }
+/** \} */
 
-  return ret;
-}
+/* -------------------------------------------------------------------- */
+/** \name Color Methods: Color Space Conversion
+ * \{ */
 
 PyDoc_STRVAR(Color_from_scene_linear_to_srgb_doc,
              ".. function:: from_scene_linear_to_srgb()\n"
@@ -190,7 +204,11 @@ static PyObject *Color_from_rec709_linear_to_scene_linear(ColorObject *self)
   return Color_CreatePyObject(col, Py_TYPE(self));
 }
 
-/* ---------------------------- Colorspace conversion -------------- */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Methods: Color Copy/Deep-Copy
+ * \{ */
 
 PyDoc_STRVAR(Color_copy_doc,
              ".. function:: copy()\n"
@@ -218,8 +236,11 @@ static PyObject *Color_deepcopy(ColorObject *self, PyObject *args)
   return Color_copy(self);
 }
 
-/* ----------------------------print object (internal)-------------- */
-/* print the object to screen */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: `__repr__` & `__str__`
+ * \{ */
 
 static PyObject *Color_repr(ColorObject *self)
 {
@@ -229,7 +250,7 @@ static PyObject *Color_repr(ColorObject *self)
     return NULL;
   }
 
-  tuple = Color_ToTupleExt(self, -1);
+  tuple = Color_to_tuple_ex(self, -1);
 
   ret = PyUnicode_FromFormat("Color(%R)", tuple);
 
@@ -255,8 +276,12 @@ static PyObject *Color_str(ColorObject *self)
 }
 #endif
 
-/* ------------------------tp_richcmpr */
-/* returns -1 exception, 0 false, 1 true */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: Rich Compare
+ * \{ */
+
 static PyObject *Color_richcmpr(PyObject *a, PyObject *b, int op)
 {
   PyObject *res;
@@ -295,6 +320,12 @@ static PyObject *Color_richcmpr(PyObject *a, PyObject *b, int op)
   return Py_INCREF_RET(res);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: Hash (`__hash__`)
+ * \{ */
+
 static Py_hash_t Color_hash(ColorObject *self)
 {
   if (BaseMath_ReadCallback(self) == -1) {
@@ -308,15 +339,19 @@ static Py_hash_t Color_hash(ColorObject *self)
   return mathutils_array_hash(self->col, COLOR_SIZE);
 }
 
-/* ---------------------SEQUENCE PROTOCOLS------------------------ */
-/* ----------------------------len(object)------------------------ */
-/* sequence length */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: Sequence & Mapping Protocols Implementation
+ * \{ */
+
+/** Sequence length: `len(object)`. */
 static int Color_len(ColorObject *UNUSED(self))
 {
   return COLOR_SIZE;
 }
-/* ----------------------------object[]--------------------------- */
-/* sequence accessor (get) */
+
+/** Sequence accessor (get): `x = object[i]`. */
 static PyObject *Color_item(ColorObject *self, int i)
 {
   if (i < 0) {
@@ -336,8 +371,8 @@ static PyObject *Color_item(ColorObject *self, int i)
 
   return PyFloat_FromDouble(self->col[i]);
 }
-/* ----------------------------object[]------------------------- */
-/* sequence accessor (set) */
+
+/** Sequence accessor (set): `object[i] = x`. */
 static int Color_ass_item(ColorObject *self, int i, PyObject *value)
 {
   float f;
@@ -373,8 +408,8 @@ static int Color_ass_item(ColorObject *self, int i, PyObject *value)
 
   return 0;
 }
-/* ----------------------------object[z:y]------------------------ */
-/* sequence slice (get) */
+
+/** Sequence slice accessor (get): `x = object[i:j]`. */
 static PyObject *Color_slice(ColorObject *self, int begin, int end)
 {
   PyObject *tuple;
@@ -398,8 +433,8 @@ static PyObject *Color_slice(ColorObject *self, int begin, int end)
 
   return tuple;
 }
-/* ----------------------------object[z:y]------------------------ */
-/* sequence slice (set) */
+
+/** Sequence slice accessor (set): `object[i:j] = x`. */
 static int Color_ass_slice(ColorObject *self, int begin, int end, PyObject *seq)
 {
   int i, size;
@@ -436,6 +471,7 @@ static int Color_ass_slice(ColorObject *self, int begin, int end, PyObject *seq)
   return 0;
 }
 
+/** Sequence generic subscript (get): `x = object[...]`. */
 static PyObject *Color_subscript(ColorObject *self, PyObject *item)
 {
   if (PyIndex_Check(item)) {
@@ -472,6 +508,7 @@ static PyObject *Color_subscript(ColorObject *self, PyObject *item)
   return NULL;
 }
 
+/** Sequence generic subscript (set): `object[...] = x`. */
 static int Color_ass_subscript(ColorObject *self, PyObject *item, PyObject *value)
 {
   if (PyIndex_Check(item)) {
@@ -504,29 +541,13 @@ static int Color_ass_subscript(ColorObject *self, PyObject *item, PyObject *valu
   return -1;
 }
 
-/* -----------------PROTCOL DECLARATIONS-------------------------- */
-static PySequenceMethods Color_SeqMethods = {
-    (lenfunc)Color_len,              /* sq_length */
-    (binaryfunc)NULL,                /* sq_concat */
-    (ssizeargfunc)NULL,              /* sq_repeat */
-    (ssizeargfunc)Color_item,        /* sq_item */
-    NULL,                            /* sq_slice, deprecated */
-    (ssizeobjargproc)Color_ass_item, /* sq_ass_item */
-    NULL,                            /* sq_ass_slice, deprecated */
-    (objobjproc)NULL,                /* sq_contains */
-    (binaryfunc)NULL,                /* sq_inplace_concat */
-    (ssizeargfunc)NULL,              /* sq_inplace_repeat */
-};
+/** \} */
 
-static PyMappingMethods Color_AsMapping = {
-    (lenfunc)Color_len,
-    (binaryfunc)Color_subscript,
-    (objobjargproc)Color_ass_subscript,
-};
-
-/* numeric */
+/* -------------------------------------------------------------------- */
+/** \name Color Type: Numeric Protocol Implementation
+ * \{ */
 
-/* addition: obj + obj */
+/** Addition: `object + object`. */
 static PyObject *Color_add(PyObject *v1, PyObject *v2)
 {
   ColorObject *color1 = NULL, *color2 = NULL;
@@ -552,7 +573,7 @@ static PyObject *Color_add(PyObject *v1, PyObject *v2)
   return Color_CreatePyObject(col, Py_TYPE(v1));
 }
 
-/* addition in-place: obj += obj */
+/** Addition in-place: `object += object`. */
 static PyObject *Color_iadd(PyObject *v1, PyObject *v2)
 {
   ColorObject *color1 = NULL, *color2 = NULL;
@@ -579,7 +600,7 @@ static PyObject *Color_iadd(PyObject *v1, PyObject *v2)
   return v1;
 }
 
-/* subtraction: obj - obj */
+/** Subtraction: `object - object`. */
 static PyObject *Color_sub(PyObject *v1, PyObject *v2)
 {
   ColorObject *color1 = NULL, *color2 = NULL;
@@ -605,7 +626,7 @@ static PyObject *Color_sub(PyObject *v1, PyObject *v2)
   return Color_CreatePyObject(col, Py_TYPE(v1));
 }
 
-/* subtraction in-place: obj -= obj */
+/** Subtraction in-place: `object -= object`. */
 static PyObject *Color_isub(PyObject *v1, PyObject *v2)
 {
   ColorObject *color1 = NULL, *color2 = NULL;
@@ -639,6 +660,7 @@ static PyObject *color_mul_float(ColorObject *color, const float scalar)
   return Color_CreatePyObject(tcol, Py_TYPE(color));
 }
 
+/** Multiplication: `object * object`. */
 static PyObject *Color_mul(PyObject *v1, PyObject *v2)
 {
   ColorObject *color1 = NULL, *color2 = NULL;
@@ -683,6 +705,7 @@ static PyObject *Color_mul(PyObject *v1, PyObject *v2)
   return NULL;
 }
 
+/** Division: `object / object`. */
 static PyObject *Color_div(PyObject *v1, PyObject *v2)
 {
   ColorObject *color1 = NULL;
@@ -716,7 +739,7 @@ static PyObject *Color_div(PyObject *v1, PyObject *v2)
   return NULL;
 }
 
-/* multiplication in-place: obj *= obj */
+/** Multiplication in-place: `object *= object`. */
 static PyObject *Color_imul(PyObject *v1, PyObject *v2)
 {
   ColorObject *color = (ColorObject *)v1;
@@ -744,7 +767,7 @@ static PyObject *Color_imul(PyObject *v1, PyObject *v2)
   return v1;
 }
 
-/* multiplication in-place: obj *= obj */
+/** Division in-place: `object *= object`. */
 static PyObject *Color_idiv(PyObject *v1, PyObject *v2)
 {
   ColorObject *color = (ColorObject *)v1;
@@ -777,8 +800,7 @@ static PyObject *Color_idiv(PyObject *v1, PyObject *v2)
   return v1;
 }
 
-/* -obj
- * returns the negative of this object */
+/** Negative (returns the negative of this object): `-object`. */
 static PyObject *Color_neg(ColorObject *self)
 {
   float tcol[COLOR_SIZE];
@@ -791,6 +813,31 @@ static PyObject *Color_neg(ColorObject *self)
   return Color_CreatePyObject(tcol, Py_TYPE(self));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: Protocol Declarations
+ * \{ */
+
+static PySequenceMethods Color_SeqMethods = {
+    (lenfunc)Color_len,              /*sq_length*/
+    (binaryfunc)NULL,                /*sq_concat*/
+    (ssizeargfunc)NULL,              /*sq_repeat*/
+    (ssizeargfunc)Color_item,        /*sq_item*/
+    NULL,                            /*sq_slice(DEPRECATED)*/
+    (ssizeobjargproc)Color_ass_item, /*sq_ass_item*/
+    NULL,                            /*sq_ass_slice(DEPRECATED)*/
+    (objobjproc)NULL,                /*sq_contains*/
+    (binaryfunc)NULL,                /*sq_inplace_concat*/
+    (ssizeargfunc)NULL,              /*sq_inplace_repeat*/
+};
+
+static PyMappingMethods Color_AsMapping = {
+    (lenfunc)Color_len,
+    (binaryfunc)Color_subscript,
+    (objobjargproc)Color_ass_subscript,
+};
+
 static PyNumberMethods Color_NumMethods = {
     (binaryfunc)Color_add, /*nb_add*/
     (binaryfunc)Color_sub, /*nb_subtract*/
@@ -811,24 +858,31 @@ static PyNumberMethods Color_NumMethods = {
     NULL,                  /*nb_int*/
     NULL,                  /*nb_reserved*/
     NULL,                  /*nb_float*/
-    Color_iadd,            /* nb_inplace_add */
-    Color_isub,            /* nb_inplace_subtract */
-    Color_imul,            /* nb_inplace_multiply */
-    NULL,                  /* nb_inplace_remainder */
-    NULL,                  /* nb_inplace_power */
-    NULL,                  /* nb_inplace_lshift */
-    NULL,                  /* nb_inplace_rshift */
-    NULL,                  /* nb_inplace_and */
-    NULL,                  /* nb_inplace_xor */
-    NULL,                  /* nb_inplace_or */
-    NULL,                  /* nb_floor_divide */
-    Color_div,             /* nb_true_divide */
-    NULL,                  /* nb_inplace_floor_divide */
-    Color_idiv,            /* nb_inplace_true_divide */
-    NULL,                  /* nb_index */
+    Color_iadd,            /*nb_inplace_add*/
+    Color_isub,            /*nb_inplace_subtract*/
+    Color_imul,            /*nb_inplace_multiply*/
+    NULL,                  /*nb_inplace_remainder*/
+    NULL,                  /*nb_inplace_power*/
+    NULL,                  /*nb_inplace_lshift*/
+    NULL,                  /*nb_inplace_rshift*/
+    NULL,                  /*nb_inplace_and*/
+    NULL,                  /*nb_inplace_xor*/
+    NULL,                  /*nb_inplace_or*/
+    NULL,                  /*nb_floor_divide*/
+    Color_div,             /*nb_true_divide*/
+    NULL,                  /*nb_inplace_floor_divide*/
+    Color_idiv,            /*nb_inplace_true_divide*/
+    NULL,                  /*nb_index*/
 };
 
-/* color channel, vector.r/g/b */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: Get/Set Item Implementation
+ * \{ */
+
+/* Color channel (RGB): `color.r/g/b`. */
+
 PyDoc_STRVAR(Color_channel_r_doc, "Red color channel.\n\n:type: float");
 PyDoc_STRVAR(Color_channel_g_doc, "Green color channel.\n\n:type: float");
 PyDoc_STRVAR(Color_channel_b_doc, "Blue color channel.\n\n:type: float");
@@ -843,7 +897,8 @@ static int Color_channel_set(ColorObject *self, PyObject *value, void *type)
   return Color_ass_item(self, POINTER_AS_INT(type), value);
 }
 
-/* color channel (HSV), color.h/s/v */
+/* Color channel (HSV): `color.h/s/v`. */
+
 PyDoc_STRVAR(Color_channel_hsv_h_doc, "HSV Hue component in [0, 1].\n\n:type: float");
 PyDoc_STRVAR(Color_channel_hsv_s_doc, "HSV Saturation component in [0, 1].\n\n:type: float");
 PyDoc_STRVAR(Color_channel_hsv_v_doc, "HSV Value component in [0, 1].\n\n:type: float");
@@ -891,8 +946,8 @@ static int Color_channel_hsv_set(ColorObject *self, PyObject *value, void *type)
   return 0;
 }
 
-/* color channel (HSV), color.h/s/v */
 PyDoc_STRVAR(Color_hsv_doc, "HSV Values in [0, 1].\n\n:type: float triplet");
+/** Color channel HSV (get): `x = color.hsv`. */
 static PyObject *Color_hsv_get(ColorObject *self, void *UNUSED(closure))
 {
   float hsv[3];
@@ -910,6 +965,7 @@ static PyObject *Color_hsv_get(ColorObject *self, void *UNUSED(closure))
   return ret;
 }
 
+/** Color channel HSV (set): `color.hsv = x`. */
 static int Color_hsv_set(ColorObject *self, PyObject *value, void *UNUSED(closure))
 {
   float hsv[3];
@@ -932,9 +988,12 @@ static int Color_hsv_set(ColorObject *self, PyObject *value, void *UNUSED(closur
   return 0;
 }
 
-/*****************************************************************************/
-/* Python attributes get/set structure:                                      */
-/*****************************************************************************/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: Get/Set Item Definitions
+ * \{ */
+
 static PyGetSetDef Color_getseters[] = {
     {"r", (getter)Color_channel_get, (setter)Color_channel_set, Color_channel_r_doc, (void *)0},
     {"g", (getter)Color_channel_get, (setter)Color_channel_set, Color_channel_g_doc, (void *)1},
@@ -977,7 +1036,12 @@ static PyGetSetDef Color_getseters[] = {
     {NULL, NULL, NULL, NULL, NULL} /* Sentinel */
 };
 
-/* -----------------------METHOD DEFINITIONS ---------------------- */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: Method Definitions
+ * \{ */
+
 static struct PyMethodDef Color_methods[] = {
     {"copy", (PyCFunction)Color_copy, METH_NOARGS, Color_copy_doc},
     {"__copy__", (PyCFunction)Color_copy, METH_NOARGS, Color_copy_doc},
@@ -1022,7 +1086,12 @@ static struct PyMethodDef Color_methods[] = {
     {NULL, NULL, 0, NULL},
 };
 
-/* ------------------PY_OBECT DEFINITION-------------------------- */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: Python Object Definition
+ * \{ */
+
 PyDoc_STRVAR(
     color_doc,
     ".. class:: Color(rgb)\n"
@@ -1087,6 +1156,12 @@ PyTypeObject color_Type = {
     NULL,                                                          /* tp_del */
 };
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Color Type: C/API Constructors
+ * \{ */
+
 PyObject *Color_CreatePyObject(const float col[3], PyTypeObject *base_type)
 {
   ColorObject *self;
@@ -1156,3 +1231,5 @@ PyObject *Color_CreatePyObject_cb(PyObject *cb_user, uchar cb_type, uchar cb_sub
 
   return (PyObject *)self;
 }
+
+/** \} */
diff --git a/source/blender/python/mathutils/mathutils_Color.h b/source/blender/python/mathutils/mathutils_Color.h
index 11a936b7bc8..b4fd2eeaa81 100644
--- a/source/blender/python/mathutils/mathutils_Color.h
+++ b/source/blender/python/mathutils/mathutils_Color.h
@@ -19,7 +19,8 @@ typedef struct {
  * be stored in py_data) or be a wrapper for data allocated through
  * Blender (stored in blend_data). This is an either/or struct not both. */
 
-/* prototypes */
+/* Prototypes. */
+
 PyObject *Color_CreatePyObject(const float col[3],
                                PyTypeObject *base_type) ATTR_WARN_UNUSED_RESULT;
 PyObject *Color_CreatePyObject_wrap(float col[3], PyTypeObject *base_type) ATTR_WARN_UNUSED_RESULT
diff --git a/source/blender/python/mathutils/mathutils_Euler.c b/source/blender/python/mathutils/mathutils_Euler.c
index 909c19e8cd2..f49868dfba7 100644
--- a/source/blender/python/mathutils/mathutils_Euler.c
+++ b/source/blender/python/mathutils/mathutils_Euler.c
@@ -19,45 +19,11 @@
 
 #define EULER_SIZE 3
 
-/* ----------------------------------mathutils.Euler() ------------------- */
-/* makes a new euler for you to play with */
-static PyObject *Euler_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
-{
-  PyObject *seq = NULL;
-  const char *order_str = NULL;
-
-  float eul[EULER_SIZE] = {0.0f, 0.0f, 0.0f};
-  short order = EULER_ORDER_XYZ;
-
-  if (kwds && PyDict_Size(kwds)) {
-    PyErr_SetString(PyExc_TypeError,
-                    "mathutils.Euler(): "
-                    "takes no keyword args");
-    return NULL;
-  }
-
-  if (!PyArg_ParseTuple(args, "|Os:mathutils.Euler", &seq, &order_str)) {
-    return NULL;
-  }
-
-  switch (PyTuple_GET_SIZE(args)) {
-    case 0:
-      break;
-    case 2:
-      if ((order = euler_order_from_string(order_str, "mathutils.Euler()")) == -1) {
-        return NULL;
-      }
-      ATTR_FALLTHROUGH;
-    case 1:
-      if (mathutils_array_parse(eul, EULER_SIZE, EULER_SIZE, seq, "mathutils.Euler()") == -1) {
-        return NULL;
-      }
-      break;
-  }
-  return Euler_CreatePyObject(eul, order, type);
-}
+/* -------------------------------------------------------------------- */
+/** \name Utilities
+ * \{ */
 
-/* internal use, assume read callback is done */
+/** Internal use, assume read callback is done. */
 static const char *euler_order_str(EulerObject *self)
 {
   static const char order[][4] = {"XYZ", "XZY", "YXZ", "YZX", "ZXY", "ZYX"};
@@ -96,8 +62,10 @@ short euler_order_from_string(const char *str, const char *error_prefix)
   return -1;
 }
 
-/* NOTE: BaseMath_ReadCallback must be called beforehand. */
-static PyObject *Euler_ToTupleExt(EulerObject *self, int ndigits)
+/**
+ * \note #BaseMath_ReadCallback must be called beforehand.
+ */
+static PyObject *Euler_to_tuple_ex(EulerObject *self, int ndigits)
 {
   PyObject *ret;
   int i;
@@ -118,8 +86,53 @@ static PyObject *Euler_ToTupleExt(EulerObject *self, int ndigits)
   return ret;
 }
 
-/* -----------------------------METHODS----------------------------
- * return a quaternion representation of the euler */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: `__new__` / `mathutils.Euler()`
+ * \{ */
+
+static PyObject *Euler_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+  PyObject *seq = NULL;
+  const char *order_str = NULL;
+
+  float eul[EULER_SIZE] = {0.0f, 0.0f, 0.0f};
+  short order = EULER_ORDER_XYZ;
+
+  if (kwds && PyDict_Size(kwds)) {
+    PyErr_SetString(PyExc_TypeError,
+                    "mathutils.Euler(): "
+                    "takes no keyword args");
+    return NULL;
+  }
+
+  if (!PyArg_ParseTuple(args, "|Os:mathutils.Euler", &seq, &order_str)) {
+    return NULL;
+  }
+
+  switch (PyTuple_GET_SIZE(args)) {
+    case 0:
+      break;
+    case 2:
+      if ((order = euler_order_from_string(order_str, "mathutils.Euler()")) == -1) {
+        return NULL;
+      }
+      ATTR_FALLTHROUGH;
+    case 1:
+      if (mathutils_array_parse(eul, EULER_SIZE, EULER_SIZE, seq, "mathutils.Euler()") == -1) {
+        return NULL;
+      }
+      break;
+  }
+  return Euler_CreatePyObject(eul, order, type);
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Methods
+ * \{ */
 
 PyDoc_STRVAR(Euler_to_quaternion_doc,
              ".. method:: to_quaternion()\n"
@@ -141,7 +154,6 @@ static PyObject *Euler_to_quaternion(EulerObject *self)
   return Quaternion_CreatePyObject(quat, NULL);
 }
 
-/* return a matrix representation of the euler */
 PyDoc_STRVAR(Euler_to_matrix_doc,
              ".. method:: to_matrix()\n"
              "\n"
@@ -279,9 +291,6 @@ static PyObject *Euler_make_compatible(EulerObject *self, PyObject *value)
   Py_RETURN_NONE;
 }
 
-/* ----------------------------Euler.rotate()-----------------------
- * return a copy of the euler */
-
 PyDoc_STRVAR(Euler_copy_doc,
              ".. function:: copy()\n"
              "\n"
@@ -308,8 +317,11 @@ static PyObject *Euler_deepcopy(EulerObject *self, PyObject *args)
   return Euler_copy(self);
 }
 
-/* ----------------------------print object (internal)--------------
- * print the object to screen */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: `__repr__` & `__str__`
+ * \{ */
 
 static PyObject *Euler_repr(EulerObject *self)
 {
@@ -319,7 +331,7 @@ static PyObject *Euler_repr(EulerObject *self)
     return NULL;
   }
 
-  tuple = Euler_ToTupleExt(self, -1);
+  tuple = Euler_to_tuple_ex(self, -1);
 
   ret = PyUnicode_FromFormat("Euler(%R, '%s')", tuple, euler_order_str(self));
 
@@ -349,6 +361,12 @@ static PyObject *Euler_str(EulerObject *self)
 }
 #endif
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: Rich Compare
+ * \{ */
+
 static PyObject *Euler_richcmpr(PyObject *a, PyObject *b, int op)
 {
   PyObject *res;
@@ -390,6 +408,12 @@ static PyObject *Euler_richcmpr(PyObject *a, PyObject *b, int op)
   return Py_INCREF_RET(res);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: Hash (`__hash__`)
+ * \{ */
+
 static Py_hash_t Euler_hash(EulerObject *self)
 {
   if (BaseMath_ReadCallback(self) == -1) {
@@ -403,15 +427,19 @@ static Py_hash_t Euler_hash(EulerObject *self)
   return mathutils_array_hash(self->eul, EULER_SIZE);
 }
 
-/* ---------------------SEQUENCE PROTOCOLS------------------------ */
-/* ----------------------------len(object)------------------------ */
-/* sequence length */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: Sequence Protocol
+ * \{ */
+
+/** Sequence length: `len(object)`. */
 static int Euler_len(EulerObject *UNUSED(self))
 {
   return EULER_SIZE;
 }
-/* ----------------------------object[]--------------------------- */
-/* sequence accessor (get) */
+
+/** Sequence accessor (get): `x = object[i]`. */
 static PyObject *Euler_item(EulerObject *self, int i)
 {
   if (i < 0) {
@@ -431,8 +459,8 @@ static PyObject *Euler_item(EulerObject *self, int i)
 
   return PyFloat_FromDouble(self->eul[i]);
 }
-/* ----------------------------object[]------------------------- */
-/* sequence accessor (set) */
+
+/** Sequence accessor (set): `object[i] = x`. */
 static int Euler_ass_item(EulerObject *self, int i, PyObject *value)
 {
   float f;
@@ -468,8 +496,8 @@ static int Euler_ass_item(EulerObject *self, int i, PyObject *value)
 
   return 0;
 }
-/* ----------------------------object[z:y]------------------------ */
-/* sequence slice (get) */
+
+/** Sequence slice accessor (get): `x = object[i:j]`. */
 static PyObject *Euler_slice(EulerObject *self, int begin, int end)
 {
   PyObject *tuple;
@@ -493,8 +521,8 @@ static PyObject *Euler_slice(EulerObject *self, int begin, int end)
 
   return tuple;
 }
-/* ----------------------------object[z:y]------------------------ */
-/* sequence slice (set) */
+
+/** Sequence slice accessor (set): `object[i:j] = x`. */
 static int Euler_ass_slice(EulerObject *self, int begin, int end, PyObject *seq)
 {
   int i, size;
@@ -531,6 +559,7 @@ static int Euler_ass_slice(EulerObject *self, int begin, int end, PyObject *seq)
   return 0;
 }
 
+/** Sequence generic subscript (get): `x = object[...]`. */
 static PyObject *Euler_subscript(EulerObject *self, PyObject *item)
 {
   if (PyIndex_Check(item)) {
@@ -567,6 +596,7 @@ static PyObject *Euler_subscript(EulerObject *self, PyObject *item)
   return NULL;
 }
 
+/** Sequence generic subscript (set): `object[...] = x`. */
 static int Euler_ass_subscript(EulerObject *self, PyObject *item, PyObject *value)
 {
   if (PyIndex_Check(item)) {
@@ -599,18 +629,23 @@ static int Euler_ass_subscript(EulerObject *self, PyObject *item, PyObject *valu
   return -1;
 }
 
-/* -----------------PROTCOL DECLARATIONS-------------------------- */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: Sequence & Mapping Protocol Declarations
+ * \{ */
+
 static PySequenceMethods Euler_SeqMethods = {
-    (lenfunc)Euler_len,              /* sq_length */
-    (binaryfunc)NULL,                /* sq_concat */
-    (ssizeargfunc)NULL,              /* sq_repeat */
-    (ssizeargfunc)Euler_item,        /* sq_item */
-    (ssizessizeargfunc)NULL,         /* sq_slice (deprecated) */
-    (ssizeobjargproc)Euler_ass_item, /* sq_ass_item */
-    (ssizessizeobjargproc)NULL,      /* sq_ass_slice (deprecated) */
-    (objobjproc)NULL,                /* sq_contains */
-    (binaryfunc)NULL,                /* sq_inplace_concat */
-    (ssizeargfunc)NULL,              /* sq_inplace_repeat */
+    (lenfunc)Euler_len,              /*sq_length*/
+    (binaryfunc)NULL,                /*sq_concat*/
+    (ssizeargfunc)NULL,              /*sq_repeat*/
+    (ssizeargfunc)Euler_item,        /*sq_item*/
+    (ssizessizeargfunc)NULL,         /*sq_slice(DEPRECATED)*/
+    (ssizeobjargproc)Euler_ass_item, /*sq_ass_item*/
+    (ssizessizeobjargproc)NULL,      /*sq_ass_slice(DEPRECATED)*/
+    (objobjproc)NULL,                /*sq_contains*/
+    (binaryfunc)NULL,                /*sq_inplace_concat*/
+    (ssizeargfunc)NULL,              /*sq_inplace_repeat*/
 };
 
 static PyMappingMethods Euler_AsMapping = {
@@ -619,7 +654,13 @@ static PyMappingMethods Euler_AsMapping = {
     (objobjargproc)Euler_ass_subscript,
 };
 
-/* euler axis, euler.x/y/z */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: Get/Set Item Implementation
+ * \{ */
+
+/* Euler axis: `euler.x/y/z`. */
 
 PyDoc_STRVAR(Euler_axis_doc, "Euler axis angle in radians.\n\n:type: float");
 static PyObject *Euler_axis_get(EulerObject *self, void *type)
@@ -632,7 +673,7 @@ static int Euler_axis_set(EulerObject *self, PyObject *value, void *type)
   return Euler_ass_item(self, POINTER_AS_INT(type), value);
 }
 
-/* rotation order */
+/* Euler rotation order: `euler.order`. */
 
 PyDoc_STRVAR(
     Euler_order_doc,
@@ -666,9 +707,12 @@ static int Euler_order_set(EulerObject *self, PyObject *value, void *UNUSED(clos
   return 0;
 }
 
-/*****************************************************************************/
-/* Python attributes get/set structure:                                      */
-/*****************************************************************************/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: Get/Set Item Definitions
+ * \{ */
+
 static PyGetSetDef Euler_getseters[] = {
     {"x", (getter)Euler_axis_get, (setter)Euler_axis_set, Euler_axis_doc, (void *)0},
     {"y", (getter)Euler_axis_get, (setter)Euler_axis_set, Euler_axis_doc, (void *)1},
@@ -694,7 +738,12 @@ static PyGetSetDef Euler_getseters[] = {
     {NULL, NULL, NULL, NULL, NULL} /* Sentinel */
 };
 
-/* -----------------------METHOD DEFINITIONS ---------------------- */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: Method Definitions
+ * \{ */
+
 static struct PyMethodDef Euler_methods[] = {
     {"zero", (PyCFunction)Euler_zero, METH_NOARGS, Euler_zero_doc},
     {"to_matrix", (PyCFunction)Euler_to_matrix, METH_NOARGS, Euler_to_matrix_doc},
@@ -711,7 +760,12 @@ static struct PyMethodDef Euler_methods[] = {
     {NULL, NULL, 0, NULL},
 };
 
-/* ------------------PY_OBECT DEFINITION-------------------------- */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: Python Object Definition
+ * \{ */
+
 PyDoc_STRVAR(
     euler_doc,
     ".. class:: Euler(angles, order='XYZ')\n"
@@ -776,6 +830,12 @@ PyTypeObject euler_Type = {
     NULL,                                                          /* tp_del */
 };
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Euler Type: C/API Constructors
+ * \{ */
+
 PyObject *Euler_CreatePyObject(const float eul[3], const short order, PyTypeObject *base_type)
 {
   EulerObject *self;
@@ -849,3 +909,5 @@ PyObject *Euler_CreatePyObject_cb(PyObject *cb_user,
 
   return (PyObject *)self;
 }
+
+/** \} */
diff --git a/source/blender/python/mathutils/mathutils_Euler.h b/source/blender/python/mathutils/mathutils_Euler.h
index 4438ee380af..ca2ca26c90a 100644
--- a/source/blender/python/mathutils/mathutils_Euler.h
+++ b/source/blender/python/mathutils/mathutils_Euler.h
@@ -22,6 +22,7 @@ typedef struct {
  * blender (stored in blend_data). This is an either/or struct not both */
 
 /* prototypes */
+
 PyObject *Euler_CreatePyObject(const float eul[3],
                                short order,
                                PyTypeObject *base_type) ATTR_WARN_UNUSED_RESULT;
diff --git a/source/blender/python/mathutils/mathutils_Matrix.c b/source/blender/python/mathutils/mathutils_Matrix.c
index 76b5424711f..8cd7a5c7d87 100644
--- a/source/blender/python/mathutils/mathutils_Matrix.c
+++ b/source/blender/python/mathutils/mathutils_Matrix.c
@@ -32,6 +32,10 @@ static PyObject *matrix__apply_to_copy(PyObject *(*matrix_func)(MatrixObject *),
                                        MatrixObject *self);
 static PyObject *MatrixAccess_CreatePyObject(MatrixObject *matrix, const eMatrixAccess_t type);
 
+/* -------------------------------------------------------------------- */
+/** \name Utilities
+ * \{ */
+
 static int matrix_row_vector_check(MatrixObject *mat, VectorObject *vec, int row)
 {
   if ((vec->vec_num != mat->col_num) || (row >= mat->row_num)) {
@@ -56,9 +60,242 @@ static int matrix_col_vector_check(MatrixObject *mat, VectorObject *vec, int col
   return 1;
 }
 
-/* ----------------------------------------------------------------------------
- * matrix row callbacks
- * this is so you can do matrix[i][j] = val OR matrix.row[i][j] = val */
+/** When a matrix is 4x4 size but initialized as a 3x3, re-assign values for 4x4. */
+static void matrix_3x3_as_4x4(float mat[16])
+{
+  mat[10] = mat[8];
+  mat[9] = mat[7];
+  mat[8] = mat[6];
+  mat[7] = 0.0f;
+  mat[6] = mat[5];
+  mat[5] = mat[4];
+  mat[4] = mat[3];
+  mat[3] = 0.0f;
+}
+
+void matrix_as_3x3(float mat[3][3], MatrixObject *self)
+{
+  copy_v3_v3(mat[0], MATRIX_COL_PTR(self, 0));
+  copy_v3_v3(mat[1], MATRIX_COL_PTR(self, 1));
+  copy_v3_v3(mat[2], MATRIX_COL_PTR(self, 2));
+}
+
+static void matrix_copy(MatrixObject *mat_dst, const MatrixObject *mat_src)
+{
+  BLI_assert((mat_dst->col_num == mat_src->col_num) && (mat_dst->row_num == mat_src->row_num));
+  BLI_assert(mat_dst != mat_src);
+
+  memcpy(mat_dst->matrix, mat_src->matrix, sizeof(float) * (mat_dst->col_num * mat_dst->row_num));
+}
+
+static void matrix_unit_internal(MatrixObject *self)
+{
+  const int mat_size = sizeof(float) * (self->col_num * self->row_num);
+  memset(self->matrix, 0x0, mat_size);
+  const int col_row_max = min_ii(self->col_num, self->row_num);
+  const int row_num = self->row_num;
+  for (int col = 0; col < col_row_max; col++) {
+    self->matrix[(col * row_num) + col] = 1.0f;
+  }
+}
+
+/** Transposes memory layout, row/columns don't have to match. */
+static void matrix_transpose_internal(float mat_dst_fl[], const MatrixObject *mat_src)
+{
+  ushort col, row;
+  uint i = 0;
+
+  for (row = 0; row < mat_src->row_num; row++) {
+    for (col = 0; col < mat_src->col_num; col++) {
+      mat_dst_fl[i++] = MATRIX_ITEM(mat_src, row, col);
+    }
+  }
+}
+
+/** Assumes `rowsize == colsize` is checked and the read callback has run. */
+static float matrix_determinant_internal(const MatrixObject *self)
+{
+  if (self->col_num == 2) {
+    return determinant_m2(MATRIX_ITEM(self, 0, 0),
+                          MATRIX_ITEM(self, 0, 1),
+                          MATRIX_ITEM(self, 1, 0),
+                          MATRIX_ITEM(self, 1, 1));
+  }
+  if (self->col_num == 3) {
+    return determinant_m3(MATRIX_ITEM(self, 0, 0),
+                          MATRIX_ITEM(self, 0, 1),
+                          MATRIX_ITEM(self, 0, 2),
+                          MATRIX_ITEM(self, 1, 0),
+                          MATRIX_ITEM(self, 1, 1),
+                          MATRIX_ITEM(self, 1, 2),
+                          MATRIX_ITEM(self, 2, 0),
+                          MATRIX_ITEM(self, 2, 1),
+                          MATRIX_ITEM(self, 2, 2));
+  }
+
+  return determinant_m4((const float(*)[4])self->matrix);
+}
+
+static void adjoint_matrix_n(float *mat_dst, const float *mat_src, const ushort dim)
+{
+  /* calculate the classical adjoint */
+  switch (dim) {
+    case 2: {
+      adjoint_m2_m2((float(*)[2])mat_dst, (const float(*)[2])mat_src);
+      break;
+    }
+    case 3: {
+      adjoint_m3_m3((float(*)[3])mat_dst, (const float(*)[3])mat_src);
+      break;
+    }
+    case 4: {
+      adjoint_m4_m4((float(*)[4])mat_dst, (const float(*)[4])mat_src);
+      break;
+    }
+    default:
+      BLI_assert_unreachable();
+      break;
+  }
+}
+
+static void matrix_invert_with_det_n_internal(float *mat_dst,
+                                              const float *mat_src,
+                                              const float det,
+                                              const ushort dim)
+{
+  float mat[MATRIX_MAX_DIM * MATRIX_MAX_DIM];
+  ushort i, j, k;
+
+  BLI_assert(det != 0.0f);
+
+  adjoint_matrix_n(mat, mat_src, dim);
+
+  /* divide by determinant & set values */
+  k = 0;
+  for (i = 0; i < dim; i++) {   /* col_num */
+    for (j = 0; j < dim; j++) { /* row_num */
+      mat_dst[MATRIX_ITEM_INDEX_NUMROW(dim, j, i)] = mat[k++] / det;
+    }
+  }
+}
+
+/**
+ * \param r_mat: can be from `self->matrix` or not.
+ */
+static bool matrix_invert_internal(const MatrixObject *self, float *r_mat)
+{
+  float det;
+  BLI_assert(self->col_num == self->row_num);
+  det = matrix_determinant_internal(self);
+
+  if (det != 0.0f) {
+    matrix_invert_with_det_n_internal(r_mat, self->matrix, det, self->col_num);
+    return true;
+  }
+
+  return false;
+}
+
+/**
+ * Similar to `matrix_invert_internal` but should never error.
+ * \param r_mat: can be from `self->matrix` or not.
+ */
+static void matrix_invert_safe_internal(const MatrixObject *self, float *r_mat)
+{
+  float det;
+  float *in_mat = self->matrix;
+  BLI_assert(self->col_num == self->row_num);
+  det = matrix_determinant_internal(self);
+
+  if (det == 0.0f) {
+    const float eps = PSEUDOINVERSE_EPSILON;
+
+    /* We will copy self->matrix into r_mat (if needed),
+     * and modify it in place to add diagonal epsilon. */
+    in_mat = r_mat;
+
+    switch (self->col_num) {
+      case 2: {
+        float(*mat)[2] = (float(*)[2])in_mat;
+
+        if (in_mat != self->matrix) {
+          copy_m2_m2(mat, (const float(*)[2])self->matrix);
+        }
+        mat[0][0] += eps;
+        mat[1][1] += eps;
+
+        if (UNLIKELY((det = determinant_m2(mat[0][0], mat[0][1], mat[1][0], mat[1][1])) == 0.0f)) {
+          unit_m2(mat);
+          det = 1.0f;
+        }
+        break;
+      }
+      case 3: {
+        float(*mat)[3] = (float(*)[3])in_mat;
+
+        if (in_mat != self->matrix) {
+          copy_m3_m3(mat, (const float(*)[3])self->matrix);
+        }
+        mat[0][0] += eps;
+        mat[1][1] += eps;
+        mat[2][2] += eps;
+
+        if (UNLIKELY((det = determinant_m3_array(mat)) == 0.0f)) {
+          unit_m3(mat);
+          det = 1.0f;
+        }
+        break;
+      }
+      case 4: {
+        float(*mat)[4] = (float(*)[4])in_mat;
+
+        if (in_mat != self->matrix) {
+          copy_m4_m4(mat, (const float(*)[4])self->matrix);
+        }
+        mat[0][0] += eps;
+        mat[1][1] += eps;
+        mat[2][2] += eps;
+        mat[3][3] += eps;
+
+        if (UNLIKELY(det = determinant_m4(mat)) == 0.0f) {
+          unit_m4(mat);
+          det = 1.0f;
+        }
+        break;
+      }
+      default:
+        BLI_assert_unreachable();
+    }
+  }
+
+  matrix_invert_with_det_n_internal(r_mat, in_mat, det, self->col_num);
+}
+
+static PyObject *matrix__apply_to_copy(PyObject *(*matrix_func)(MatrixObject *),
+                                       MatrixObject *self)
+{
+  PyObject *ret = Matrix_copy(self);
+  if (ret) {
+    PyObject *ret_dummy = matrix_func((MatrixObject *)ret);
+    if (ret_dummy) {
+      Py_DECREF(ret_dummy);
+      return ret;
+    }
+    /* error */
+    Py_DECREF(ret);
+    return NULL;
+  }
+
+  /* copy may fail if the read callback errors out */
+  return NULL;
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Row Callbacks
+ * This is so you can do `matrix[i][j] = val` or `matrix.row[i][j] = val`.
+ * \{ */
 
 uchar mathutils_matrix_row_cb_index = -1;
 
@@ -147,9 +384,12 @@ Mathutils_Callback mathutils_matrix_row_cb = {
     mathutils_matrix_row_set_index,
 };
 
-/* ----------------------------------------------------------------------------
- * matrix row callbacks
- * this is so you can do matrix.col[i][j] = val */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Column Callbacks
+ * This is so you can do `matrix.col[i][j] = val`.
+ * \{ */
 
 uchar mathutils_matrix_col_cb_index = -1;
 
@@ -246,10 +486,14 @@ Mathutils_Callback mathutils_matrix_col_cb = {
     mathutils_matrix_col_set_index,
 };
 
-/* ----------------------------------------------------------------------------
- * matrix row callbacks
- * this is so you can do matrix.translation = val
- * NOTE: this is _exactly like matrix.col except the 4th component is always omitted. */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Translation Callbacks
+ * This is so you can do `matrix.translation = val`.
+ *
+ * \note this is _exactly like matrix.col except the 4th component is always omitted.
+ * \{ */
 
 uchar mathutils_matrix_translation_cb_index = -1;
 
@@ -326,11 +570,12 @@ Mathutils_Callback mathutils_matrix_translation_cb = {
     mathutils_matrix_translation_set_index,
 };
 
-/* matrix column callbacks, this is so you can do `matrix.translation = Vector()`. */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: `__new__` / `mathutils.Matrix()`
+ * \{ */
 
-/* ----------------------------------mathutils.Matrix() ----------------- */
-/* mat is a 1D array of floats - row[0][0], row[0][1], row[1][0], etc. */
-/* create a new matrix type */
 static PyObject *Matrix_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
 {
   if (kwds && PyDict_Size(kwds)) {
@@ -379,41 +624,13 @@ static PyObject *Matrix_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
   return NULL;
 }
 
-static PyObject *matrix__apply_to_copy(PyObject *(*matrix_func)(MatrixObject *),
-                                       MatrixObject *self)
-{
-  PyObject *ret = Matrix_copy(self);
-  if (ret) {
-    PyObject *ret_dummy = matrix_func((MatrixObject *)ret);
-    if (ret_dummy) {
-      Py_DECREF(ret_dummy);
-      return ret;
-    }
-    /* error */
-    Py_DECREF(ret);
-    return NULL;
-  }
-
-  /* copy may fail if the read callback errors out */
-  return NULL;
-}
-
-/* when a matrix is 4x4 size but initialized as a 3x3, re-assign values for 4x4 */
-static void matrix_3x3_as_4x4(float mat[16])
-{
-  mat[10] = mat[8];
-  mat[9] = mat[7];
-  mat[8] = mat[6];
-  mat[7] = 0.0f;
-  mat[6] = mat[5];
-  mat[5] = mat[4];
-  mat[4] = mat[3];
-  mat[3] = 0.0f;
-}
+/** \} */
 
-/*-----------------------CLASS-METHODS----------------------------*/
+/* -------------------------------------------------------------------- */
+/** \name Matrix Class Methods
+ * \{ */
 
-/* mat is a 1D array of floats - row[0][0], row[0][1], row[1][0], etc. */
+/** Identity constructor: `mathutils.Matrix.Identity()`. */
 PyDoc_STRVAR(C_Matrix_Identity_doc,
              ".. classmethod:: Identity(size)\n"
              "\n"
@@ -441,6 +658,7 @@ static PyObject *C_Matrix_Identity(PyObject *cls, PyObject *args)
   return Matrix_CreatePyObject(NULL, matSize, matSize, (PyTypeObject *)cls);
 }
 
+/** Rotation constructor: `mathutils.Matrix.Rotation()`. */
 PyDoc_STRVAR(C_Matrix_Rotation_doc,
              ".. classmethod:: Rotation(angle, size, axis)\n"
              "\n"
@@ -460,25 +678,8 @@ static PyObject *C_Matrix_Rotation(PyObject *cls, PyObject *args)
   PyObject *vec = NULL;
   const char *axis = NULL;
   int matSize;
-  double angle; /* use double because of precision problems at high values */
-  float mat[16] = {
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      1.0f,
-  };
+  double angle; /* Use double because of precision problems at high values. */
+  float mat[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1};
 
   if (!PyArg_ParseTuple(args, "di|O:Matrix.Rotation", &angle, &matSize, &vec)) {
     return NULL;
@@ -545,6 +746,7 @@ static PyObject *C_Matrix_Rotation(PyObject *cls, PyObject *args)
   return Matrix_CreatePyObject(mat, matSize, matSize, (PyTypeObject *)cls);
 }
 
+/** Translation constructor: `mathutils.Matrix.Translation()`. */
 PyDoc_STRVAR(C_Matrix_Translation_doc,
              ".. classmethod:: Translation(vector)\n"
              "\n"
@@ -567,7 +769,7 @@ static PyObject *C_Matrix_Translation(PyObject *cls, PyObject *value)
 
   return Matrix_CreatePyObject(&mat[0][0], 4, 4, (PyTypeObject *)cls);
 }
-/* ----------------------------------mathutils.Matrix.Diagonal() ------------- */
+
 PyDoc_STRVAR(C_Matrix_Diagonal_doc,
              ".. classmethod:: Diagonal(vector)\n"
              "\n"
@@ -577,6 +779,7 @@ PyDoc_STRVAR(C_Matrix_Diagonal_doc,
              "   :type vector: :class:`Vector`\n"
              "   :return: A diagonal matrix.\n"
              "   :rtype: :class:`Matrix`\n");
+/** Diagonal constructor: `mathutils.Matrix.Diagonal()`. */
 static PyObject *C_Matrix_Diagonal(PyObject *cls, PyObject *value)
 {
   float mat[16] = {0.0f};
@@ -595,8 +798,8 @@ static PyObject *C_Matrix_Diagonal(PyObject *cls, PyObject *value)
 
   return Matrix_CreatePyObject(mat, size, size, (PyTypeObject *)cls);
 }
-/* ----------------------------------mathutils.Matrix.Scale() ------------- */
-/* mat is a 1D array of floats - row[0][0], row[0][1], row[1][0], etc. */
+
+/** Scale constructor: `mathutils.Matrix.Scale()`. */
 PyDoc_STRVAR(C_Matrix_Scale_doc,
              ".. classmethod:: Scale(factor, size, axis)\n"
              "\n"
@@ -617,24 +820,7 @@ static PyObject *C_Matrix_Scale(PyObject *cls, PyObject *args)
   float tvec[3];
   float factor;
   int matSize;
-  float mat[16] = {
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      1.0f,
-  };
+  float mat[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1};
 
   if (!PyArg_ParseTuple(args, "fi|O:Matrix.Scale", &factor, &matSize, &vec)) {
     return NULL;
@@ -700,8 +886,7 @@ static PyObject *C_Matrix_Scale(PyObject *cls, PyObject *args)
   /* pass to matrix creation */
   return Matrix_CreatePyObject(mat, matSize, matSize, (PyTypeObject *)cls);
 }
-/* ----------------------------------mathutils.Matrix.OrthoProjection() --- */
-/* mat is a 1D array of floats - row[0][0], row[0][1], row[1][0], etc. */
+/** Orthographic projection constructor: `mathutils.Matrix.OrthoProjection()`. */
 PyDoc_STRVAR(C_Matrix_OrthoProjection_doc,
              ".. classmethod:: OrthoProjection(axis, size)\n"
              "\n"
@@ -721,24 +906,7 @@ static PyObject *C_Matrix_OrthoProjection(PyObject *cls, PyObject *args)
 
   int matSize, x;
   float norm = 0.0f;
-  float mat[16] = {
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      1.0f,
-  };
+  float mat[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1};
 
   if (!PyArg_ParseTuple(args, "Oi:Matrix.OrthoProjection", &axis, &matSize)) {
     return NULL;
@@ -837,6 +1005,7 @@ static PyObject *C_Matrix_OrthoProjection(PyObject *cls, PyObject *args)
   return Matrix_CreatePyObject(mat, matSize, matSize, (PyTypeObject *)cls);
 }
 
+/** Shear constructor: `mathutils.Matrix.Shear()`. */
 PyDoc_STRVAR(C_Matrix_Shear_doc,
              ".. classmethod:: Shear(plane, size, factor)\n"
              "\n"
@@ -857,24 +1026,7 @@ static PyObject *C_Matrix_Shear(PyObject *cls, PyObject *args)
   int matSize;
   const char *plane;
   PyObject *fac;
-  float mat[16] = {
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      0.0f,
-      1.0f,
-  };
+  float mat[16] = {0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1};
 
   if (!PyArg_ParseTuple(args, "siO:Matrix.Shear", &plane, &matSize, &fac)) {
     return NULL;
@@ -1051,205 +1203,12 @@ static PyObject *C_Matrix_LocRotScale(PyObject *cls, PyObject *args)
   return Matrix_CreatePyObject(&mat[0][0], 4, 4, (PyTypeObject *)cls);
 }
 
-void matrix_as_3x3(float mat[3][3], MatrixObject *self)
-{
-  copy_v3_v3(mat[0], MATRIX_COL_PTR(self, 0));
-  copy_v3_v3(mat[1], MATRIX_COL_PTR(self, 1));
-  copy_v3_v3(mat[2], MATRIX_COL_PTR(self, 2));
-}
-
-static void matrix_copy(MatrixObject *mat_dst, const MatrixObject *mat_src)
-{
-  BLI_assert((mat_dst->col_num == mat_src->col_num) && (mat_dst->row_num == mat_src->row_num));
-  BLI_assert(mat_dst != mat_src);
+/** \} */
 
-  memcpy(mat_dst->matrix, mat_src->matrix, sizeof(float) * (mat_dst->col_num * mat_dst->row_num));
-}
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: To Quaternion
+ * \{ */
 
-static void matrix_unit_internal(MatrixObject *self)
-{
-  const int mat_size = sizeof(float) * (self->col_num * self->row_num);
-  memset(self->matrix, 0x0, mat_size);
-  const int col_row_max = min_ii(self->col_num, self->row_num);
-  const int row_num = self->row_num;
-  for (int col = 0; col < col_row_max; col++) {
-    self->matrix[(col * row_num) + col] = 1.0f;
-  }
-}
-
-/* transposes memory layout, rol/col's don't have to match */
-static void matrix_transpose_internal(float mat_dst_fl[], const MatrixObject *mat_src)
-{
-  ushort col, row;
-  uint i = 0;
-
-  for (row = 0; row < mat_src->row_num; row++) {
-    for (col = 0; col < mat_src->col_num; col++) {
-      mat_dst_fl[i++] = MATRIX_ITEM(mat_src, row, col);
-    }
-  }
-}
-
-/* assumes rowsize == colsize is checked and the read callback has run */
-static float matrix_determinant_internal(const MatrixObject *self)
-{
-  if (self->col_num == 2) {
-    return determinant_m2(MATRIX_ITEM(self, 0, 0),
-                          MATRIX_ITEM(self, 0, 1),
-                          MATRIX_ITEM(self, 1, 0),
-                          MATRIX_ITEM(self, 1, 1));
-  }
-  if (self->col_num == 3) {
-    return determinant_m3(MATRIX_ITEM(self, 0, 0),
-                          MATRIX_ITEM(self, 0, 1),
-                          MATRIX_ITEM(self, 0, 2),
-                          MATRIX_ITEM(self, 1, 0),
-                          MATRIX_ITEM(self, 1, 1),
-                          MATRIX_ITEM(self, 1, 2),
-                          MATRIX_ITEM(self, 2, 0),
-                          MATRIX_ITEM(self, 2, 1),
-                          MATRIX_ITEM(self, 2, 2));
-  }
-
-  return determinant_m4((const float(*)[4])self->matrix);
-}
-
-static void adjoint_matrix_n(float *mat_dst, const float *mat_src, const ushort dim)
-{
-  /* calculate the classical adjoint */
-  switch (dim) {
-    case 2: {
-      adjoint_m2_m2((float(*)[2])mat_dst, (const float(*)[2])mat_src);
-      break;
-    }
-    case 3: {
-      adjoint_m3_m3((float(*)[3])mat_dst, (const float(*)[3])mat_src);
-      break;
-    }
-    case 4: {
-      adjoint_m4_m4((float(*)[4])mat_dst, (const float(*)[4])mat_src);
-      break;
-    }
-    default:
-      BLI_assert_unreachable();
-      break;
-  }
-}
-
-static void matrix_invert_with_det_n_internal(float *mat_dst,
-                                              const float *mat_src,
-                                              const float det,
-                                              const ushort dim)
-{
-  float mat[MATRIX_MAX_DIM * MATRIX_MAX_DIM];
-  ushort i, j, k;
-
-  BLI_assert(det != 0.0f);
-
-  adjoint_matrix_n(mat, mat_src, dim);
-
-  /* divide by determinant & set values */
-  k = 0;
-  for (i = 0; i < dim; i++) {   /* col_num */
-    for (j = 0; j < dim; j++) { /* row_num */
-      mat_dst[MATRIX_ITEM_INDEX_NUMROW(dim, j, i)] = mat[k++] / det;
-    }
-  }
-}
-
-/**
- * \param r_mat: can be from `self->matrix` or not.
- */
-static bool matrix_invert_internal(const MatrixObject *self, float *r_mat)
-{
-  float det;
-  BLI_assert(self->col_num == self->row_num);
-  det = matrix_determinant_internal(self);
-
-  if (det != 0.0f) {
-    matrix_invert_with_det_n_internal(r_mat, self->matrix, det, self->col_num);
-    return true;
-  }
-
-  return false;
-}
-
-/**
- * Similar to `matrix_invert_internal` but should never error.
- * \param r_mat: can be from `self->matrix` or not.
- */
-static void matrix_invert_safe_internal(const MatrixObject *self, float *r_mat)
-{
-  float det;
-  float *in_mat = self->matrix;
-  BLI_assert(self->col_num == self->row_num);
-  det = matrix_determinant_internal(self);
-
-  if (det == 0.0f) {
-    const float eps = PSEUDOINVERSE_EPSILON;
-
-    /* We will copy self->matrix into r_mat (if needed),
-     * and modify it in place to add diagonal epsilon. */
-    in_mat = r_mat;
-
-    switch (self->col_num) {
-      case 2: {
-        float(*mat)[2] = (float(*)[2])in_mat;
-
-        if (in_mat != self->matrix) {
-          copy_m2_m2(mat, (const float(*)[2])self->matrix);
-        }
-        mat[0][0] += eps;
-        mat[1][1] += eps;
-
-        if (UNLIKELY((det = determinant_m2(mat[0][0], mat[0][1], mat[1][0], mat[1][1])) == 0.0f)) {
-          unit_m2(mat);
-          det = 1.0f;
-        }
-        break;
-      }
-      case 3: {
-        float(*mat)[3] = (float(*)[3])in_mat;
-
-        if (in_mat != self->matrix) {
-          copy_m3_m3(mat, (const float(*)[3])self->matrix);
-        }
-        mat[0][0] += eps;
-        mat[1][1] += eps;
-        mat[2][2] += eps;
-
-        if (UNLIKELY((det = determinant_m3_array(mat)) == 0.0f)) {
-          unit_m3(mat);
-          det = 1.0f;
-        }
-        break;
-      }
-      case 4: {
-        float(*mat)[4] = (float(*)[4])in_mat;
-
-        if (in_mat != self->matrix) {
-          copy_m4_m4(mat, (const float(*)[4])self->matrix);
-        }
-        mat[0][0] += eps;
-        mat[1][1] += eps;
-        mat[2][2] += eps;
-        mat[3][3] += eps;
-
-        if (UNLIKELY(det = determinant_m4(mat)) == 0.0f) {
-          unit_m4(mat);
-          det = 1.0f;
-        }
-        break;
-      }
-      default:
-        BLI_assert_unreachable();
-    }
-  }
-
-  matrix_invert_with_det_n_internal(r_mat, in_mat, det, self->col_num);
-}
-
-/*-----------------------------METHODS----------------------------*/
 PyDoc_STRVAR(Matrix_to_quaternion_doc,
              ".. method:: to_quaternion()\n"
              "\n"
@@ -1282,7 +1241,12 @@ static PyObject *Matrix_to_quaternion(MatrixObject *self)
   return Quaternion_CreatePyObject(quat, NULL);
 }
 
-/*---------------------------matrix.toEuler() --------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: To Euler
+ * \{ */
+
 PyDoc_STRVAR(Matrix_to_euler_doc,
              ".. method:: to_euler(order, euler_compat)\n"
              "\n"
@@ -1367,6 +1331,12 @@ static PyObject *Matrix_to_euler(MatrixObject *self, PyObject *args)
   return Euler_CreatePyObject(eul, order, NULL);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Resize
+ * \{ */
+
 PyDoc_STRVAR(Matrix_resize_4x4_doc,
              ".. method:: resize_4x4()\n"
              "\n"
@@ -1411,6 +1381,12 @@ static PyObject *Matrix_resize_4x4(MatrixObject *self)
   Py_RETURN_NONE;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: To NxN
+ * \{ */
+
 static PyObject *Matrix_to_NxN(MatrixObject *self, const int col_num, const int row_num)
 {
   const int mat_size = sizeof(float) * (col_num * row_num);
@@ -1480,6 +1456,12 @@ static PyObject *Matrix_to_4x4(MatrixObject *self)
   return Matrix_to_NxN(self, 4, 4);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: To Translation/Scale
+ * \{ */
+
 PyDoc_STRVAR(Matrix_to_translation_doc,
              ".. method:: to_translation()\n"
              "\n"
@@ -1539,9 +1521,13 @@ static PyObject *Matrix_to_scale(MatrixObject *self)
   return Vector_CreatePyObject(size, 3, NULL);
 }
 
-/*---------------------------matrix.invert() ---------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Invert
+ * \{ */
 
-/* re-usable checks for invert */
+/** Re-usable checks for invert. */
 static bool matrix_invert_is_compat(const MatrixObject *self)
 {
   if (self->col_num != self->row_num) {
@@ -1763,7 +1749,12 @@ static PyObject *Matrix_inverted_safe(MatrixObject *self)
   return Matrix_copy_notest(self, mat);
 }
 
-/*---------------------------matrix.adjugate() ---------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Adjugate
+ * \{ */
+
 PyDoc_STRVAR(
     Matrix_adjugate_doc,
     ".. method:: adjugate()\n"
@@ -1852,7 +1843,12 @@ static PyObject *Matrix_rotate(MatrixObject *self, PyObject *value)
   Py_RETURN_NONE;
 }
 
-/*---------------------------matrix.decompose() ---------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Decompose
+ * \{ */
+
 PyDoc_STRVAR(Matrix_decompose_doc,
              ".. method:: decompose()\n"
              "\n"
@@ -1890,6 +1886,12 @@ static PyObject *Matrix_decompose(MatrixObject *self)
   return ret;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Linear Interpolate (lerp)
+ * \{ */
+
 PyDoc_STRVAR(Matrix_lerp_doc,
              ".. function:: lerp(other, factor)\n"
              "\n"
@@ -1947,7 +1949,6 @@ static PyObject *Matrix_lerp(MatrixObject *self, PyObject *args)
   return Matrix_CreatePyObject(mat, self->col_num, self->row_num, Py_TYPE(self));
 }
 
-/*---------------------------matrix.determinant() ----------------*/
 PyDoc_STRVAR(
     Matrix_determinant_doc,
     ".. method:: determinant()\n"
@@ -1973,7 +1974,13 @@ static PyObject *Matrix_determinant(MatrixObject *self)
 
   return PyFloat_FromDouble((double)matrix_determinant_internal(self));
 }
-/*---------------------------matrix.transpose() ------------------*/
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Transpose
+ * \{ */
+
 PyDoc_STRVAR(
     Matrix_transpose_doc,
     ".. method:: transpose()\n"
@@ -2022,7 +2029,12 @@ static PyObject *Matrix_transposed(MatrixObject *self)
   return matrix__apply_to_copy(Matrix_transpose, self);
 }
 
-/*---------------------------matrix.normalize() ------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Normalize
+ * \{ */
+
 PyDoc_STRVAR(Matrix_normalize_doc,
              ".. method:: normalize()\n"
              "\n"
@@ -2068,7 +2080,12 @@ static PyObject *Matrix_normalized(MatrixObject *self)
   return matrix__apply_to_copy(Matrix_normalize, self);
 }
 
-/*---------------------------matrix.zero() -----------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Zero
+ * \{ */
+
 PyDoc_STRVAR(Matrix_zero_doc,
              ".. method:: zero()\n"
              "\n"
@@ -2089,7 +2106,13 @@ static PyObject *Matrix_zero(MatrixObject *self)
 
   Py_RETURN_NONE;
 }
-/*---------------------------matrix.identity(() ------------------*/
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Set Identity
+ * \{ */
+
 static void matrix_identity_internal(MatrixObject *self)
 {
   BLI_assert((self->col_num == self->row_num) && (self->row_num <= 4));
@@ -2137,8 +2160,13 @@ static PyObject *Matrix_identity(MatrixObject *self)
   Py_RETURN_NONE;
 }
 
-/*---------------------------Matrix.copy() ------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Methods: Copy
+ * \{ */
 
+/** Copy `Matrix.copy()` */
 static PyObject *Matrix_copy_notest(MatrixObject *self, const float *matrix)
 {
   return Matrix_CreatePyObject((const float *)matrix, self->col_num, self->row_num, Py_TYPE(self));
@@ -2159,6 +2187,8 @@ static PyObject *Matrix_copy(MatrixObject *self)
 
   return Matrix_copy_notest(self, self->matrix);
 }
+
+/** Deep-copy `Matrix.deepcopy()` */
 static PyObject *Matrix_deepcopy(MatrixObject *self, PyObject *args)
 {
   if (!PyC_CheckArgs_DeepCopy(args)) {
@@ -2167,8 +2197,12 @@ static PyObject *Matrix_deepcopy(MatrixObject *self, PyObject *args)
   return Matrix_copy(self);
 }
 
-/*----------------------------print object (internal)-------------*/
-/* print the object to screen */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: `__repr__` & `__str__`
+ * \{ */
+
 static PyObject *Matrix_repr(MatrixObject *self)
 {
   int col, row;
@@ -2257,6 +2291,12 @@ static PyObject *Matrix_str(MatrixObject *self)
 }
 #endif
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: Rich Compare
+ * \{ */
+
 static PyObject *Matrix_richcmpr(PyObject *a, PyObject *b, int op)
 {
   PyObject *res;
@@ -2298,6 +2338,12 @@ static PyObject *Matrix_richcmpr(PyObject *a, PyObject *b, int op)
   return Py_INCREF_RET(res);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: Hash (`__hash__`)
+ * \{ */
+
 static Py_hash_t Matrix_hash(MatrixObject *self)
 {
   float mat[MATRIX_MAX_DIM * MATRIX_MAX_DIM];
@@ -2315,16 +2361,22 @@ static Py_hash_t Matrix_hash(MatrixObject *self)
   return mathutils_array_hash(mat, self->row_num * self->col_num);
 }
 
-/*---------------------SEQUENCE PROTOCOLS------------------------
- * ----------------------------len(object)------------------------
- * sequence length */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: Sequence & Mapping Protocol Implementation
+ * \{ */
+
+/** Sequence length: `len(object)`. */
 static int Matrix_len(MatrixObject *self)
 {
   return self->row_num;
 }
-/*----------------------------object[]---------------------------
- * sequence accessor (get)
- * the wrapped vector gives direct access to the matrix data */
+
+/**
+ * Sequence accessor (get): `x = object[i]`.
+ * \note the wrapped vector gives direct access to the matrix data.
+ */
 static PyObject *Matrix_item_row(MatrixObject *self, int row)
 {
   if (BaseMath_ReadCallback_ForWrite(self) == -1) {
@@ -2340,7 +2392,10 @@ static PyObject *Matrix_item_row(MatrixObject *self, int row)
   return Vector_CreatePyObject_cb(
       (PyObject *)self, self->col_num, mathutils_matrix_row_cb_index, row);
 }
-/* same but column access */
+/**
+ * Sequence accessor (get): `x = object.col[i]`.
+ * \note the wrapped vector gives direct access to the matrix data.
+ */
 static PyObject *Matrix_item_col(MatrixObject *self, int col)
 {
   if (BaseMath_ReadCallback_ForWrite(self) == -1) {
@@ -2357,9 +2412,7 @@ static PyObject *Matrix_item_col(MatrixObject *self, int col)
       (PyObject *)self, self->row_num, mathutils_matrix_col_cb_index, col);
 }
 
-/*----------------------------object[]-------------------------
- * sequence accessor (set) */
-
+/** Sequence accessor (set): `object[i] = x`. */
 static int Matrix_ass_item_row(MatrixObject *self, int row, PyObject *value)
 {
   int col;
@@ -2386,6 +2439,8 @@ static int Matrix_ass_item_row(MatrixObject *self, int row, PyObject *value)
   (void)BaseMath_WriteCallback(self);
   return 0;
 }
+
+/** Sequence accessor (set): `object.col[i] = x`. */
 static int Matrix_ass_item_col(MatrixObject *self, int col, PyObject *value)
 {
   int row;
@@ -2413,8 +2468,7 @@ static int Matrix_ass_item_col(MatrixObject *self, int col, PyObject *value)
   return 0;
 }
 
-/*----------------------------object[z:y]------------------------
- * Sequence slice (get). */
+/** Sequence slice accessor (get): `x = object[i:j]`. */
 static PyObject *Matrix_slice(MatrixObject *self, int begin, int end)
 {
 
@@ -2439,8 +2493,8 @@ static PyObject *Matrix_slice(MatrixObject *self, int begin, int end)
 
   return tuple;
 }
-/*----------------------------object[z:y]------------------------
- * Sequence slice (set). */
+
+/** Sequence slice accessor (set): `object[i:j] = x`. */
 static int Matrix_ass_slice(MatrixObject *self, int begin, int end, PyObject *value)
 {
   PyObject *value_fast;
@@ -2500,8 +2554,84 @@ static int Matrix_ass_slice(MatrixObject *self, int begin, int end, PyObject *va
   (void)BaseMath_WriteCallback(self);
   return 0;
 }
-/*------------------------NUMERIC PROTOCOLS----------------------
- *------------------------obj + obj------------------------------*/
+
+/** Sequence generic subscript (get): `x = object[...]`. */
+static PyObject *Matrix_subscript(MatrixObject *self, PyObject *item)
+{
+  if (PyIndex_Check(item)) {
+    Py_ssize_t i;
+    i = PyNumber_AsSsize_t(item, PyExc_IndexError);
+    if (i == -1 && PyErr_Occurred()) {
+      return NULL;
+    }
+    if (i < 0) {
+      i += self->row_num;
+    }
+    return Matrix_item_row(self, i);
+  }
+  if (PySlice_Check(item)) {
+    Py_ssize_t start, stop, step, slicelength;
+
+    if (PySlice_GetIndicesEx(item, self->row_num, &start, &stop, &step, &slicelength) < 0) {
+      return NULL;
+    }
+
+    if (slicelength <= 0) {
+      return PyTuple_New(0);
+    }
+    if (step == 1) {
+      return Matrix_slice(self, start, stop);
+    }
+
+    PyErr_SetString(PyExc_IndexError, "slice steps not supported with matrices");
+    return NULL;
+  }
+
+  PyErr_Format(
+      PyExc_TypeError, "matrix indices must be integers, not %.200s", Py_TYPE(item)->tp_name);
+  return NULL;
+}
+
+/** Sequence generic subscript (set): `object[...] = x`. */
+static int Matrix_ass_subscript(MatrixObject *self, PyObject *item, PyObject *value)
+{
+  if (PyIndex_Check(item)) {
+    Py_ssize_t i = PyNumber_AsSsize_t(item, PyExc_IndexError);
+    if (i == -1 && PyErr_Occurred()) {
+      return -1;
+    }
+    if (i < 0) {
+      i += self->row_num;
+    }
+    return Matrix_ass_item_row(self, i, value);
+  }
+  if (PySlice_Check(item)) {
+    Py_ssize_t start, stop, step, slicelength;
+
+    if (PySlice_GetIndicesEx(item, self->row_num, &start, &stop, &step, &slicelength) < 0) {
+      return -1;
+    }
+
+    if (step == 1) {
+      return Matrix_ass_slice(self, start, stop, value);
+    }
+
+    PyErr_SetString(PyExc_IndexError, "slice steps not supported with matrices");
+    return -1;
+  }
+
+  PyErr_Format(
+      PyExc_TypeError, "matrix indices must be integers, not %.200s", Py_TYPE(item)->tp_name);
+  return -1;
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: Numeric Protocol Implementation
+ * \{ */
+
+/** Addition: `object + object`. */
 static PyObject *Matrix_add(PyObject *m1, PyObject *m2)
 {
   float mat[MATRIX_MAX_DIM * MATRIX_MAX_DIM];
@@ -2534,8 +2664,8 @@ static PyObject *Matrix_add(PyObject *m1, PyObject *m2)
 
   return Matrix_CreatePyObject(mat, mat1->col_num, mat1->row_num, Py_TYPE(mat1));
 }
-/*------------------------obj - obj------------------------------
- * subtraction */
+
+/** Subtraction: `object - object`. */
 static PyObject *Matrix_sub(PyObject *m1, PyObject *m2)
 {
   float mat[MATRIX_MAX_DIM * MATRIX_MAX_DIM];
@@ -2568,8 +2698,8 @@ static PyObject *Matrix_sub(PyObject *m1, PyObject *m2)
 
   return Matrix_CreatePyObject(mat, mat1->col_num, mat1->row_num, Py_TYPE(mat1));
 }
-/*------------------------obj * obj------------------------------
- * element-wise multiplication */
+
+/** Multiplication (element-wise): `object * object`. */
 static PyObject *matrix_mul_float(MatrixObject *mat, const float scalar)
 {
   float tmat[MATRIX_MAX_DIM * MATRIX_MAX_DIM];
@@ -2631,8 +2761,8 @@ static PyObject *Matrix_mul(PyObject *m1, PyObject *m2)
                Py_TYPE(m2)->tp_name);
   return NULL;
 }
-/*------------------------obj *= obj------------------------------
- * In place element-wise multiplication */
+
+/** Multiplication in-place (element-wise): `object *= object`. */
 static PyObject *Matrix_imul(PyObject *m1, PyObject *m2)
 {
   float scalar;
@@ -2680,8 +2810,8 @@ static PyObject *Matrix_imul(PyObject *m1, PyObject *m2)
   Py_INCREF(m1);
   return m1;
 }
-/*------------------------obj @ obj------------------------------
- * matrix multiplication */
+
+/** Multiplication (matrix multiply): `object @ object`. */
 static PyObject *Matrix_matmul(PyObject *m1, PyObject *m2)
 {
   int vec_num;
@@ -2756,8 +2886,8 @@ static PyObject *Matrix_matmul(PyObject *m1, PyObject *m2)
                Py_TYPE(m2)->tp_name);
   return NULL;
 }
-/*------------------------obj @= obj------------------------------
- * In place matrix multiplication */
+
+/** Multiplication in-place (matrix multiply): `object @= object`. */
 static PyObject *Matrix_imatmul(PyObject *m1, PyObject *m2)
 {
   MatrixObject *mat1 = NULL, *mat2 = NULL;
@@ -2816,87 +2946,24 @@ static PyObject *Matrix_imatmul(PyObject *m1, PyObject *m2)
   return m1;
 }
 
-/*-----------------PROTOCOL DECLARATIONS--------------------------*/
-static PySequenceMethods Matrix_SeqMethods = {
-    (lenfunc)Matrix_len,                  /* sq_length */
-    (binaryfunc)NULL,                     /* sq_concat */
-    (ssizeargfunc)NULL,                   /* sq_repeat */
-    (ssizeargfunc)Matrix_item_row,        /* sq_item */
-    (ssizessizeargfunc)NULL,              /* sq_slice, deprecated */
-    (ssizeobjargproc)Matrix_ass_item_row, /* sq_ass_item */
-    (ssizessizeobjargproc)NULL,           /* sq_ass_slice, deprecated */
-    (objobjproc)NULL,                     /* sq_contains */
-    (binaryfunc)NULL,                     /* sq_inplace_concat */
-    (ssizeargfunc)NULL,                   /* sq_inplace_repeat */
-};
-
-static PyObject *Matrix_subscript(MatrixObject *self, PyObject *item)
-{
-  if (PyIndex_Check(item)) {
-    Py_ssize_t i;
-    i = PyNumber_AsSsize_t(item, PyExc_IndexError);
-    if (i == -1 && PyErr_Occurred()) {
-      return NULL;
-    }
-    if (i < 0) {
-      i += self->row_num;
-    }
-    return Matrix_item_row(self, i);
-  }
-  if (PySlice_Check(item)) {
-    Py_ssize_t start, stop, step, slicelength;
-
-    if (PySlice_GetIndicesEx(item, self->row_num, &start, &stop, &step, &slicelength) < 0) {
-      return NULL;
-    }
-
-    if (slicelength <= 0) {
-      return PyTuple_New(0);
-    }
-    if (step == 1) {
-      return Matrix_slice(self, start, stop);
-    }
-
-    PyErr_SetString(PyExc_IndexError, "slice steps not supported with matrices");
-    return NULL;
-  }
-
-  PyErr_Format(
-      PyExc_TypeError, "matrix indices must be integers, not %.200s", Py_TYPE(item)->tp_name);
-  return NULL;
-}
-
-static int Matrix_ass_subscript(MatrixObject *self, PyObject *item, PyObject *value)
-{
-  if (PyIndex_Check(item)) {
-    Py_ssize_t i = PyNumber_AsSsize_t(item, PyExc_IndexError);
-    if (i == -1 && PyErr_Occurred()) {
-      return -1;
-    }
-    if (i < 0) {
-      i += self->row_num;
-    }
-    return Matrix_ass_item_row(self, i, value);
-  }
-  if (PySlice_Check(item)) {
-    Py_ssize_t start, stop, step, slicelength;
+/** \} */
 
-    if (PySlice_GetIndicesEx(item, self->row_num, &start, &stop, &step, &slicelength) < 0) {
-      return -1;
-    }
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: Protocol Declarations
+ * \{ */
 
-    if (step == 1) {
-      return Matrix_ass_slice(self, start, stop, value);
-    }
-
-    PyErr_SetString(PyExc_IndexError, "slice steps not supported with matrices");
-    return -1;
-  }
-
-  PyErr_Format(
-      PyExc_TypeError, "matrix indices must be integers, not %.200s", Py_TYPE(item)->tp_name);
-  return -1;
-}
+static PySequenceMethods Matrix_SeqMethods = {
+    (lenfunc)Matrix_len,                  /*sq_length*/
+    (binaryfunc)NULL,                     /*sq_concat*/
+    (ssizeargfunc)NULL,                   /*sq_repeat*/
+    (ssizeargfunc)Matrix_item_row,        /*sq_item*/
+    (ssizessizeargfunc)NULL,              /*sq_slice(DEPRECATED)*/
+    (ssizeobjargproc)Matrix_ass_item_row, /*sq_ass_item*/
+    (ssizessizeobjargproc)NULL,           /*sq_ass_slice(DEPRECATED)*/
+    (objobjproc)NULL,                     /*sq_contains*/
+    (binaryfunc)NULL,                     /*sq_inplace_concat*/
+    (ssizeargfunc)NULL,                   /*sq_inplace_repeat*/
+};
 
 static PyMappingMethods Matrix_AsMapping = {
     (lenfunc)Matrix_len,
@@ -2924,25 +2991,31 @@ static PyNumberMethods Matrix_NumMethods = {
     NULL,                              /*nb_int*/
     NULL,                              /*nb_reserved*/
     NULL,                              /*nb_float*/
-    NULL,                              /* nb_inplace_add */
-    NULL,                              /* nb_inplace_subtract */
-    (binaryfunc)Matrix_imul,           /* nb_inplace_multiply */
-    NULL,                              /* nb_inplace_remainder */
-    NULL,                              /* nb_inplace_power */
-    NULL,                              /* nb_inplace_lshift */
-    NULL,                              /* nb_inplace_rshift */
-    NULL,                              /* nb_inplace_and */
-    NULL,                              /* nb_inplace_xor */
-    NULL,                              /* nb_inplace_or */
-    NULL,                              /* nb_floor_divide */
-    NULL,                              /* nb_true_divide */
-    NULL,                              /* nb_inplace_floor_divide */
-    NULL,                              /* nb_inplace_true_divide */
-    NULL,                              /* nb_index */
-    (binaryfunc)Matrix_matmul,         /* nb_matrix_multiply */
-    (binaryfunc)Matrix_imatmul,        /* nb_inplace_matrix_multiply */
+    NULL,                              /*nb_inplace_add*/
+    NULL,                              /*nb_inplace_subtract*/
+    (binaryfunc)Matrix_imul,           /*nb_inplace_multiply*/
+    NULL,                              /*nb_inplace_remainder*/
+    NULL,                              /*nb_inplace_power*/
+    NULL,                              /*nb_inplace_lshift*/
+    NULL,                              /*nb_inplace_rshift*/
+    NULL,                              /*nb_inplace_and*/
+    NULL,                              /*nb_inplace_xor*/
+    NULL,                              /*nb_inplace_or*/
+    NULL,                              /*nb_floor_divide*/
+    NULL,                              /*nb_true_divide*/
+    NULL,                              /*nb_inplace_floor_divide*/
+    NULL,                              /*nb_inplace_true_divide*/
+    NULL,                              /*nb_index*/
+    (binaryfunc)Matrix_matmul,         /*nb_matrix_multiply*/
+    (binaryfunc)Matrix_imatmul,        /*nb_inplace_matrix_multiply*/
 };
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: Get/Set Item Implementation
+ * \{ */
+
 PyDoc_STRVAR(Matrix_translation_doc, "The translation component of the matrix.\n\n:type: Vector");
 static PyObject *Matrix_translation_get(MatrixObject *self, void *UNUSED(closure))
 {
@@ -3099,9 +3172,12 @@ static PyObject *Matrix_is_orthogonal_axis_vectors_get(MatrixObject *self, void
   return NULL;
 }
 
-/*****************************************************************************/
-/* Python attributes get/set structure:                                      */
-/*****************************************************************************/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: Get/Set Item Definitions
+ * \{ */
+
 static PyGetSetDef Matrix_getseters[] = {
     {"median_scale", (getter)Matrix_median_scale_get, (setter)NULL, Matrix_median_scale_doc, NULL},
     {"translation",
@@ -3141,7 +3217,12 @@ static PyGetSetDef Matrix_getseters[] = {
     {NULL, NULL, NULL, NULL, NULL} /* Sentinel */
 };
 
-/*-----------------------METHOD DEFINITIONS ----------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: Method Definitions
+ * \{ */
+
 static struct PyMethodDef Matrix_methods[] = {
     /* Derived values. */
     {"determinant", (PyCFunction)Matrix_determinant, METH_NOARGS, Matrix_determinant_doc},
@@ -3205,7 +3286,12 @@ static struct PyMethodDef Matrix_methods[] = {
     {NULL, NULL, 0, NULL},
 };
 
-/*------------------PY_OBECT DEFINITION--------------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: Python Object Definition
+ * \{ */
+
 PyDoc_STRVAR(
     matrix_doc,
     ".. class:: Matrix([rows])\n"
@@ -3268,6 +3354,12 @@ PyTypeObject matrix_Type = {
     NULL,                                                          /*tp_del*/
 };
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: C/API Constructors
+ * \{ */
+
 PyObject *Matrix_CreatePyObject(const float *mat,
                                 const ushort col_num,
                                 const ushort row_num,
@@ -3380,6 +3472,12 @@ PyObject *Matrix_CreatePyObject_alloc(float *mat,
   return (PyObject *)self;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix Type: C/API Parse Utilities
+ * \{ */
+
 /**
  * Use with PyArg_ParseTuple's "O&" formatting.
  */
@@ -3460,8 +3558,11 @@ int Matrix_Parse4x4(PyObject *o, void *p)
   return 1;
 }
 
-/* ----------------------------------------------------------------------------
- * special type for alternate access */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix-Access Type: Struct & Internal Functions
+ * \{ */
 
 typedef struct {
   PyObject_HEAD /* Required Python macro. */
@@ -3491,7 +3592,11 @@ static void MatrixAccess_dealloc(MatrixAccessObject *self)
   Py_TYPE(self)->tp_free(self);
 }
 
-/* sequence access */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix-Access Type: Sequence Protocol
+ * \{ */
 
 static int MatrixAccess_len(MatrixAccessObject *self)
 {
@@ -3609,13 +3714,13 @@ static int MatrixAccess_ass_subscript(MatrixAccessObject *self, PyObject *item,
 
 static PyObject *MatrixAccess_iter(MatrixAccessObject *self)
 {
-  /* Try get values from a collection */
+  /* Try get values from a collection. */
   PyObject *ret;
   PyObject *iter = NULL;
   ret = MatrixAccess_slice(self, 0, MATRIX_MAX_DIM);
 
-  /* we know this is a tuple so no need to PyIter_Check
-   * otherwise it could be NULL (unlikely) if conversion failed */
+  /* We know this is a tuple so no need to #PyIter_Check
+   * otherwise it could be NULL (unlikely) if conversion failed. */
   if (ret) {
     iter = PyObject_GetIter(ret);
     Py_DECREF(ret);
@@ -3630,6 +3735,12 @@ static PyMappingMethods MatrixAccess_AsMapping = {
     (objobjargproc)MatrixAccess_ass_subscript,
 };
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix-Access Type: Python Object Definition
+ * \{ */
+
 PyTypeObject matrix_access_Type = {
     PyVarObject_HEAD_INIT(NULL, 0) "MatrixAccess",           /*tp_name*/
     sizeof(MatrixAccessObject),                              /*tp_basicsize*/
@@ -3658,6 +3769,12 @@ PyTypeObject matrix_access_Type = {
     (getiterfunc)MatrixAccess_iter,                          /* getiterfunc tp_iter; */
 };
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Matrix-Access Type: C/API Constructor
+ * \{ */
+
 static PyObject *MatrixAccess_CreatePyObject(MatrixObject *matrix, const eMatrixAccess_t type)
 {
   MatrixAccessObject *matrix_access = (MatrixAccessObject *)PyObject_GC_New(MatrixObject,
@@ -3671,5 +3788,4 @@ static PyObject *MatrixAccess_CreatePyObject(MatrixObject *matrix, const eMatrix
   return (PyObject *)matrix_access;
 }
 
-/* end special access
- * -------------------------------------------------------------------------- */
+/** \} */
diff --git a/source/blender/python/mathutils/mathutils_Matrix.h b/source/blender/python/mathutils/mathutils_Matrix.h
index bc596ce6ac8..c8c207c13f3 100644
--- a/source/blender/python/mathutils/mathutils_Matrix.h
+++ b/source/blender/python/mathutils/mathutils_Matrix.h
@@ -45,7 +45,8 @@ typedef struct {
  * be stored in py_data) or be a wrapper for data allocated through
  * blender (stored in blend_data). This is an either/or struct not both */
 
-/* prototypes */
+/* Prototypes. */
+
 PyObject *Matrix_CreatePyObject(const float *mat,
                                 ushort col_num,
                                 ushort row_num,
@@ -70,6 +71,7 @@ PyObject *Matrix_CreatePyObject_alloc(float *mat,
                                       PyTypeObject *base_type) ATTR_WARN_UNUSED_RESULT;
 
 /* PyArg_ParseTuple's "O&" formatting helpers. */
+
 int Matrix_ParseAny(PyObject *o, void *p);
 int Matrix_Parse2x2(PyObject *o, void *p);
 int Matrix_Parse3x3(PyObject *o, void *p);
diff --git a/source/blender/python/mathutils/mathutils_Quaternion.c b/source/blender/python/mathutils/mathutils_Quaternion.c
index 7b51154f0d0..6994a313237 100644
--- a/source/blender/python/mathutils/mathutils_Quaternion.c
+++ b/source/blender/python/mathutils/mathutils_Quaternion.c
@@ -26,9 +26,49 @@ static void quat__axis_angle_sanitize(float axis[3], float *angle);
 static PyObject *Quaternion_copy(QuaternionObject *self);
 static PyObject *Quaternion_deepcopy(QuaternionObject *self, PyObject *args);
 
-/* -----------------------------METHODS------------------------------ */
+/* -------------------------------------------------------------------- */
+/** \name Utilities
+ * \{ */
 
-/* NOTE: BaseMath_ReadCallback must be called beforehand. */
+static PyObject *quat__apply_to_copy(PyObject *(*quat_func)(QuaternionObject *),
+                                     QuaternionObject *self)
+{
+  PyObject *ret = Quaternion_copy(self);
+  PyObject *ret_dummy = quat_func((QuaternionObject *)ret);
+  if (ret_dummy) {
+    Py_DECREF(ret_dummy);
+    return ret;
+  }
+  /* error */
+  Py_DECREF(ret);
+  return NULL;
+}
+
+/** Axis vector suffers from precision errors, use this function to ensure. */
+static void quat__axis_angle_sanitize(float axis[3], float *angle)
+{
+  if (axis) {
+    if (is_zero_v3(axis) || !isfinite(axis[0]) || !isfinite(axis[1]) || !isfinite(axis[2])) {
+      axis[0] = 1.0f;
+      axis[1] = 0.0f;
+      axis[2] = 0.0f;
+    }
+    else if (EXPP_FloatsAreEqual(axis[0], 0.0f, 10) && EXPP_FloatsAreEqual(axis[1], 0.0f, 10) &&
+             EXPP_FloatsAreEqual(axis[2], 0.0f, 10)) {
+      axis[0] = 1.0f;
+    }
+  }
+
+  if (angle) {
+    if (!isfinite(*angle)) {
+      *angle = 0.0f;
+    }
+  }
+}
+
+/**
+ * \note #BaseMath_ReadCallback must be called beforehand.
+ */
 static PyObject *Quaternion_to_tuple_ext(QuaternionObject *self, int ndigits)
 {
   PyObject *ret;
@@ -50,6 +90,72 @@ static PyObject *Quaternion_to_tuple_ext(QuaternionObject *self, int ndigits)
   return ret;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: `__new__` / `mathutils.Quaternion()`
+ * \{ */
+
+static PyObject *Quaternion_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
+{
+  PyObject *seq = NULL;
+  double angle = 0.0f;
+  float quat[QUAT_SIZE];
+  unit_qt(quat);
+
+  if (kwds && PyDict_Size(kwds)) {
+    PyErr_SetString(PyExc_TypeError,
+                    "mathutils.Quaternion(): "
+                    "takes no keyword args");
+    return NULL;
+  }
+
+  if (!PyArg_ParseTuple(args, "|Od:mathutils.Quaternion", &seq, &angle)) {
+    return NULL;
+  }
+
+  switch (PyTuple_GET_SIZE(args)) {
+    case 0:
+      break;
+    case 1: {
+      int size;
+
+      if ((size = mathutils_array_parse(quat, 3, QUAT_SIZE, seq, "mathutils.Quaternion()")) ==
+          -1) {
+        return NULL;
+      }
+
+      if (size == 4) {
+        /* 4d: Quaternion (common case) */
+      }
+      else {
+        /* 3d: Interpret as exponential map */
+        BLI_assert(size == 3);
+        expmap_to_quat(quat, quat);
+      }
+
+      break;
+    }
+    case 2: {
+      float axis[3];
+      if (mathutils_array_parse(axis, 3, 3, seq, "mathutils.Quaternion()") == -1) {
+        return NULL;
+      }
+      angle = angle_wrap_rad(angle); /* clamp because of precision issues */
+      axis_angle_to_quat(quat, axis, angle);
+      break;
+      /* PyArg_ParseTuple assures no more than 2 */
+    }
+  }
+  return Quaternion_CreatePyObject(quat, type);
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: To Euler
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_to_euler_doc,
              ".. method:: to_euler(order, euler_compat)\n"
              "\n"
@@ -114,6 +220,12 @@ static PyObject *Quaternion_to_euler(QuaternionObject *self, PyObject *args)
   return Euler_CreatePyObject(eul, order, NULL);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: To Matrix
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_to_matrix_doc,
              ".. method:: to_matrix()\n"
              "\n"
@@ -133,6 +245,12 @@ static PyObject *Quaternion_to_matrix(QuaternionObject *self)
   return Matrix_CreatePyObject(mat, 3, 3, NULL);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: To Axis/Angle
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_to_axis_angle_doc,
              ".. method:: to_axis_angle()\n"
              "\n"
@@ -163,6 +281,12 @@ static PyObject *Quaternion_to_axis_angle(QuaternionObject *self)
   return ret;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: To Swing/Twist
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_to_swing_twist_doc,
              ".. method:: to_swing_twist(axis)\n"
              "\n"
@@ -207,6 +331,12 @@ static PyObject *Quaternion_to_swing_twist(QuaternionObject *self, PyObject *axi
   return ret;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: To Exponential Map
+ * \{ */
+
 PyDoc_STRVAR(
     Quaternion_to_exponential_map_doc,
     ".. method:: to_exponential_map()\n"
@@ -232,6 +362,12 @@ static PyObject *Quaternion_to_exponential_map(QuaternionObject *self)
   return Vector_CreatePyObject(expmap, 3, NULL);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Cross Product
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_cross_doc,
              ".. method:: cross(other)\n"
              "\n"
@@ -259,6 +395,12 @@ static PyObject *Quaternion_cross(QuaternionObject *self, PyObject *value)
   return Quaternion_CreatePyObject(quat, Py_TYPE(self));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Dot Product
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_dot_doc,
              ".. method:: dot(other)\n"
              "\n"
@@ -285,6 +427,12 @@ static PyObject *Quaternion_dot(QuaternionObject *self, PyObject *value)
   return PyFloat_FromDouble(dot_qtqt(self->quat, tquat));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Rotation Difference
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_rotation_difference_doc,
              ".. function:: rotation_difference(other)\n"
              "\n"
@@ -315,6 +463,12 @@ static PyObject *Quaternion_rotation_difference(QuaternionObject *self, PyObject
   return Quaternion_CreatePyObject(quat, Py_TYPE(self));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Spherical Interpolation (slerp)
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_slerp_doc,
              ".. function:: slerp(other, factor)\n"
              "\n"
@@ -360,6 +514,12 @@ static PyObject *Quaternion_slerp(QuaternionObject *self, PyObject *args)
   return Quaternion_CreatePyObject(quat, Py_TYPE(self));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Rotate
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_rotate_doc,
              ".. method:: rotate(other)\n"
              "\n"
@@ -423,9 +583,15 @@ static PyObject *Quaternion_make_compatible(QuaternionObject *self, PyObject *va
   Py_RETURN_NONE;
 }
 
-/* ----------------------------Quaternion.normalize()---------------- */
-/* Normalize the quaternion. This may change the angle as well as the
- * rotation axis, as all of (w, x, y, z) are scaled. */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Normalize
+ *
+ * Normalize the quaternion. This may change the angle as well as the
+ * rotation axis, as all of (w, x, y, z) are scaled.
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_normalize_doc,
              ".. function:: normalize()\n"
              "\n"
@@ -453,6 +619,15 @@ static PyObject *Quaternion_normalized(QuaternionObject *self)
   return quat__apply_to_copy(Quaternion_normalize, self);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Invert
+ *
+ * Normalize the quaternion. This may change the angle as well as the
+ * rotation axis, as all of (w, x, y, z) are scaled.
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_invert_doc,
              ".. function:: invert()\n"
              "\n"
@@ -480,6 +655,12 @@ static PyObject *Quaternion_inverted(QuaternionObject *self)
   return quat__apply_to_copy(Quaternion_invert, self);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Set Identity
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_identity_doc,
              ".. function:: identity()\n"
              "\n"
@@ -498,6 +679,12 @@ static PyObject *Quaternion_identity(QuaternionObject *self)
   Py_RETURN_NONE;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Negate
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_negate_doc,
              ".. function:: negate()\n"
              "\n"
@@ -516,6 +703,12 @@ static PyObject *Quaternion_negate(QuaternionObject *self)
   Py_RETURN_NONE;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Conjugate
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_conjugate_doc,
              ".. function:: conjugate()\n"
              "\n"
@@ -543,6 +736,12 @@ static PyObject *Quaternion_conjugated(QuaternionObject *self)
   return quat__apply_to_copy(Quaternion_conjugate, self);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Methods: Copy/Deep-Copy
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_copy_doc,
              ".. function:: copy()\n"
              "\n"
@@ -569,7 +768,12 @@ static PyObject *Quaternion_deepcopy(QuaternionObject *self, PyObject *args)
   return Quaternion_copy(self);
 }
 
-/* print the object to screen */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: `__repr__` & `__str__`
+ * \{ */
+
 static PyObject *Quaternion_repr(QuaternionObject *self)
 {
   PyObject *ret, *tuple;
@@ -608,6 +812,12 @@ static PyObject *Quaternion_str(QuaternionObject *self)
 }
 #endif
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: Rich Compare
+ * \{ */
+
 static PyObject *Quaternion_richcmpr(PyObject *a, PyObject *b, int op)
 {
   PyObject *res;
@@ -646,6 +856,12 @@ static PyObject *Quaternion_richcmpr(PyObject *a, PyObject *b, int op)
   return Py_INCREF_RET(res);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: Hash (`__hash__`)
+ * \{ */
+
 static Py_hash_t Quaternion_hash(QuaternionObject *self)
 {
   if (BaseMath_ReadCallback(self) == -1) {
@@ -659,15 +875,19 @@ static Py_hash_t Quaternion_hash(QuaternionObject *self)
   return mathutils_array_hash(self->quat, QUAT_SIZE);
 }
 
-/* ---------------------SEQUENCE PROTOCOLS------------------------ */
-/* ----------------------------len(object)------------------------ */
-/* sequence length */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: Sequence & Mapping Protocols Implementation
+ * \{ */
+
+/** Sequence length: `len(object)`. */
 static int Quaternion_len(QuaternionObject *UNUSED(self))
 {
   return QUAT_SIZE;
 }
-/* ----------------------------object[]--------------------------- */
-/* sequence accessor (get) */
+
+/** Sequence accessor (get): `x = object[i]`. */
 static PyObject *Quaternion_item(QuaternionObject *self, int i)
 {
   if (i < 0) {
@@ -687,8 +907,8 @@ static PyObject *Quaternion_item(QuaternionObject *self, int i)
 
   return PyFloat_FromDouble(self->quat[i]);
 }
-/* ----------------------------object[]------------------------- */
-/* sequence accessor (set) */
+
+/** Sequence accessor (set): `object[i] = x`. */
 static int Quaternion_ass_item(QuaternionObject *self, int i, PyObject *ob)
 {
   float f;
@@ -724,8 +944,8 @@ static int Quaternion_ass_item(QuaternionObject *self, int i, PyObject *ob)
 
   return 0;
 }
-/* ----------------------------object[z:y]------------------------ */
-/* sequence slice (get) */
+
+/** Sequence slice accessor (get): `x = object[i:j]`. */
 static PyObject *Quaternion_slice(QuaternionObject *self, int begin, int end)
 {
   PyObject *tuple;
@@ -749,8 +969,8 @@ static PyObject *Quaternion_slice(QuaternionObject *self, int begin, int end)
 
   return tuple;
 }
-/* ----------------------------object[z:y]------------------------ */
-/* sequence slice (set) */
+
+/** Sequence slice accessor (set): `object[i:j] = x`. */
 static int Quaternion_ass_slice(QuaternionObject *self, int begin, int end, PyObject *seq)
 {
   int i, size;
@@ -779,7 +999,7 @@ static int Quaternion_ass_slice(QuaternionObject *self, int begin, int end, PyOb
     return -1;
   }
 
-  /* parsed well - now set in vector */
+  /* Parsed well, now set in vector. */
   for (i = 0; i < size; i++) {
     self->quat[begin + i] = quat[i];
   }
@@ -788,6 +1008,7 @@ static int Quaternion_ass_slice(QuaternionObject *self, int begin, int end, PyOb
   return 0;
 }
 
+/** Sequence generic subscript (get): `x = object[...]`. */
 static PyObject *Quaternion_subscript(QuaternionObject *self, PyObject *item)
 {
   if (PyIndex_Check(item)) {
@@ -824,6 +1045,7 @@ static PyObject *Quaternion_subscript(QuaternionObject *self, PyObject *item)
   return NULL;
 }
 
+/** Sequence generic subscript (set): `object[...] = x`. */
 static int Quaternion_ass_subscript(QuaternionObject *self, PyObject *item, PyObject *value)
 {
   if (PyIndex_Check(item)) {
@@ -856,9 +1078,13 @@ static int Quaternion_ass_subscript(QuaternionObject *self, PyObject *item, PyOb
   return -1;
 }
 
-/* ------------------------NUMERIC PROTOCOLS---------------------- */
-/* ------------------------obj + obj------------------------------ */
-/* addition */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: Numeric Protocol Implementation
+ * \{ */
+
+/** Addition: `object + object`. */
 static PyObject *Quaternion_add(PyObject *q1, PyObject *q2)
 {
   float quat[QUAT_SIZE];
@@ -882,8 +1108,8 @@ static PyObject *Quaternion_add(PyObject *q1, PyObject *q2)
   add_qt_qtqt(quat, quat1->quat, quat2->quat, 1.0f);
   return Quaternion_CreatePyObject(quat, Py_TYPE(q1));
 }
-/* ------------------------obj - obj------------------------------ */
-/* subtraction */
+
+/** Subtraction: `object - object`. */
 static PyObject *Quaternion_sub(PyObject *q1, PyObject *q2)
 {
   int x;
@@ -921,8 +1147,7 @@ static PyObject *quat_mul_float(QuaternionObject *quat, const float scalar)
   return Quaternion_CreatePyObject(tquat, Py_TYPE(quat));
 }
 
-/*------------------------obj * obj------------------------------
- * multiplication */
+/** Multiplication (element-wise or scalar): `object * object`. */
 static PyObject *Quaternion_mul(PyObject *q1, PyObject *q2)
 {
   float scalar;
@@ -965,8 +1190,8 @@ static PyObject *Quaternion_mul(PyObject *q1, PyObject *q2)
                Py_TYPE(q2)->tp_name);
   return NULL;
 }
-/*------------------------obj *= obj------------------------------
- * in-place multiplication */
+
+/** Multiplication in-place (element-wise or scalar): `object *= object`. */
 static PyObject *Quaternion_imul(PyObject *q1, PyObject *q2)
 {
   float scalar;
@@ -1005,8 +1230,8 @@ static PyObject *Quaternion_imul(PyObject *q1, PyObject *q2)
   Py_INCREF(q1);
   return q1;
 }
-/*------------------------obj @ obj------------------------------
- * quaternion multiplication */
+
+/** Multiplication (quaternion multiply): `object @ object`. */
 static PyObject *Quaternion_matmul(PyObject *q1, PyObject *q2)
 {
   float quat[QUAT_SIZE];
@@ -1060,8 +1285,8 @@ static PyObject *Quaternion_matmul(PyObject *q1, PyObject *q2)
                Py_TYPE(q2)->tp_name);
   return NULL;
 }
-/*------------------------obj @= obj------------------------------
- * in-place quaternion multiplication */
+
+/** Multiplication in-place (quaternion multiply): `object @= object`. */
 static PyObject *Quaternion_imatmul(PyObject *q1, PyObject *q2)
 {
   float quat[QUAT_SIZE];
@@ -1098,8 +1323,7 @@ static PyObject *Quaternion_imatmul(PyObject *q1, PyObject *q2)
   return q1;
 }
 
-/* -obj
- * Returns the negative of this object. */
+/** Negative (returns the negative of this object): `-object`. */
 static PyObject *Quaternion_neg(QuaternionObject *self)
 {
   float tquat[QUAT_SIZE];
@@ -1112,18 +1336,23 @@ static PyObject *Quaternion_neg(QuaternionObject *self)
   return Quaternion_CreatePyObject(tquat, Py_TYPE(self));
 }
 
-/* -----------------PROTOCOL DECLARATIONS-------------------------- */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: Protocol Declarations
+ * \{ */
+
 static PySequenceMethods Quaternion_SeqMethods = {
-    (lenfunc)Quaternion_len,              /* sq_length */
-    (binaryfunc)NULL,                     /* sq_concat */
-    (ssizeargfunc)NULL,                   /* sq_repeat */
-    (ssizeargfunc)Quaternion_item,        /* sq_item */
-    (ssizessizeargfunc)NULL,              /* sq_slice, deprecated */
-    (ssizeobjargproc)Quaternion_ass_item, /* sq_ass_item */
-    (ssizessizeobjargproc)NULL,           /* sq_ass_slice, deprecated */
-    (objobjproc)NULL,                     /* sq_contains */
-    (binaryfunc)NULL,                     /* sq_inplace_concat */
-    (ssizeargfunc)NULL,                   /* sq_inplace_repeat */
+    (lenfunc)Quaternion_len,              /*sq_length*/
+    (binaryfunc)NULL,                     /*sq_concat*/
+    (ssizeargfunc)NULL,                   /*sq_repeat*/
+    (ssizeargfunc)Quaternion_item,        /*sq_item*/
+    (ssizessizeargfunc)NULL,              /*sq_slice(deprecated)*/
+    (ssizeobjargproc)Quaternion_ass_item, /*sq_ass_item*/
+    (ssizessizeobjargproc)NULL,           /*sq_ass_slice(deprecated)*/
+    (objobjproc)NULL,                     /*sq_contains*/
+    (binaryfunc)NULL,                     /*sq_inplace_concat*/
+    (ssizeargfunc)NULL,                   /*sq_inplace_repeat*/
 };
 
 static PyMappingMethods Quaternion_AsMapping = {
@@ -1152,25 +1381,31 @@ static PyNumberMethods Quaternion_NumMethods = {
     NULL,                           /*nb_int*/
     NULL,                           /*nb_reserved*/
     NULL,                           /*nb_float*/
-    NULL,                           /* nb_inplace_add */
-    NULL,                           /* nb_inplace_subtract */
-    (binaryfunc)Quaternion_imul,    /* nb_inplace_multiply */
-    NULL,                           /* nb_inplace_remainder */
-    NULL,                           /* nb_inplace_power */
-    NULL,                           /* nb_inplace_lshift */
-    NULL,                           /* nb_inplace_rshift */
-    NULL,                           /* nb_inplace_and */
-    NULL,                           /* nb_inplace_xor */
-    NULL,                           /* nb_inplace_or */
-    NULL,                           /* nb_floor_divide */
-    NULL,                           /* nb_true_divide */
-    NULL,                           /* nb_inplace_floor_divide */
-    NULL,                           /* nb_inplace_true_divide */
-    NULL,                           /* nb_index */
-    (binaryfunc)Quaternion_matmul,  /* nb_matrix_multiply */
-    (binaryfunc)Quaternion_imatmul, /* nb_inplace_matrix_multiply */
+    NULL,                           /*nb_inplace_add*/
+    NULL,                           /*nb_inplace_subtract*/
+    (binaryfunc)Quaternion_imul,    /*nb_inplace_multiply*/
+    NULL,                           /*nb_inplace_remainder*/
+    NULL,                           /*nb_inplace_power*/
+    NULL,                           /*nb_inplace_lshift*/
+    NULL,                           /*nb_inplace_rshift*/
+    NULL,                           /*nb_inplace_and*/
+    NULL,                           /*nb_inplace_xor*/
+    NULL,                           /*nb_inplace_or*/
+    NULL,                           /*nb_floor_divide*/
+    NULL,                           /*nb_true_divide*/
+    NULL,                           /*nb_inplace_floor_divide*/
+    NULL,                           /*nb_inplace_true_divide*/
+    NULL,                           /*nb_index*/
+    (binaryfunc)Quaternion_matmul,  /*nb_matrix_multiply*/
+    (binaryfunc)Quaternion_imatmul, /*nb_inplace_matrix_multiply*/
 };
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: Get/Set Item Implementation
+ * \{ */
+
 PyDoc_STRVAR(Quaternion_axis_doc, "Quaternion axis value.\n\n:type: float");
 static PyObject *Quaternion_axis_get(QuaternionObject *self, void *type)
 {
@@ -1300,98 +1535,69 @@ static int Quaternion_axis_vector_set(QuaternionObject *self,
   return 0;
 }
 
-/* ----------------------------------mathutils.Quaternion() -------------- */
-static PyObject *Quaternion_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
-{
-  PyObject *seq = NULL;
-  double angle = 0.0f;
-  float quat[QUAT_SIZE];
-  unit_qt(quat);
-
-  if (kwds && PyDict_Size(kwds)) {
-    PyErr_SetString(PyExc_TypeError,
-                    "mathutils.Quaternion(): "
-                    "takes no keyword args");
-    return NULL;
-  }
-
-  if (!PyArg_ParseTuple(args, "|Od:mathutils.Quaternion", &seq, &angle)) {
-    return NULL;
-  }
-
-  switch (PyTuple_GET_SIZE(args)) {
-    case 0:
-      break;
-    case 1: {
-      int size;
-
-      if ((size = mathutils_array_parse(quat, 3, QUAT_SIZE, seq, "mathutils.Quaternion()")) ==
-          -1) {
-        return NULL;
-      }
+/** \} */
 
-      if (size == 4) {
-        /* 4d: Quaternion (common case) */
-      }
-      else {
-        /* 3d: Interpret as exponential map */
-        BLI_assert(size == 3);
-        expmap_to_quat(quat, quat);
-      }
-
-      break;
-    }
-    case 2: {
-      float axis[3];
-      if (mathutils_array_parse(axis, 3, 3, seq, "mathutils.Quaternion()") == -1) {
-        return NULL;
-      }
-      angle = angle_wrap_rad(angle); /* clamp because of precision issues */
-      axis_angle_to_quat(quat, axis, angle);
-      break;
-      /* PyArg_ParseTuple assures no more than 2 */
-    }
-  }
-  return Quaternion_CreatePyObject(quat, type);
-}
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: Get/Set Item Definitions
+ * \{ */
 
-static PyObject *quat__apply_to_copy(PyObject *(*quat_func)(QuaternionObject *),
-                                     QuaternionObject *self)
-{
-  PyObject *ret = Quaternion_copy(self);
-  PyObject *ret_dummy = quat_func((QuaternionObject *)ret);
-  if (ret_dummy) {
-    Py_DECREF(ret_dummy);
-    return ret;
-  }
-  /* error */
-  Py_DECREF(ret);
-  return NULL;
-}
+static PyGetSetDef Quaternion_getseters[] = {
+    {"w",
+     (getter)Quaternion_axis_get,
+     (setter)Quaternion_axis_set,
+     Quaternion_axis_doc,
+     (void *)0},
+    {"x",
+     (getter)Quaternion_axis_get,
+     (setter)Quaternion_axis_set,
+     Quaternion_axis_doc,
+     (void *)1},
+    {"y",
+     (getter)Quaternion_axis_get,
+     (setter)Quaternion_axis_set,
+     Quaternion_axis_doc,
+     (void *)2},
+    {"z",
+     (getter)Quaternion_axis_get,
+     (setter)Quaternion_axis_set,
+     Quaternion_axis_doc,
+     (void *)3},
+    {"magnitude", (getter)Quaternion_magnitude_get, (setter)NULL, Quaternion_magnitude_doc, NULL},
+    {"angle",
+     (getter)Quaternion_angle_get,
+     (setter)Quaternion_angle_set,
+     Quaternion_angle_doc,
+     NULL},
+    {"axis",
+     (getter)Quaternion_axis_vector_get,
+     (setter)Quaternion_axis_vector_set,
+     Quaternion_axis_vector_doc,
+     NULL},
+    {"is_wrapped",
+     (getter)BaseMathObject_is_wrapped_get,
+     (setter)NULL,
+     BaseMathObject_is_wrapped_doc,
+     NULL},
+    {"is_frozen",
+     (getter)BaseMathObject_is_frozen_get,
+     (setter)NULL,
+     BaseMathObject_is_frozen_doc,
+     NULL},
+    {"is_valid",
+     (getter)BaseMathObject_is_valid_get,
+     (setter)NULL,
+     BaseMathObject_is_valid_doc,
+     NULL},
+    {"owner", (getter)BaseMathObject_owner_get, (setter)NULL, BaseMathObject_owner_doc, NULL},
+    {NULL, NULL, NULL, NULL, NULL} /* Sentinel */
+};
 
-/* axis vector suffers from precision errors, use this function to ensure */
-static void quat__axis_angle_sanitize(float axis[3], float *angle)
-{
-  if (axis) {
-    if (is_zero_v3(axis) || !isfinite(axis[0]) || !isfinite(axis[1]) || !isfinite(axis[2])) {
-      axis[0] = 1.0f;
-      axis[1] = 0.0f;
-      axis[2] = 0.0f;
-    }
-    else if (EXPP_FloatsAreEqual(axis[0], 0.0f, 10) && EXPP_FloatsAreEqual(axis[1], 0.0f, 10) &&
-             EXPP_FloatsAreEqual(axis[2], 0.0f, 10)) {
-      axis[0] = 1.0f;
-    }
-  }
+/** \} */
 
-  if (angle) {
-    if (!isfinite(*angle)) {
-      *angle = 0.0f;
-    }
-  }
-}
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: Method Definitions
+ * \{ */
 
-/* -----------------------METHOD DEFINITIONS ---------------------- */
 static struct PyMethodDef Quaternion_methods[] = {
     /* In place only. */
     {"identity", (PyCFunction)Quaternion_identity, METH_NOARGS, Quaternion_identity_doc},
@@ -1446,61 +1652,12 @@ static struct PyMethodDef Quaternion_methods[] = {
     {NULL, NULL, 0, NULL},
 };
 
-/*****************************************************************************/
-/* Python attributes get/set structure:                                      */
-/*****************************************************************************/
-static PyGetSetDef Quaternion_getseters[] = {
-    {"w",
-     (getter)Quaternion_axis_get,
-     (setter)Quaternion_axis_set,
-     Quaternion_axis_doc,
-     (void *)0},
-    {"x",
-     (getter)Quaternion_axis_get,
-     (setter)Quaternion_axis_set,
-     Quaternion_axis_doc,
-     (void *)1},
-    {"y",
-     (getter)Quaternion_axis_get,
-     (setter)Quaternion_axis_set,
-     Quaternion_axis_doc,
-     (void *)2},
-    {"z",
-     (getter)Quaternion_axis_get,
-     (setter)Quaternion_axis_set,
-     Quaternion_axis_doc,
-     (void *)3},
-    {"magnitude", (getter)Quaternion_magnitude_get, (setter)NULL, Quaternion_magnitude_doc, NULL},
-    {"angle",
-     (getter)Quaternion_angle_get,
-     (setter)Quaternion_angle_set,
-     Quaternion_angle_doc,
-     NULL},
-    {"axis",
-     (getter)Quaternion_axis_vector_get,
-     (setter)Quaternion_axis_vector_set,
-     Quaternion_axis_vector_doc,
-     NULL},
-    {"is_wrapped",
-     (getter)BaseMathObject_is_wrapped_get,
-     (setter)NULL,
-     BaseMathObject_is_wrapped_doc,
-     NULL},
-    {"is_frozen",
-     (getter)BaseMathObject_is_frozen_get,
-     (setter)NULL,
-     BaseMathObject_is_frozen_doc,
-     NULL},
-    {"is_valid",
-     (getter)BaseMathObject_is_valid_get,
-     (setter)NULL,
-     BaseMathObject_is_valid_doc,
-     NULL},
-    {"owner", (getter)BaseMathObject_owner_get, (setter)NULL, BaseMathObject_owner_doc, NULL},
-    {NULL, NULL, NULL, NULL, NULL} /* Sentinel */
-};
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: Python Object Definition
+ * \{ */
 
-/* ------------------PY_OBECT DEFINITION-------------------------- */
 PyDoc_STRVAR(quaternion_doc,
              ".. class:: Quaternion([seq, [angle]])\n"
              "\n"
@@ -1577,6 +1734,12 @@ PyTypeObject quaternion_Type = {
     NULL,                                                          /* tp_del */
 };
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Quaternion Type: C/API Constructors
+ * \{ */
+
 PyObject *Quaternion_CreatePyObject(const float quat[4], PyTypeObject *base_type)
 {
   QuaternionObject *self;
@@ -1643,3 +1806,5 @@ PyObject *Quaternion_CreatePyObject_cb(PyObject *cb_user, uchar cb_type, uchar c
 
   return (PyObject *)self;
 }
+
+/** \} */
diff --git a/source/blender/python/mathutils/mathutils_Quaternion.h b/source/blender/python/mathutils/mathutils_Quaternion.h
index 96e3d2e0055..c45b0a98a7b 100644
--- a/source/blender/python/mathutils/mathutils_Quaternion.h
+++ b/source/blender/python/mathutils/mathutils_Quaternion.h
@@ -20,7 +20,8 @@ typedef struct {
  * be stored in py_data) or be a wrapper for data allocated through
  * blender (stored in blend_data). This is an either/or struct not both */
 
-/* prototypes */
+/* Prototypes. */
+
 PyObject *Quaternion_CreatePyObject(const float quat[4],
                                     PyTypeObject *base_type) ATTR_WARN_UNUSED_RESULT;
 PyObject *Quaternion_CreatePyObject_wrap(float quat[4],
diff --git a/source/blender/python/mathutils/mathutils_Vector.c b/source/blender/python/mathutils/mathutils_Vector.c
index ffeb121b3d1..0c9cbd6ccfa 100644
--- a/source/blender/python/mathutils/mathutils_Vector.c
+++ b/source/blender/python/mathutils/mathutils_Vector.c
@@ -24,19 +24,108 @@
  */
 #define MAX_DIMENSIONS 4
 
-/* Swizzle axes get packed into a single value that is used as a closure. Each
+/**
+ * Swizzle axes get packed into a single value that is used as a closure. Each
  * axis uses SWIZZLE_BITS_PER_AXIS bits. The first bit (SWIZZLE_VALID_AXIS) is
- * used as a sentinel: if it is unset, the axis is not valid. */
+ * used as a sentinel: if it is unset, the axis is not valid.
+ */
 #define SWIZZLE_BITS_PER_AXIS 3
 #define SWIZZLE_VALID_AXIS 0x4
 #define SWIZZLE_AXIS 0x3
 
 static PyObject *Vector_copy(VectorObject *self);
 static PyObject *Vector_deepcopy(VectorObject *self, PyObject *args);
-static PyObject *Vector_to_tuple_ex(VectorObject *self, int ndigits);
+
+/* -------------------------------------------------------------------- */
+/** \name Utilities
+ * \{ */
+
+/**
+ * Row vector multiplication - (Vector * Matrix)
+ * <pre>
+ * [x][y][z] * [1][4][7]
+ *             [2][5][8]
+ *             [3][6][9]
+ * </pre>
+ * \note vector/matrix multiplication is not commutative.
+ */
 static int row_vector_multiplication(float r_vec[MAX_DIMENSIONS],
                                      VectorObject *vec,
-                                     MatrixObject *mat);
+                                     MatrixObject *mat)
+{
+  float vec_cpy[MAX_DIMENSIONS];
+  int row, col, z = 0, vec_num = vec->vec_num;
+
+  if (mat->row_num != vec_num) {
+    if (mat->row_num == 4 && vec_num == 3) {
+      vec_cpy[3] = 1.0f;
+    }
+    else {
+      PyErr_SetString(PyExc_ValueError,
+                      "vector * matrix: matrix column size "
+                      "and the vector size must be the same");
+      return -1;
+    }
+  }
+
+  if (BaseMath_ReadCallback(vec) == -1 || BaseMath_ReadCallback(mat) == -1) {
+    return -1;
+  }
+
+  memcpy(vec_cpy, vec->vec, vec_num * sizeof(float));
+
+  r_vec[3] = 1.0f;
+  /* Multiplication. */
+  for (col = 0; col < mat->col_num; col++) {
+    double dot = 0.0;
+    for (row = 0; row < mat->row_num; row++) {
+      dot += (double)(MATRIX_ITEM(mat, row, col) * vec_cpy[row]);
+    }
+    r_vec[z++] = (float)dot;
+  }
+  return 0;
+}
+
+static PyObject *vec__apply_to_copy(PyObject *(*vec_func)(VectorObject *), VectorObject *self)
+{
+  PyObject *ret = Vector_copy(self);
+  PyObject *ret_dummy = vec_func((VectorObject *)ret);
+  if (ret_dummy) {
+    Py_DECREF(ret_dummy);
+    return (PyObject *)ret;
+  }
+  /* error */
+  Py_DECREF(ret);
+  return NULL;
+}
+
+/** \note #BaseMath_ReadCallback must be called beforehand. */
+static PyObject *Vector_to_tuple_ex(VectorObject *self, int ndigits)
+{
+  PyObject *ret;
+  int i;
+
+  ret = PyTuple_New(self->vec_num);
+
+  if (ndigits >= 0) {
+    for (i = 0; i < self->vec_num; i++) {
+      PyTuple_SET_ITEM(ret, i, PyFloat_FromDouble(double_round((double)self->vec[i], ndigits)));
+    }
+  }
+  else {
+    for (i = 0; i < self->vec_num; i++) {
+      PyTuple_SET_ITEM(ret, i, PyFloat_FromDouble(self->vec[i]));
+    }
+  }
+
+  return ret;
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: `__new__` / `mathutils.Vector()`
+ * \{ */
 
 /**
  * Supports 2D, 3D, and 4D vector objects both int and float values
@@ -82,20 +171,12 @@ static PyObject *Vector_new(PyTypeObject *type, PyObject *args, PyObject *kwds)
   return Vector_CreatePyObject_alloc(vec, vec_num, type);
 }
 
-static PyObject *vec__apply_to_copy(PyObject *(*vec_func)(VectorObject *), VectorObject *self)
-{
-  PyObject *ret = Vector_copy(self);
-  PyObject *ret_dummy = vec_func((VectorObject *)ret);
-  if (ret_dummy) {
-    Py_DECREF(ret_dummy);
-    return (PyObject *)ret;
-  }
-  /* error */
-  Py_DECREF(ret);
-  return NULL;
-}
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Class Methods
+ * \{ */
 
-/*-----------------------CLASS-METHODS----------------------------*/
 PyDoc_STRVAR(C_Vector_Fill_doc,
              ".. classmethod:: Fill(size, fill=0.0)\n"
              "\n"
@@ -311,7 +392,12 @@ static PyObject *C_Vector_Repeat(PyObject *cls, PyObject *args)
   return Vector_CreatePyObject_alloc(vec, vec_num, (PyTypeObject *)cls);
 }
 
-/*-----------------------------METHODS---------------------------- */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Zero
+ * \{ */
+
 PyDoc_STRVAR(Vector_zero_doc,
              ".. method:: zero()\n"
              "\n"
@@ -331,6 +417,12 @@ static PyObject *Vector_zero(VectorObject *self)
   Py_RETURN_NONE;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Normalize
+ * \{ */
+
 PyDoc_STRVAR(Vector_normalize_doc,
              ".. method:: normalize()\n"
              "\n"
@@ -364,6 +456,12 @@ static PyObject *Vector_normalized(VectorObject *self)
   return vec__apply_to_copy(Vector_normalize, self);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Resize
+ * \{ */
+
 PyDoc_STRVAR(Vector_resize_doc,
              ".. method:: resize(size=3)\n"
              "\n"
@@ -553,6 +651,13 @@ static PyObject *Vector_resize_4d(VectorObject *self)
   self->vec_num = 4;
   Py_RETURN_NONE;
 }
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: To N-dimensions
+ * \{ */
+
 PyDoc_STRVAR(Vector_to_2d_doc,
              ".. method:: to_2d()\n"
              "\n"
@@ -605,6 +710,12 @@ static PyObject *Vector_to_4d(VectorObject *self)
   return Vector_CreatePyObject(tvec, 4, Py_TYPE(self));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: To Tuple
+ * \{ */
+
 PyDoc_STRVAR(Vector_to_tuple_doc,
              ".. method:: to_tuple(precision=-1)\n"
              "\n"
@@ -614,28 +725,6 @@ PyDoc_STRVAR(Vector_to_tuple_doc,
              "   :type precision: int\n"
              "   :return: the values of the vector rounded by *precision*\n"
              "   :rtype: tuple\n");
-/* NOTE: BaseMath_ReadCallback must be called beforehand. */
-static PyObject *Vector_to_tuple_ex(VectorObject *self, int ndigits)
-{
-  PyObject *ret;
-  int i;
-
-  ret = PyTuple_New(self->vec_num);
-
-  if (ndigits >= 0) {
-    for (i = 0; i < self->vec_num; i++) {
-      PyTuple_SET_ITEM(ret, i, PyFloat_FromDouble(double_round((double)self->vec[i], ndigits)));
-    }
-  }
-  else {
-    for (i = 0; i < self->vec_num; i++) {
-      PyTuple_SET_ITEM(ret, i, PyFloat_FromDouble(self->vec[i]));
-    }
-  }
-
-  return ret;
-}
-
 static PyObject *Vector_to_tuple(VectorObject *self, PyObject *args)
 {
   int ndigits = 0;
@@ -662,6 +751,12 @@ static PyObject *Vector_to_tuple(VectorObject *self, PyObject *args)
   return Vector_to_tuple_ex(self, ndigits);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: To Track Quaternion
+ * \{ */
+
 PyDoc_STRVAR(Vector_to_track_quat_doc,
              ".. method:: to_track_quat(track, up)\n"
              "\n"
@@ -781,6 +876,12 @@ static PyObject *Vector_to_track_quat(VectorObject *self, PyObject *args)
   return Quaternion_CreatePyObject(quat, NULL);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Orthogonal
+ * \{ */
+
 PyDoc_STRVAR(
     Vector_orthogonal_doc,
     ".. method:: orthogonal()\n"
@@ -816,12 +917,15 @@ static PyObject *Vector_orthogonal(VectorObject *self)
   return Vector_CreatePyObject(vec, self->vec_num, Py_TYPE(self));
 }
 
-/**
- * Vector.reflect(mirror): return a reflected vector on the mirror normal.
- * <pre>
- * vec - ((2 * dot(vec, mirror)) * mirror)
- * </pre>
- */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Reflect
+ *
+ * `Vector.reflect(mirror)`: return a reflected vector on the mirror normal:
+ * `vec - ((2 * dot(vec, mirror)) * mirror)`.
+ * \{ */
+
 PyDoc_STRVAR(Vector_reflect_doc,
              ".. method:: reflect(mirror)\n"
              "\n"
@@ -866,6 +970,12 @@ static PyObject *Vector_reflect(VectorObject *self, PyObject *value)
   return Vector_CreatePyObject(reflect, self->vec_num, Py_TYPE(self));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Cross Product
+ * \{ */
+
 PyDoc_STRVAR(Vector_cross_doc,
              ".. method:: cross(other)\n"
              "\n"
@@ -908,6 +1018,12 @@ static PyObject *Vector_cross(VectorObject *self, PyObject *value)
   return ret;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Dot Product
+ * \{ */
+
 PyDoc_STRVAR(Vector_dot_doc,
              ".. method:: dot(other)\n"
              "\n"
@@ -936,6 +1052,12 @@ static PyObject *Vector_dot(VectorObject *self, PyObject *value)
   return ret;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Angle
+ * \{ */
+
 PyDoc_STRVAR(
     Vector_angle_doc,
     ".. function:: angle(other, fallback=None)\n"
@@ -1001,6 +1123,12 @@ static PyObject *Vector_angle(VectorObject *self, PyObject *args)
   return PyFloat_FromDouble(saacos(dot / (sqrt(dot_self) * sqrt(dot_other))));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Angle Signed
+ * \{ */
+
 PyDoc_STRVAR(
     Vector_angle_signed_doc,
     ".. function:: angle_signed(other, fallback)\n"
@@ -1055,6 +1183,12 @@ static PyObject *Vector_angle_signed(VectorObject *self, PyObject *args)
   return PyFloat_FromDouble(angle_signed_v2v2(self->vec, tvec));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Rotation Difference
+ * \{ */
+
 PyDoc_STRVAR(Vector_rotation_difference_doc,
              ".. function:: rotation_difference(other)\n"
              "\n"
@@ -1096,6 +1230,12 @@ static PyObject *Vector_rotation_difference(VectorObject *self, PyObject *value)
   return Quaternion_CreatePyObject(quat, NULL);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Project
+ * \{ */
+
 PyDoc_STRVAR(Vector_project_doc,
              ".. function:: project(other)\n"
              "\n"
@@ -1134,6 +1274,12 @@ static PyObject *Vector_project(VectorObject *self, PyObject *value)
   return Vector_CreatePyObject_alloc(tvec, vec_num, Py_TYPE(self));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Linear Interpolation
+ * \{ */
+
 PyDoc_STRVAR(Vector_lerp_doc,
              ".. function:: lerp(other, factor)\n"
              "\n"
@@ -1170,6 +1316,12 @@ static PyObject *Vector_lerp(VectorObject *self, PyObject *args)
   return Vector_CreatePyObject_alloc(tvec, vec_num, Py_TYPE(self));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Spherical Interpolation
+ * \{ */
+
 PyDoc_STRVAR(Vector_slerp_doc,
              ".. function:: slerp(other, factor, fallback=None)\n"
              "\n"
@@ -1256,6 +1408,12 @@ static PyObject *Vector_slerp(VectorObject *self, PyObject *args)
   return Vector_CreatePyObject(ret_vec, vec_num, Py_TYPE(self));
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Rotate
+ * \{ */
+
 PyDoc_STRVAR(
     Vector_rotate_doc,
     ".. function:: rotate(other)\n"
@@ -1297,6 +1455,34 @@ static PyObject *Vector_rotate(VectorObject *self, PyObject *value)
   Py_RETURN_NONE;
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Negate
+ * \{ */
+
+PyDoc_STRVAR(Vector_negate_doc,
+             ".. method:: negate()\n"
+             "\n"
+             "   Set all values to their negative.\n");
+static PyObject *Vector_negate(VectorObject *self)
+{
+  if (BaseMath_ReadCallback(self) == -1) {
+    return NULL;
+  }
+
+  negate_vn(self->vec, self->vec_num);
+
+  (void)BaseMath_WriteCallback(self); /* already checked for error */
+  Py_RETURN_NONE;
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Methods: Copy/Deep-Copy
+ * \{ */
+
 PyDoc_STRVAR(Vector_copy_doc,
              ".. function:: copy()\n"
              "\n"
@@ -1323,6 +1509,12 @@ static PyObject *Vector_deepcopy(VectorObject *self, PyObject *args)
   return Vector_copy(self);
 }
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: `__repr__` & `__str__`
+ * \{ */
+
 static PyObject *Vector_repr(VectorObject *self)
 {
   PyObject *ret, *tuple;
@@ -1362,13 +1554,124 @@ static PyObject *Vector_str(VectorObject *self)
 }
 #endif
 
-/* Sequence Protocol */
-/* sequence length len(vector) */
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: Rich Compare
+ * \{ */
+
+static PyObject *Vector_richcmpr(PyObject *objectA, PyObject *objectB, int comparison_type)
+{
+  VectorObject *vecA = NULL, *vecB = NULL;
+  int result = 0;
+  const double epsilon = 0.000001f;
+  double lenA, lenB;
+
+  if (!VectorObject_Check(objectA) || !VectorObject_Check(objectB)) {
+    if (comparison_type == Py_NE) {
+      Py_RETURN_TRUE;
+    }
+
+    Py_RETURN_FALSE;
+  }
+  vecA = (VectorObject *)objectA;
+  vecB = (VectorObject *)objectB;
+
+  if (BaseMath_ReadCallback(vecA) == -1 || BaseMath_ReadCallback(vecB) == -1) {
+    return NULL;
+  }
+
+  if (vecA->vec_num != vecB->vec_num) {
+    if (comparison_type == Py_NE) {
+      Py_RETURN_TRUE;
+    }
+
+    Py_RETURN_FALSE;
+  }
+
+  switch (comparison_type) {
+    case Py_LT:
+      lenA = len_squared_vn(vecA->vec, vecA->vec_num);
+      lenB = len_squared_vn(vecB->vec, vecB->vec_num);
+      if (lenA < lenB) {
+        result = 1;
+      }
+      break;
+    case Py_LE:
+      lenA = len_squared_vn(vecA->vec, vecA->vec_num);
+      lenB = len_squared_vn(vecB->vec, vecB->vec_num);
+      if (lenA < lenB) {
+        result = 1;
+      }
+      else {
+        result = (((lenA + epsilon) > lenB) && ((lenA - epsilon) < lenB));
+      }
+      break;
+    case Py_EQ:
+      result = EXPP_VectorsAreEqual(vecA->vec, vecB->vec, vecA->vec_num, 1);
+      break;
+    case Py_NE:
+      result = !EXPP_VectorsAreEqual(vecA->vec, vecB->vec, vecA->vec_num, 1);
+      break;
+    case Py_GT:
+      lenA = len_squared_vn(vecA->vec, vecA->vec_num);
+      lenB = len_squared_vn(vecB->vec, vecB->vec_num);
+      if (lenA > lenB) {
+        result = 1;
+      }
+      break;
+    case Py_GE:
+      lenA = len_squared_vn(vecA->vec, vecA->vec_num);
+      lenB = len_squared_vn(vecB->vec, vecB->vec_num);
+      if (lenA > lenB) {
+        result = 1;
+      }
+      else {
+        result = (((lenA + epsilon) > lenB) && ((lenA - epsilon) < lenB));
+      }
+      break;
+    default:
+      printf("The result of the comparison could not be evaluated");
+      break;
+  }
+  if (result == 1) {
+    Py_RETURN_TRUE;
+  }
+
+  Py_RETURN_FALSE;
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: Hash (`__hash__`)
+ * \{ */
+
+static Py_hash_t Vector_hash(VectorObject *self)
+{
+  if (BaseMath_ReadCallback(self) == -1) {
+    return -1;
+  }
+
+  if (BaseMathObject_Prepare_ForHash(self) == -1) {
+    return -1;
+  }
+
+  return mathutils_array_hash(self->vec, self->vec_num);
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: Sequence & Mapping Protocols Implementation
+ * \{ */
+
+/** Sequence length: `len(object)`. */
 static int Vector_len(VectorObject *self)
 {
   return self->vec_num;
 }
-/* sequence accessor (get): vector[index] */
+
 static PyObject *vector_item_internal(VectorObject *self, int i, const bool is_attr)
 {
   if (i < 0) {
@@ -1395,11 +1698,12 @@ static PyObject *vector_item_internal(VectorObject *self, int i, const bool is_a
   return PyFloat_FromDouble(self->vec[i]);
 }
 
+/** Sequence accessor (get): `x = object[i]`. */
 static PyObject *Vector_item(VectorObject *self, int i)
 {
   return vector_item_internal(self, i, false);
 }
-/* sequence accessor (set): vector[index] = value */
+
 static int vector_ass_item_internal(VectorObject *self, int i, PyObject *value, const bool is_attr)
 {
   float scalar;
@@ -1442,12 +1746,13 @@ static int vector_ass_item_internal(VectorObject *self, int i, PyObject *value,
   return 0;
 }
 
+/** Sequence accessor (set): `object[i] = x`. */
 static int Vector_ass_item(VectorObject *self, int i, PyObject *value)
 {
   return vector_ass_item_internal(self, i, value, false);
 }
 
-/* sequence slice (get): vector[a:b] */
+/** Sequence slice accessor (get): `x = object[i:j]`. */
 static PyObject *Vector_slice(VectorObject *self, int begin, int end)
 {
   PyObject *tuple;
@@ -1471,7 +1776,8 @@ static PyObject *Vector_slice(VectorObject *self, int begin, int end)
 
   return tuple;
 }
-/* sequence slice (set): vector[a:b] = value */
+
+/** Sequence slice accessor (set): `object[i:j] = x`. */
 static int Vector_ass_slice(VectorObject *self, int begin, int end, PyObject *seq)
 {
   int vec_num = 0;
@@ -1509,8 +1815,83 @@ static int Vector_ass_slice(VectorObject *self, int begin, int end, PyObject *se
   return 0;
 }
 
-/* Numeric Protocols */
-/* addition: obj + obj */
+/** Sequence generic subscript (get): `x = object[...]`. */
+static PyObject *Vector_subscript(VectorObject *self, PyObject *item)
+{
+  if (PyIndex_Check(item)) {
+    Py_ssize_t i;
+    i = PyNumber_AsSsize_t(item, PyExc_IndexError);
+    if (i == -1 && PyErr_Occurred()) {
+      return NULL;
+    }
+    if (i < 0) {
+      i += self->vec_num;
+    }
+    return Vector_item(self, i);
+  }
+  if (PySlice_Check(item)) {
+    Py_ssize_t start, stop, step, slicelength;
+
+    if (PySlice_GetIndicesEx(item, self->vec_num, &start, &stop, &step, &slicelength) < 0) {
+      return NULL;
+    }
+
+    if (slicelength <= 0) {
+      return PyTuple_New(0);
+    }
+    if (step == 1) {
+      return Vector_slice(self, start, stop);
+    }
+
+    PyErr_SetString(PyExc_IndexError, "slice steps not supported with vectors");
+    return NULL;
+  }
+
+  PyErr_Format(
+      PyExc_TypeError, "vector indices must be integers, not %.200s", Py_TYPE(item)->tp_name);
+  return NULL;
+}
+
+/** Sequence generic subscript (set): `object[...] = x`. */
+static int Vector_ass_subscript(VectorObject *self, PyObject *item, PyObject *value)
+{
+  if (PyIndex_Check(item)) {
+    Py_ssize_t i = PyNumber_AsSsize_t(item, PyExc_IndexError);
+    if (i == -1 && PyErr_Occurred()) {
+      return -1;
+    }
+    if (i < 0) {
+      i += self->vec_num;
+    }
+    return Vector_ass_item(self, i, value);
+  }
+  if (PySlice_Check(item)) {
+    Py_ssize_t start, stop, step, slicelength;
+
+    if (PySlice_GetIndicesEx(item, self->vec_num, &start, &stop, &step, &slicelength) < 0) {
+      return -1;
+    }
+
+    if (step == 1) {
+      return Vector_ass_slice(self, start, stop, value);
+    }
+
+    PyErr_SetString(PyExc_IndexError, "slice steps not supported with vectors");
+    return -1;
+  }
+
+  PyErr_Format(
+      PyExc_TypeError, "vector indices must be integers, not %.200s", Py_TYPE(item)->tp_name);
+  return -1;
+}
+
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: Numeric Protocol Implementation
+ * \{ */
+
+/** Addition: `object + object`. */
 static PyObject *Vector_add(PyObject *v1, PyObject *v2)
 {
   VectorObject *vec1 = NULL, *vec2 = NULL;
@@ -1552,7 +1933,7 @@ static PyObject *Vector_add(PyObject *v1, PyObject *v2)
   return Vector_CreatePyObject_alloc(vec, vec1->vec_num, Py_TYPE(v1));
 }
 
-/* addition in-place: obj += obj */
+/** Addition in-place: `object += object`. */
 static PyObject *Vector_iadd(PyObject *v1, PyObject *v2)
 {
   VectorObject *vec1 = NULL, *vec2 = NULL;
@@ -1586,7 +1967,7 @@ static PyObject *Vector_iadd(PyObject *v1, PyObject *v2)
   return v1;
 }
 
-/* subtraction: obj - obj */
+/** Subtraction: `object - object`. */
 static PyObject *Vector_sub(PyObject *v1, PyObject *v2)
 {
   VectorObject *vec1 = NULL, *vec2 = NULL;
@@ -1627,7 +2008,7 @@ static PyObject *Vector_sub(PyObject *v1, PyObject *v2)
   return Vector_CreatePyObject_alloc(vec, vec1->vec_num, Py_TYPE(v1));
 }
 
-/* subtraction in-place: obj -= obj */
+/** Subtraction in-place: `object -= object`. */
 static PyObject *Vector_isub(PyObject *v1, PyObject *v2)
 {
   VectorObject *vec1 = NULL, *vec2 = NULL;
@@ -1661,8 +2042,7 @@ static PyObject *Vector_isub(PyObject *v1, PyObject *v2)
   return v1;
 }
 
-/*------------------------obj * obj------------------------------
- * multiplication */
+/* Multiply internal implementation `object * object`, `object *= object`. */
 
 int column_vector_multiplication(float r_vec[MAX_DIMENSIONS], VectorObject *vec, MatrixObject *mat)
 {
@@ -1725,6 +2105,7 @@ static PyObject *vector_mul_vec(VectorObject *vec1, VectorObject *vec2)
   return Vector_CreatePyObject_alloc(tvec, vec1->vec_num, Py_TYPE(vec1));
 }
 
+/** Multiplication (element-wise or scalar): `object * object`. */
 static PyObject *Vector_mul(PyObject *v1, PyObject *v2)
 {
   VectorObject *vec1 = NULL, *vec2 = NULL;
@@ -1776,7 +2157,7 @@ static PyObject *Vector_mul(PyObject *v1, PyObject *v2)
   return NULL;
 }
 
-/* multiplication in-place: obj *= obj */
+/** Multiplication in-place (element-wise or scalar): `object *= object`. */
 static PyObject *Vector_imul(PyObject *v1, PyObject *v2)
 {
   VectorObject *vec1 = NULL, *vec2 = NULL;
@@ -1830,6 +2211,7 @@ static PyObject *Vector_imul(PyObject *v1, PyObject *v2)
   return v1;
 }
 
+/** Multiplication (matrix multiply): `object @ object`. */
 static PyObject *Vector_matmul(PyObject *v1, PyObject *v2)
 {
   VectorObject *vec1 = NULL, *vec2 = NULL;
@@ -1893,6 +2275,7 @@ static PyObject *Vector_matmul(PyObject *v1, PyObject *v2)
   return NULL;
 }
 
+/** Multiplication in-place (matrix multiply): `object @= object`. */
 static PyObject *Vector_imatmul(PyObject *v1, PyObject *v2)
 {
   PyErr_Format(PyExc_TypeError,
@@ -1903,7 +2286,7 @@ static PyObject *Vector_imatmul(PyObject *v1, PyObject *v2)
   return NULL;
 }
 
-/* divide: obj / obj */
+/** Division: `object / object`. */
 static PyObject *Vector_div(PyObject *v1, PyObject *v2)
 {
   float *vec = NULL, scalar;
@@ -1950,7 +2333,7 @@ static PyObject *Vector_div(PyObject *v1, PyObject *v2)
   return Vector_CreatePyObject_alloc(vec, vec1->vec_num, Py_TYPE(v1));
 }
 
-/* divide in-place: obj /= obj */
+/** Division in-place: `object /= object`. */
 static PyObject *Vector_idiv(PyObject *v1, PyObject *v2)
 {
   float scalar;
@@ -1983,8 +2366,7 @@ static PyObject *Vector_idiv(PyObject *v1, PyObject *v2)
   return v1;
 }
 
-/* -obj
- * Returns the negative of this object. */
+/** Negative (returns the negative of this object): `-object`. */
 static PyObject *Vector_neg(VectorObject *self)
 {
   float *tvec;
@@ -1998,184 +2380,25 @@ static PyObject *Vector_neg(VectorObject *self)
   return Vector_CreatePyObject_alloc(tvec, self->vec_num, Py_TYPE(self));
 }
 
-/*------------------------tp_richcmpr
- * returns -1 exception, 0 false, 1 true */
-static PyObject *Vector_richcmpr(PyObject *objectA, PyObject *objectB, int comparison_type)
-{
-  VectorObject *vecA = NULL, *vecB = NULL;
-  int result = 0;
-  const double epsilon = 0.000001f;
-  double lenA, lenB;
-
-  if (!VectorObject_Check(objectA) || !VectorObject_Check(objectB)) {
-    if (comparison_type == Py_NE) {
-      Py_RETURN_TRUE;
-    }
-
-    Py_RETURN_FALSE;
-  }
-  vecA = (VectorObject *)objectA;
-  vecB = (VectorObject *)objectB;
-
-  if (BaseMath_ReadCallback(vecA) == -1 || BaseMath_ReadCallback(vecB) == -1) {
-    return NULL;
-  }
+/** \} */
 
-  if (vecA->vec_num != vecB->vec_num) {
-    if (comparison_type == Py_NE) {
-      Py_RETURN_TRUE;
-    }
-
-    Py_RETURN_FALSE;
-  }
-
-  switch (comparison_type) {
-    case Py_LT:
-      lenA = len_squared_vn(vecA->vec, vecA->vec_num);
-      lenB = len_squared_vn(vecB->vec, vecB->vec_num);
-      if (lenA < lenB) {
-        result = 1;
-      }
-      break;
-    case Py_LE:
-      lenA = len_squared_vn(vecA->vec, vecA->vec_num);
-      lenB = len_squared_vn(vecB->vec, vecB->vec_num);
-      if (lenA < lenB) {
-        result = 1;
-      }
-      else {
-        result = (((lenA + epsilon) > lenB) && ((lenA - epsilon) < lenB));
-      }
-      break;
-    case Py_EQ:
-      result = EXPP_VectorsAreEqual(vecA->vec, vecB->vec, vecA->vec_num, 1);
-      break;
-    case Py_NE:
-      result = !EXPP_VectorsAreEqual(vecA->vec, vecB->vec, vecA->vec_num, 1);
-      break;
-    case Py_GT:
-      lenA = len_squared_vn(vecA->vec, vecA->vec_num);
-      lenB = len_squared_vn(vecB->vec, vecB->vec_num);
-      if (lenA > lenB) {
-        result = 1;
-      }
-      break;
-    case Py_GE:
-      lenA = len_squared_vn(vecA->vec, vecA->vec_num);
-      lenB = len_squared_vn(vecB->vec, vecB->vec_num);
-      if (lenA > lenB) {
-        result = 1;
-      }
-      else {
-        result = (((lenA + epsilon) > lenB) && ((lenA - epsilon) < lenB));
-      }
-      break;
-    default:
-      printf("The result of the comparison could not be evaluated");
-      break;
-  }
-  if (result == 1) {
-    Py_RETURN_TRUE;
-  }
-
-  Py_RETURN_FALSE;
-}
-
-static Py_hash_t Vector_hash(VectorObject *self)
-{
-  if (BaseMath_ReadCallback(self) == -1) {
-    return -1;
-  }
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: Protocol Declarations
+ * \{ */
 
-  if (BaseMathObject_Prepare_ForHash(self) == -1) {
-    return -1;
-  }
-
-  return mathutils_array_hash(self->vec, self->vec_num);
-}
-
-/*-----------------PROTCOL DECLARATIONS--------------------------*/
 static PySequenceMethods Vector_SeqMethods = {
-    (lenfunc)Vector_len,              /* sq_length */
-    (binaryfunc)NULL,                 /* sq_concat */
-    (ssizeargfunc)NULL,               /* sq_repeat */
-    (ssizeargfunc)Vector_item,        /* sq_item */
-    NULL,                             /* py3 deprecated slice func */
-    (ssizeobjargproc)Vector_ass_item, /* sq_ass_item */
-    NULL,                             /* py3 deprecated slice assign func */
-    (objobjproc)NULL,                 /* sq_contains */
-    (binaryfunc)NULL,                 /* sq_inplace_concat */
-    (ssizeargfunc)NULL,               /* sq_inplace_repeat */
+    (lenfunc)Vector_len,              /*sq_length*/
+    (binaryfunc)NULL,                 /*sq_concat*/
+    (ssizeargfunc)NULL,               /*sq_repeat*/
+    (ssizeargfunc)Vector_item,        /*sq_item*/
+    NULL,                             /*sq_slice(DEPRECATED)*/
+    (ssizeobjargproc)Vector_ass_item, /*sq_ass_item*/
+    NULL,                             /*sq_ass_slice(DEPRECATED)*/
+    (objobjproc)NULL,                 /*sq_contains*/
+    (binaryfunc)NULL,                 /*sq_inplace_concat */
+    (ssizeargfunc)NULL,               /*sq_inplace_repeat */
 };
 
-static PyObject *Vector_subscript(VectorObject *self, PyObject *item)
-{
-  if (PyIndex_Check(item)) {
-    Py_ssize_t i;
-    i = PyNumber_AsSsize_t(item, PyExc_IndexError);
-    if (i == -1 && PyErr_Occurred()) {
-      return NULL;
-    }
-    if (i < 0) {
-      i += self->vec_num;
-    }
-    return Vector_item(self, i);
-  }
-  if (PySlice_Check(item)) {
-    Py_ssize_t start, stop, step, slicelength;
-
-    if (PySlice_GetIndicesEx(item, self->vec_num, &start, &stop, &step, &slicelength) < 0) {
-      return NULL;
-    }
-
-    if (slicelength <= 0) {
-      return PyTuple_New(0);
-    }
-    if (step == 1) {
-      return Vector_slice(self, start, stop);
-    }
-
-    PyErr_SetString(PyExc_IndexError, "slice steps not supported with vectors");
-    return NULL;
-  }
-
-  PyErr_Format(
-      PyExc_TypeError, "vector indices must be integers, not %.200s", Py_TYPE(item)->tp_name);
-  return NULL;
-}
-
-static int Vector_ass_subscript(VectorObject *self, PyObject *item, PyObject *value)
-{
-  if (PyIndex_Check(item)) {
-    Py_ssize_t i = PyNumber_AsSsize_t(item, PyExc_IndexError);
-    if (i == -1 && PyErr_Occurred()) {
-      return -1;
-    }
-    if (i < 0) {
-      i += self->vec_num;
-    }
-    return Vector_ass_item(self, i, value);
-  }
-  if (PySlice_Check(item)) {
-    Py_ssize_t start, stop, step, slicelength;
-
-    if (PySlice_GetIndicesEx(item, self->vec_num, &start, &stop, &step, &slicelength) < 0) {
-      return -1;
-    }
-
-    if (step == 1) {
-      return Vector_ass_slice(self, start, stop, value);
-    }
-
-    PyErr_SetString(PyExc_IndexError, "slice steps not supported with vectors");
-    return -1;
-  }
-
-  PyErr_Format(
-      PyExc_TypeError, "vector indices must be integers, not %.200s", Py_TYPE(item)->tp_name);
-  return -1;
-}
-
 static PyMappingMethods Vector_AsMapping = {
     (lenfunc)Vector_len,
     (binaryfunc)Vector_subscript,
@@ -2202,28 +2425,32 @@ static PyNumberMethods Vector_NumMethods = {
     NULL,                       /*nb_int*/
     NULL,                       /*nb_reserved*/
     NULL,                       /*nb_float*/
-    Vector_iadd,                /* nb_inplace_add */
-    Vector_isub,                /* nb_inplace_subtract */
-    Vector_imul,                /* nb_inplace_multiply */
-    NULL,                       /* nb_inplace_remainder */
-    NULL,                       /* nb_inplace_power */
-    NULL,                       /* nb_inplace_lshift */
-    NULL,                       /* nb_inplace_rshift */
-    NULL,                       /* nb_inplace_and */
-    NULL,                       /* nb_inplace_xor */
-    NULL,                       /* nb_inplace_or */
-    NULL,                       /* nb_floor_divide */
-    Vector_div,                 /* nb_true_divide */
-    NULL,                       /* nb_inplace_floor_divide */
-    Vector_idiv,                /* nb_inplace_true_divide */
-    NULL,                       /* nb_index */
-    (binaryfunc)Vector_matmul,  /* nb_matrix_multiply */
-    (binaryfunc)Vector_imatmul, /* nb_inplace_matrix_multiply */
+    Vector_iadd,                /*nb_inplace_add*/
+    Vector_isub,                /*nb_inplace_subtract*/
+    Vector_imul,                /*nb_inplace_multiply*/
+    NULL,                       /*nb_inplace_remainder*/
+    NULL,                       /*nb_inplace_power*/
+    NULL,                       /*nb_inplace_lshift*/
+    NULL,                       /*nb_inplace_rshift*/
+    NULL,                       /*nb_inplace_and*/
+    NULL,                       /*nb_inplace_xor*/
+    NULL,                       /*nb_inplace_or*/
+    NULL,                       /*nb_floor_divide*/
+    Vector_div,                 /*nb_true_divide*/
+    NULL,                       /*nb_inplace_floor_divide*/
+    Vector_idiv,                /*nb_inplace_true_divide*/
+    NULL,                       /*nb_index*/
+    (binaryfunc)Vector_matmul,  /*nb_matrix_multiply*/
+    (binaryfunc)Vector_imatmul, /*nb_inplace_matrix_multiply*/
 };
 
-/*------------------PY_OBECT DEFINITION--------------------------*/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: Get/Set Item Implementation
+ * \{ */
 
-/* vector axis, vector.x/y/z/w */
+/* Vector axis: `vector.x/y/z/w`. */
 
 PyDoc_STRVAR(Vector_axis_x_doc, "Vector X axis.\n\n:type: float");
 PyDoc_STRVAR(Vector_axis_y_doc, "Vector Y axis.\n\n:type: float");
@@ -2240,7 +2467,7 @@ static int Vector_axis_set(VectorObject *self, PyObject *value, void *type)
   return vector_ass_item_internal(self, POINTER_AS_INT(type), value, true);
 }
 
-/* vector.length */
+/* `Vector.length`. */
 
 PyDoc_STRVAR(Vector_length_doc, "Vector Length.\n\n:type: float");
 static PyObject *Vector_length_get(VectorObject *self, void *UNUSED(closure))
@@ -2296,7 +2523,7 @@ static int Vector_length_set(VectorObject *self, PyObject *value)
   return 0;
 }
 
-/* vector.length_squared */
+/* `Vector.length_squared`. */
 PyDoc_STRVAR(Vector_length_squared_doc, "Vector length squared (v.dot(v)).\n\n:type: float");
 static PyObject *Vector_length_squared_get(VectorObject *self, void *UNUSED(closure))
 {
@@ -2503,9 +2730,12 @@ static int Vector_swizzle_set(VectorObject *self, PyObject *value, void *closure
 #define SWIZZLE3(a, b, c) POINTER_FROM_INT(_SWIZZLE3(a, b, c))
 #define SWIZZLE4(a, b, c, d) POINTER_FROM_INT(_SWIZZLE4(a, b, c, d))
 
-/*****************************************************************************/
-/* Python attributes get/set structure:                                      */
-/*****************************************************************************/
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: Get/Set Item Definitions
+ * \{ */
+
 static PyGetSetDef Vector_getseters[] = {
     {"x", (getter)Vector_axis_get, (setter)Vector_axis_set, Vector_axis_x_doc, (void *)0},
     {"y", (getter)Vector_axis_get, (setter)Vector_axis_set, Vector_axis_y_doc, (void *)1},
@@ -2886,68 +3116,11 @@ static PyGetSetDef Vector_getseters[] = {
     {NULL, NULL, NULL, NULL, NULL} /* Sentinel */
 };
 
-/**
- * Row vector multiplication - (Vector * Matrix)
- * <pre>
- * [x][y][z] * [1][4][7]
- *             [2][5][8]
- *             [3][6][9]
- * </pre>
- * \note vector/matrix multiplication is not commutative.
- */
-static int row_vector_multiplication(float r_vec[MAX_DIMENSIONS],
-                                     VectorObject *vec,
-                                     MatrixObject *mat)
-{
-  float vec_cpy[MAX_DIMENSIONS];
-  int row, col, z = 0, vec_num = vec->vec_num;
-
-  if (mat->row_num != vec_num) {
-    if (mat->row_num == 4 && vec_num == 3) {
-      vec_cpy[3] = 1.0f;
-    }
-    else {
-      PyErr_SetString(PyExc_ValueError,
-                      "vector * matrix: matrix column size "
-                      "and the vector size must be the same");
-      return -1;
-    }
-  }
-
-  if (BaseMath_ReadCallback(vec) == -1 || BaseMath_ReadCallback(mat) == -1) {
-    return -1;
-  }
-
-  memcpy(vec_cpy, vec->vec, vec_num * sizeof(float));
-
-  r_vec[3] = 1.0f;
-  /* Multiplication. */
-  for (col = 0; col < mat->col_num; col++) {
-    double dot = 0.0;
-    for (row = 0; row < mat->row_num; row++) {
-      dot += (double)(MATRIX_ITEM(mat, row, col) * vec_cpy[row]);
-    }
-    r_vec[z++] = (float)dot;
-  }
-  return 0;
-}
-
-/*----------------------------Vector.negate() -------------------- */
-PyDoc_STRVAR(Vector_negate_doc,
-             ".. method:: negate()\n"
-             "\n"
-             "   Set all values to their negative.\n");
-static PyObject *Vector_negate(VectorObject *self)
-{
-  if (BaseMath_ReadCallback(self) == -1) {
-    return NULL;
-  }
-
-  negate_vn(self->vec, self->vec_num);
+/** \} */
 
-  (void)BaseMath_WriteCallback(self); /* already checked for error */
-  Py_RETURN_NONE;
-}
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: Method Definitions
+ * \{ */
 
 static struct PyMethodDef Vector_methods[] = {
     /* Class Methods */
@@ -3000,11 +3173,15 @@ static struct PyMethodDef Vector_methods[] = {
     {NULL, NULL, 0, NULL},
 };
 
-/**
- * NOTE: #Py_TPFLAGS_CHECKTYPES allows us to avoid casting all types to Vector when coercing
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: Python Object Definition
+ *
+ * \note #Py_TPFLAGS_CHECKTYPES allows us to avoid casting all types to Vector when coercing
  * but this means for eg that (vec * mat) and (mat * vec)
  * both get sent to Vector_mul and it needs to sort out the order
- */
+ * \{ */
 
 PyDoc_STRVAR(vector_doc,
              ".. class:: Vector(seq)\n"
@@ -3098,6 +3275,12 @@ PyTypeObject vector_Type = {
     NULL,
 };
 
+/** \} */
+
+/* -------------------------------------------------------------------- */
+/** \name Vector Type: C/API Constructors
+ * \{ */
+
 PyObject *Vector_CreatePyObject(const float *vec, const int vec_num, PyTypeObject *base_type)
 {
   VectorObject *self;
@@ -3190,3 +3373,5 @@ PyObject *Vector_CreatePyObject_alloc(float *vec, const int vec_num, PyTypeObjec
 
   return (PyObject *)self;
 }
+
+/** \} */
diff --git a/source/blender/python/mathutils/mathutils_Vector.h b/source/blender/python/mathutils/mathutils_Vector.h
index 3bc4e9d6b6f..7006ece6bf0 100644
--- a/source/blender/python/mathutils/mathutils_Vector.h
+++ b/source/blender/python/mathutils/mathutils_Vector.h
@@ -18,7 +18,8 @@ typedef struct {
   int vec_num;
 } VectorObject;
 
-/*prototypes*/
+/* Prototypes. */
+
 PyObject *Vector_CreatePyObject(const float *vec,
                                 int vec_num,
                                 PyTypeObject *base_type) ATTR_WARN_UNUSED_RESULT;