{ "cells": [ { "cell_type": "code", "execution_count": 1, "id": "78b84bf1-1786-45bd-bd5d-c173f27e59ab", "metadata": { "tags": [ "hide-input" ] }, "outputs": [], "source": [ "# Copyright (C) 2020-2021 Mael Kerbiriou \n", "#\n", "# This file may be distributed under the terms of the GNU GPLv3 license.\n", "from contextlib import contextmanager\n", "import numpy as np\n", "from matplotlib import pyplot as plt\n", "from matplotlib_inline.backend_inline import set_matplotlib_formats\n", "%matplotlib inline\n", "set_matplotlib_formats('svg')\n", "plt.style.use(['dark_background'])\n", "\n", "plt_kwargs = dict(\n", " toolhead=dict(color=(.5,.5,1), zorder=20),\n", " X=dict(color='red'),\n", " Y=dict(color=(.0,.8,.0)),\n", " A=dict(color='cyan'),\n", " B=dict(color='yellow'),\n", " )\n", "for k in plt_kwargs: plt_kwargs[k]['linewidth'] = 1.25\n", "plt_kwargs['toolhead']['linewidth'] = 3\n", "\n", "@contextmanager\n", "def fig_tmpl(title, ordinate):\n", " f = plt.figure(figsize=(10,5))\n", " plt.grid(False)\n", " plt.title(title)\n", " plt.xlabel('Move orientation (°)')\n", " plt.ylabel('Max %s (relative)' % ordinate)\n", " plt.xticks(np.linspace(0,180,4*3, endpoint=False))\n", " try:\n", " yield f\n", " finally:\n", " _, ymax = plt.ylim()\n", " plt.ylim((0,ymax))\n", " plt.xlim((0,180))\n", "\n", " plt.legend(loc='lower right')\n", " plt.tight_layout()\n", " \n", " \n", "@contextmanager\n", "def fig_tmpl(title, ordinate):\n", " f = plt.figure(figsize=(7,7))\n", " ax = f.add_subplot(projection='polar')\n", " \n", " #ax.set_thetamin(0)\n", " #ax.set_thetamax(180)\n", "\n", " plt.grid(c='gray', linewidth=0.5, linestyle='--', alpha=0.5)\n", " plt.title(title)\n", " plt.xlabel('Move orientation (°)')\n", " #plt.ylabel('Max %s (relative)' % ordinate)\n", " #plt.xticks(np.linspace(0,180,4*3, endpoint=False))\n", " try:\n", " yield f\n", " finally:\n", " _, ymax = plt.ylim()\n", " plt.ylim((0,ymax))\n", " #plt.xlim((0,180))\n", "\n", " plt.legend(loc='lower right')\n", " plt.tight_layout()\n", " \n", "def gen_moves(theta):\n", " d = np.random.rand(*theta.shape)*10+1 # Random move lengths\n", " return d*np.cos(theta), d*np.sin(theta)\n", "\n", "theta = np.linspace(0,2*np.pi, 127)\n", "theta_deg = (180/np.pi) * theta\n", "x,y = gen_moves(theta)\n", "d = np.hypot(x,y)" ] }, { "cell_type": "markdown", "id": "9f5a44c2", "metadata": {}, "source": [ "## Cartesian\n", "\n", "### Velocitiy limits\n", "\n", "When the toolhead on a cartesian printer moves at a constant speed, each independant axis move at a lower speed depending on the orientation. The following plot in polar coordiantes shows axes velocities when the toolhead moves at a constant velocity in every possible orientation:" ] }, { "cell_type": "code", "execution_count": 2, "id": "a71926c9", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n" ], "text/plain": [ "

" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "v_x, v_y = x / d, y / d\n", "\n", "with fig_tmpl('Cartesian axes velocities', 'velocity'):\n", " plt.plot(theta, x*0.0+1.0, label='toolhead', **plt_kwargs['toolhead'])\n", " plt.plot(theta, np.abs(v_x), label='$|\\dot X|$', **plt_kwargs['X'])\n", " plt.plot(theta, np.abs(v_y), label='$|\\dot Y|$', **plt_kwargs['Y'])" ] }, { "cell_type": "markdown", "id": "58e82efb", "metadata": {}, "source": [ "Red and green lines are the velocities of the axes at that constant velocity limit. All the plot will use the convention that a value of 1 represents the target limit.\n", "\n", "One of the goal of speed limits is to avoid reaching the step frequencies where the pull-out torque of steppers drops drastically. Notice that, at 45° and 135° degree, the axes moves at $1/\\sqrt 2 = 70.7\\%$ the speed of the toolhead. This is well bellow the speed limit of the motor, here represented by the value 1.\n", "\n", "This next polar plot shows what happens if we limit the velocities of the axes instead of the toolhead. The geometrical procedure for acheiving that requires clipping the X and Y velocity curves to a circle, and then scaling up to our limit $v_\\text{max}=1$:" ] }, { "cell_type": "code", "execution_count": 3, "id": "fd4d348e", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n" ], "text/plain": [ "

" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "max_vel = 1\n", "#max_vel = np.minimum(max_vel*d/np.abs(x), max_vel*d/np.abs(y))\n", "max_vel = max_vel * d / np.maximum(np.abs(x), np.abs(y))\n", "v_x, v_y = max_vel * x / d, max_vel * y / d\n", "\n", "with fig_tmpl('Cartesian per-axis velocity limiting', 'velocity'):\n", " plt.plot(theta, max_vel, label='toolhead', **plt_kwargs['toolhead'])\n", " plt.plot(theta, np.abs(v_x), label='$|\\dot X|$', **plt_kwargs['X'])\n", " plt.plot(theta, np.abs(v_y), label='$|\\dot Y|$', **plt_kwargs['Y'])" ] }, { "cell_type": "markdown", "id": "86855071-5782-4448-8551-d61c083bda28", "metadata": {}, "source": [ "When the toolhead is traveling at the speed limit (blue square), there is always one of the two steppers working at its speed limit, characterized by the circle of radius $v_\\text{max}=1$.\n", "\n", "At an angle of 45°, both motors run at $v_\\text{max}$, yielding a combined toolhead velocity of $\\sqrt{2}\\cdot v_\\text{max} \\simeq 141\\% \\cdot v_\\text{max}$." ] }, { "cell_type": "markdown", "id": "47ac8fb0", "metadata": {}, "source": [ "### Acceleration limits\n", "\n", "On most printers, like bed slingers and h-bots, one axis is heavier and require more torque to accelerate. To work closer to the machine limits, it is desirable to restrict the acceleration on the heaviest axis (often Y).\n", "\n", "Given a normed move orientation vector $ [ r_X \\; r_Y ] = [\\Delta X \\; \\Delta Y] / \\sqrt{\\Delta X^2 + \\Delta Y^2} $, we can rewrite the acceleration components as a function of the toolhead acceleration magnitude $a$:\n", "\n", "$$\n", "\\begin{cases}\n", "\\ddot X = a \\cdot r_X \\\\ \n", "\\ddot Y = a \\cdot r_Y \\\\ \n", "\\end{cases}\n", "$$\n", "\n", "Then, for each move, we must scale $a_\\text{max}$ to ensure that the acceleration of each steppers stays bounded, observing $ | \\ddot X | \\leq \\ddot X_\\text{max} $ and $| \\ddot Y | \\leq \\ddot Y _\\text{max} $:\n", "\n", "$$\n", "a_\\text{max} = \\min \\left( \\frac{\\ddot X_\\text{max}}{|r_X|}, \\, \\frac{\\ddot Y_\\text{max}}{|r_Y|} \\right)\n", "$$\n", "\n", "In polar coordinate, $a_\\text{max}$ will show up as a rectangle. For example with $X_\\text{max} = 1$ and $Y_\\text{max} = 3/4$:" ] }, { "cell_type": "code", "execution_count": 4, "id": "09bae3a5", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n" ], "text/plain": [ "

" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "max_acc_x = 1\n", "max_acc_y = 0.75\n", "max_acc = d/np.maximum(np.abs(x)/max_acc_x, np.abs(y)/max_acc_y)\n", "v_x, v_y = max_acc * x / d, max_acc * y / d\n", "\n", "with fig_tmpl('Cartesian per-axis acceleration limiting', 'acceleration'):\n", " plt.plot(theta, max_acc, label='toolhead', **plt_kwargs['toolhead'])\n", " plt.plot(theta, np.abs(v_x), label='$|\\ddot X|$', **plt_kwargs['X'])\n", " plt.plot(theta, np.abs(v_y), label='$|\\ddot Y|$', **plt_kwargs['Y'])" ] }, { "cell_type": "markdown", "id": "86c609f3", "metadata": {}, "source": [ "The angle with maximum acceleration is given by $\\arctan \\frac{\\ddot Y_\\text{max}}{\\ddot X_\\text{max}}$ ($\\simeq 37°$ here).\n", "\n", "Some printers have a torquier motor on Y, often acheiving a lower top speed due to higher back EMF and inductance. This could be a motivation to use different velocity limits on X and Y axes. The velocity profile will be a rectangle, exactly like the acceleration above.\n", "\n", "## CoreXY\n", "\n", "[CoreXY](http://corexy.com/theory.html) uses a different belt path allowing the motors to stay stationnary. The forward and reverse kinematics are:\n", "\n", "$$\n", "\\begin{cases}\n", "2X = A + B \\\\ \n", "2Y = A - B \\\\ \n", "\\end{cases}\n", "$$\n", "\n", "$$\n", "\\begin{cases}\n", "A = X + Y \\\\ \n", "B = X - Y \\\\ \n", "\\end{cases}\n", "$$ $A$ and $B$ designate the motor position. In order to avoid confusion with angular position A and B will be called \"belt\" positions.\n", "\n", "Since these relations are linear, the equations works equaly for position, velocity and acceleration.\n", "\n", "### Velocity limits\n", "\n", "A constant toolhead velocity $v$ requires belt velocites of up to $\\sqrt{2}\\cdot v$:" ] }, { "cell_type": "code", "execution_count": 10, "id": "c8626369", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n" ], "text/plain": [ "

" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "v_x, v_y = x / d, y / d\n", "v_a, v_b = v_x + v_y, v_x - v_y\n", "\n", "with fig_tmpl('CoreXY belts and axes velocities', 'velocity'):\n", " plt.plot(theta, np.hypot(v_x, v_y), label='toolhead', **plt_kwargs['toolhead'])\n", " plt.plot(theta, np.abs(v_a), label='$|\\dot A|$', **plt_kwargs['A'])\n", " plt.plot(theta, np.abs(v_b), label='$|\\dot B|$', **plt_kwargs['B'])\n", " plt.plot(theta, np.abs(v_x), label='$|\\dot X|$', **plt_kwargs['X'])\n", " plt.plot(theta, np.abs(v_y), label='$|\\dot Y|$', **plt_kwargs['Y'])" ] }, { "cell_type": "markdown", "id": "d1b2fcf3", "metadata": {}, "source": [ "X and Y velocities are identical to the cartesian plot. However, we now have to look at the belt velocities (or OD velocities of motor pulleys if you want), A and B. They form bigger pairs of circles, rotated by 45° compared to X and Y. This is caused by the factor 2 and the diagonal projection in the kinematics equation.\n", "\n", "For accelerations, Klipper's constant limiting will look exactly the same, thanks to the linearity of the kinematic equation.\n", "\n", "Limiting the belt velocities, or equivalently, clipping the A and B curves to a circle, decreases the max toolhead velocities on diagonals, with minimums at 45° and 135°:" ] }, { "cell_type": "code", "execution_count": 6, "id": "8e1d6700", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n" ], "text/plain": [ "

" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "max_vel = 1\n", "max_vel *= d/np.maximum(np.abs(x+y), np.abs(x-y))\n", "v_x, v_y = max_vel * x / d, max_vel * y / d\n", "v_a, v_b = v_x + v_y, v_x - v_y\n", "\n", "\n", "with fig_tmpl('CoreXY limiting belts velocities', 'velocity'):\n", " plt.plot(theta, max_vel, label='toolhead', **plt_kwargs['toolhead'])\n", " plt.plot(theta, np.abs(v_a), label='$|\\dot A|$', **plt_kwargs['A'])\n", " plt.plot(theta, np.abs(v_b), label='$|\\dot B|$', **plt_kwargs['B'])\n", " plt.plot(theta, np.abs(v_x), label='$|\\dot X|$', **plt_kwargs['X'])\n", " plt.plot(theta, np.abs(v_y), label='$|\\dot Y|$', **plt_kwargs['Y'])" ] }, { "cell_type": "markdown", "id": "e34093e0", "metadata": { "tags": [] }, "source": [ "Unlike cartesian, the motor velocities on CoreXY match the toolhead velocity when the max speed is acheived (at 0° and 90° orientation). However, in practive limiting motor velocity yields, again, a potential 41% increase of top speed, compared to just limiting the toolhead velocity.\n", "\n", "### Torque limiting\n", "\n", "On cartesian printers, acceleration is a direct proxy for torque since each moving mass accelerate on it own independant axis. This is no longer true for CoreXY where the effects of inertia felt by the motor working together or against each other depend on the orientation. The forces exerced on the toolhead are related to the accelerations and the masses involved:\n", "\n", "$$\n", "\\begin{cases}\n", "F_X = m_X \\ddot X \\\\ \n", "F_Y = m_Y \\ddot Y \\\\ \n", "\\end{cases}\n", "$$\n", "\n", "where\n", "\n", "$$\n", "\\begin{cases}\n", "m_X = m_\\text{toolhead} \\\\\n", "m_Y = m_\\text{toolhead} + m_\\text{gantry}\n", "\\end{cases}\n", "$$\n", "\n", "On the belts, at the motors, we get the forces:\n", "\n", "$$\n", "\\begin{cases}\n", "F_A = m_X \\ddot X + m_Y \\ddot Y \\\\ \n", "F_B = m_Y \\ddot X - m_Y \\ddot Y \\\\ \n", "\\end{cases}\n", "$$\n", "\n", "The following plot shows the forces for masses $m_\\text{toolhead}=1$, $m_\\text{gantry}=2/3$ and unitary acceleration:" ] }, { "cell_type": "code", "execution_count": 7, "id": "7c835384", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n" ], "text/plain": [ "

" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "m_toolhead = 1\n", "m_gantry = 2/3\n", "m_x = m_toolhead\n", "m_y = m_toolhead + m_gantry\n", "\n", "a_x = x / d\n", "a_y = y / d\n", "\n", "f_x = m_x * a_x\n", "f_y = m_y * a_y\n", "f = np.hypot(f_x, f_y)\n", "\n", "f_a = f_x + f_y\n", "f_b = f_x - f_y\n", "\n", "with fig_tmpl('CoreXY forces on belts and axes', 'force'):\n", " plt.plot(theta, f, label=r'$|F_\\mathrm{toolhead}|$', **plt_kwargs['toolhead'])\n", " plt.plot(theta, np.abs(f_a), label='$|F_A|$', **plt_kwargs['A'])\n", " plt.plot(theta, np.abs(f_b), label='$|F_B|$', **plt_kwargs['B'])\n", " plt.plot(theta, np.abs(f_x), label='$|F_X|$', **plt_kwargs['X'])\n", " plt.plot(theta, np.abs(f_y), label='$|F_Y|$', **plt_kwargs['Y'])" ] }, { "cell_type": "markdown", "id": "8f3a9b1a", "metadata": {}, "source": [ "If the displaced masses were equal in every direction ($m_X=m_Y$), we would again get the same profile as \"CoreXY belts and axes velocities\".\n", "Instead, the magnitude of the force on the toolhead (blue line) is squished horizontally: moving the toolhead plus gantry in the Y direction requires more force than just moving the toolhead in the X direction. As a result, the pairs of circles for the forces on A and B are rotated closer to the vertical axis.\n", "\n", "\n", "With the cartesian kinematics we substituted the axis accelerations by the magnitude of toolhead acceleration and normed move orentation components. Here, we do the same for our forces, yielding:\n", "\n", "$$\n", "\\begin{cases}\n", "F_A = a (m_X r_X + m_Y r_Y) \\\\ \n", "F_B = a (m_Y r_Y - m_Y r_Y) \\\\ \n", "\\end{cases}\n", "$$\n", "\n", "Then we derive the maximum acceleration such that we don't exceed the available torque, or equivalently, such that $F_A \\leq F_\\text{max}$ and $F_B \\leq F_\\text{max}$ (supposing that motors have equal torque):\n", "\n", "$$\n", "\\begin{align}\n", "a_\\text{max} &= \\min \\left( \\frac{F_\\text{max}}{|m_X r_X + m_Y r_Y|},\\, \\frac{F_\\text{max}}{|m_X r_X - m_Y r_Y|} \\right) \\\\\n", " &= \\frac{ F_\\text{max} }{ \\max \\left( |m_X r_X + m_Y r_Y|,\\, |m_X r_X - m_Y r_Y| \\right) } \\\\\n", "% &= \\frac{ F_\\text{max} \\sqrt{\\Delta X^2 + \\Delta Y^2} }\n", "% { \\max \\left( |m_X \\Delta X + m_Y \\Delta Y|,\\, |m_X \\Delta X - m_Y \\Delta Y| \\right) }\n", "\\end{align}\n", "$$\n", "\n", "Again, with Klipper's constant acceleration limiting, the acceleration plot shares the profile of \"CoreXY belts and axes velocities\". Instead, when limiting the forces at the belts to $F_\\text{max}=1$, we are now getting these accelerations:" ] }, { "cell_type": "code", "execution_count": 8, "id": "42824736-3893-42ed-9f61-16a93d93a285", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n" ], "text/plain": [ "

" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "f_max = 1\n", "m_x = 1\n", "m_y = m_x + 2/3\n", "\n", "r_x = x / d\n", "r_y = y / d\n", "\n", "#a = np.minimum(f_max / np.abs(m_x * r_x + m_y * r_y), f_max / np.abs(m_x * r_x - m_y * r_y))\n", "#a = f_max / np.maximum(np.abs(m_x * r_x + m_y * r_y), np.abs(m_x * r_x - m_y * r_y))\n", "#a = f_max * d / np.maximum(np.abs(m_x * x + m_y * y), np.abs(m_x * x - m_y * y))\n", "amax_x = f_max / m_x\n", "amax_y = f_max / m_y\n", "a = d / np.maximum(np.abs(x / amax_x + y / amax_y), np.abs(x / amax_x - y / amax_y))\n", "\n", "a_x = a * r_x\n", "a_y = a * r_y\n", "\n", "a_a = a_x + a_y\n", "a_b = a_x - a_y\n", "\n", "with fig_tmpl('CoreXY limited acceleration on belts and axes', 'acceleration') as f:\n", " plt.plot(theta, a, label='toolhead', **plt_kwargs['toolhead'])\n", " plt.plot(theta, np.abs(a_a), label='$|\\ddot A|$', **plt_kwargs['A'])\n", " plt.plot(theta, np.abs(a_b), label='$|\\ddot B|$', **plt_kwargs['B'])\n", " plt.plot(theta, np.abs(a_x), label='$|\\ddot X|$', **plt_kwargs['X'])\n", " plt.plot(theta, np.abs(a_y), label='$|\\ddot Y|$', **plt_kwargs['Y'])" ] }, { "cell_type": "markdown", "id": "e539ec54", "metadata": { "tags": [ "hide-input" ] }, "source": [ "#### Parametrization\n", "\n", "The formula above requires specifying the masses of the system and the maximum force that the motors can apply on the belts. Empirically, it is easier to measure the maximum acceleration at which the motors start skipping steps. Noticing that a force divided by an mass is an acceleration, we can write the formula in such a way that they are instead parametrized by the maximum accelerations at 0° and 90°:\n", "\n", "$$\n", "\\begin{align}\n", "a_\\text{max} &= \\frac{ F_\\text{max} }{ \\max \\left( |m_X r_X + m_Y r_Y|,\\, |m_X r_X - m_Y r_Y| \\right) } \\\\\n", "\\frac{1}{a_\\text{max}} \n", " &= \\max \\left( \n", " \\left| \\frac{m_X}{F_\\text{max}} r_X + \\frac{m_Y}{F_\\text{max}} r_Y \\right|,\\, \n", " \\left| \\frac{m_X}{F_\\text{max}} r_X - \\frac{m_Y}{F_\\text{max}} r_Y \\right|\n", " \\right) \\\\\n", " &= \\max \\left( \n", " \\left| \\frac{r_X}{\\ddot X_\\text{max}} + \\frac{r_Y}{\\ddot Y_\\text{max}} \\right|,\\, \n", " \\left| \\frac{r_X}{\\ddot X_\\text{max}} - \\frac{r_Y}{\\ddot Y_\\text{max}} \\right|\n", " \\right)\n", "\\end{align}\n", "$$\n", "\n", "#### Lowest acceleration\n", "\n", "The minimum acceleration $a_\\text{min}$ on diagonal moves is the altitude of the blue triangles. Inverse Pythagorean theorem yields:\n", "\n", "$$ a_\\text{min}^{-2} = \\ddot X_\\text{max}^{-2} + \\ddot Y_\\text{max}^{-2} $$\n", "\n", "The angle $\\beta_\\text{min}$ of this altitude is the same as the one found at the top of the blue triangle, between the 90° axis and its hypoteneuse:\n", "\n", "$$ \\beta_\\text{min} = \\arctan \\frac{\\ddot X_\\text{max}}{\\ddot Y_\\text{max}} $$\n", "\n", "For cartesian per-axis limited acceleration we got curves clamped to circles. The acceleration plot for CoreXY doesn't exactly look like that. This is because we clamped the *forces* on the belts instead of their acceleration. Then let's check that the forces are indeed clamped when applying our acceleration limits:" ] }, { "cell_type": "code", "execution_count": 9, "id": "5f63d210", "metadata": { "tags": [ "hide-input" ] }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n" ], "text/plain": [ "

" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "f_x = m_x * a_x\n", "f_y = m_y * a_y\n", "f = np.hypot(f_x, f_y)\n", "\n", "f_a = f_x + f_y\n", "f_b = f_x - f_y\n", "\n", "with fig_tmpl('CoreXY forces on belts and axes (limited acceleration)', 'force'):\n", " plt.plot(theta, f, label=r'$|F_\\mathrm{toolhead}|$', **plt_kwargs['toolhead'])\n", " plt.plot(theta, np.abs(f_a), label='$|F_A|$', **plt_kwargs['A'])\n", " plt.plot(theta, np.abs(f_b), label='$|F_B|$', **plt_kwargs['B'])\n", " plt.plot(theta, np.abs(f_x), label='$|F_X|$', **plt_kwargs['X'])\n", " plt.plot(theta, np.abs(f_y), label='$|F_Y|$', **plt_kwargs['Y'])" ] } ], "metadata": { "celltoolbar": "Tags", "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.9.7" } }, "nbformat": 4, "nbformat_minor": 5 }