Skip to content

Instantly share code, notes, and snippets.

@osamutake

osamutake/Slinky.ipynb

Last active April 25, 2019 00:49
Show Gist options

  • Select an option

    Learn more about clone URLs
  • Save osamutake/b57f0e72524c336d92e4de4fba507c43 to your computer and use it in GitHub Desktop.

Select an option

Learn more about clone URLs
Save osamutake/b57f0e72524c336d92e4de4fba507c43 to your computer and use it in GitHub Desktop.
Download ZIP
https://twitter.com/Rainmaker1973/status/1119638564269166602 を再現するためのコード。ipynb 形式にして、この系に成り立つスケール則について追記した。
Display the source blob
Display the rendered blob
Raw
Slinky.ipynb
{
"cells": [
{
"cell_type": "code",
"execution_count": 20,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"-0.3000000000004938\n"
]
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAAmgAAAGgCAYAAAAevJsNAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAAPYQAAD2EBqD+naQAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMi4zLCBodHRwOi8vbWF0cGxvdGxpYi5vcmcvIxREBQAAG2VJREFUeJzt3X+QVfV9//HXwsIq4q5FAqK4EmRRwQhCaBmNoknNlzpNaFDTanXGiYqxddI2+Q5LZjpNmnEcyRjbNJkGzVimMzRtjTBJp4nRJupKp7FDa6yJFQXKsoJLaYWI+APBPd8/+OZONrnALiy7n4XHY2ZH7vmce++58J7rc869d29DVVVVAAAoxoihPgAAAHoTaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIUpKtDefPPNPPPMM3nzzTeH+lAAAIZMUYG2fv36zJ07N+vXrz/i23jttdcG8Ig4XpgL6jEX1GMuqGew56KoQBsI77777lAfAgUyF9RjLqjHXFDPYM/FcRdoAADDnUADAChM41AfAABw4tq7d29efPHF7N+/f6gP5ZBee+21tLS0HHR95MiROeecc3LaaacNyP0JNABgSGzatCmzZ8/Onj17hvpQBsxtt92WFStWZMSIo3uRUqABAIOup6cnt9xyS8aPH59HHnkkY8aMGepDOirvvPNOnnrqqSxbtixJ8sADDxzV7Qk0AGDQdXd3p6OjI9/4xjfygQ98YKgPZ0DMnz8/SdLe3p4vfvGLR/Vypw8JAACD7n/+53+SJOeee+4QH8nAuvzyy5MkW7ZsOarbcQYNABh0PT09SZLGxv6nyLvvvpu1a9emu7s7kyZNymWXXZaRI0cO9CEekdGjRyc5+t+b5gwaADBsrFmzJlOmTsuVV16ZG264IVdeeWWmTJ2WNWvWHLP77OzszBVXXJGWlpbMnj37mN3PzxNoAMCwsGbNmlx77bXZ+iszk2Vrk7/YmSxbm22/MjPXXnvtMYu05ubm3HXXXfnGN75xTG6/HoEGABTv3XffzR/80WdSXXR1csfqZOqvJSeNTab+Wqo7VicXXZ0//PT/PaqXFu+9994sWbKkdvmnP/1pxo8fnyT5wAc+kFNOOeWoH0dfCTQAoHhr167N1q7O5DeWJb/4O8ZGjEi1sD0vb9mctWvXHvF93HrrrfnWt76Vn/70p0mSlStXZtGiRRk3btxRHPmREWgAQPG6u7sP/OHMmfV3OGtm7/2OwGmnnZZrr702f/VXf5WqqvK1r30td9555xHf3tHwKU4AoHiTJk068IdXnj/w8uYv2vZ87/2O0Kc+9al89KMfzQUXXJD3vOc9ufjii4/q9o6UM2gAQPEuu+yyTG6dkoZH7kn+/6/oqOnpScP3lufsc96byy677Kju5/zzz8/UqVOzZMmSITt7lgg0AGAYGDlyZL78Z19KnvtuGr52TbLp6eTt15NNTx+4/Nx38+f33Tsgvw/ttttuy/79+3PttdcmSd58881Mnjw51113Xf7zP/8zkydPzmc/+9mjvp9D8RInADAsLF68OA8//HD+4I8+k63LL69tn3zOe/PnDz+cxYsXD8j9PPHEE/m93/u9jBo1KkkyZsyYbN26dUBuu68EGgAwbCxevDiLFi06Jt8k8Morr+SDH/xgxo0bl0cffXQAjvbI9fslzu9+97uZM2dOZs+enQsvvDB//dd/nSTZsWNHFi5cmLa2tlx44YV56qmnatc51BoAQH+MHDkyV1xxRa6//vpcccUVA/Y1T2eeeWbWr1+ff/mXf8mpp546ILd5pPp1Bq2qqtx444158sknc9FFF6WzszPnn39+Fi9enGXLlmX+/Pn53ve+l3Xr1uVjH/tYNm/enFGjRh1yDQCA3vr9EmdDQ0PtF7jt3r07p59+epqamvLQQw9l48aNSZJ58+blzDPPTEdHR37913/9kGv17NmzJ7t3765dbmpqSlNTU78fHABQpp+d9XrnnXeG+EgG1ptvvpkkR30Sql+B1tDQkL//+7/P4sWLc8opp2TXrl1Zs2ZNXn/99ezbty9nnHFGbd8pU6akq6srr7766kHXDmbBggW9Li9dujTt7e19OsZdu3b15yFxgjAX1GMuqMdcDI6WlpacdNJJ+cIXvpA/+ZM/yejRo4f6kI7K/v37s2nTpnz2s5/NqaeemvHjx2fnzp219f5+G0G/Am3//v256667smbNmlx++eVZt25dPvrRj+bZZ5/t150eTkdHR69vi+/vGbSh+EoGymcuqMdcUI+5OPbGjRuXb3/721m0aFEeeeSRoT6cAXPFFVfk8ccfP+pfmNuvQHv22Wfzyiuv5PLLD3y0dd68eZk8eXKee+65NDY2Zvv27bUzZZ2dnWltbc3pp59+0LWDGTt2bJqbm4/0MQEAw8CHP/zhbN++PZ2dnUf1JeeD4bXXXktLS8tB10eMGJEJEybkjDPOyIhf/K7QI9CvQDv77LPT3d2dF154IRdccEE2btyYTZs25bzzzst1112XFStW5POf/3zWrVuXbdu21V6qPNQaAHDiamlpyaxZs4b6MA5r586dg3pmtV+BNnHixDzwwAP5+Mc/nhEjRqSnpydf/epX09ramuXLl+emm25KW1tbRo8enVWrVtXeIHeoNQAAeuv3pzivv/76XH/99b+0feLEiXnsscfqXudQawAA9Oa7OAEACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAArT70Dbu3dv7rzzzrS1teV973tfbrzxxiTJhg0bcskll2T69OmZN29enn/++dp1DrUGAEBv/Q60ZcuWpaGhIS+99FJ+/OMf5957702S3H777VmyZEleeumltLe35+abb65d51BrAAD01lBVVdXXnd94441MmjQpW7duTXNzc237jh07Mm3atOzcuTONjY2pqiqTJk3KP//zP6e5ufmga9OmTet1+88880zmzp2bjo6OzJ49u7a9qakpTU1NfTrGnTt3Zty4cX19SJwgzAX1mAvqMRfUM9hz0difnTdt2pRx48bl7rvvzve///2cfPLJ+fznP5/TTjstkyZNSmPjgZtraGhIa2trurq60tLSctC1Xwy0n1mwYEGvy0uXLk17e3ufjnHXrl39eUicIMwF9ZgL6jEX1HO0c9HfuOtXoO3fvz9btmzJjBkzcs899+RHP/pRrrrqqnznO9/p150eztGcQUv6/5fAicFcUI+5oB5zQT3FnkFrbW3NiBEj8ru/+7tJkosvvjjvfe97s2XLlnR3d2f//v21lzG7urrS2tqa5ubmg64dzNixY3u9hAoAcCLp14cExo8fnw996EN59NFHkySbN2/O5s2bc+mll2bOnDlZtWpVkmT16tWZPHlypk2blgkTJhx0DQCAX9avM2hJsmLFitxyyy1pb2/PiBEjcv/99+ess87K/fffn5tvvjl33313mpubs3Llytp1DrUGAEBv/Q60qVOn5oknnvil7eedd15++MMf1r3OodYAAOjNNwkAABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAURqABABRGoAEAFEagAQAU5ogDbeXKlWloaMi3vvWtJMmOHTuycOHCtLW15cILL8xTTz1V2/dQawAA9HZEgdbZ2Zmvf/3rmT9/fm3bsmXLMn/+/GzYsCErV67MDTfckH379h12DQCA3vodaD09Pbn11lvzla98JU1NTbXtDz30UD75yU8mSebNm5czzzwzHR0dh12rZ8+ePdm9e3ftZ+/evf09TACAYauxv1e47777cumll2bu3Lm1ba+++mr27duXM844o7ZtypQp6erqOuTawSxYsKDX5aVLl6a9vb1Px7dr166+PhROIOaCeswF9ZgL6jnauRg3bly/9u9XoP3kJz/J6tWrj/l7yDo6OjJ79uza5aampl5n6w6nv38JnBjMBfWYC+oxF9QzmHPRr0Bbu3ZtOjs709bWliTZvn17lixZkj/90z9NY2Njtm/fXjtT1tnZmdbW1px++ukHXTuYsWPHprm5+UgfEwDAsNav96Ddcccd6e7uTmdnZzo7OzN//vw88MADueOOO3LddddlxYoVSZJ169Zl27ZttZcqD7UGAEBv/X4P2sEsX748N910U9ra2jJ69OisWrUqo0aNOuwaAAC9HVWgPfnkk7U/T5w4MY899ljd/Q61BgBAb75JAACgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDD9CrS33347v/Vbv5Xp06dn1qxZueqqq7Jx48YkyY4dO7Jw4cK0tbXlwgsvzFNPPVW73qHWAADord9n0JYsWZIXX3wx//Ef/5FFixbl1ltvTZIsW7Ys8+fPz4YNG7Jy5crccMMN2bdv32HXAADorV+BdtJJJ+Xqq69OQ0NDkmT+/Pnp7OxMkjz00EP55Cc/mSSZN29ezjzzzHR0dBx2rZ49e/Zk9+7dtZ+9e/f2+4EBAAxXjUdz5S9/+ctZtGhRXn311ezbty9nnHFGbW3KlCnp6uo65NrBLFiwoNflpUuXpr29vU/HtGvXrn4+Ck4E5oJ6zAX1mAvqOdq5GDduXL/2P+JAu/vuu7Nx48b84Ac/yFtvvXWkN1NXR0dHZs+eXbvc1NSUpqamPl+/v38JnBjMBfWYC+oxF9QzmHNxRIF27733Zs2aNfn+97+fMWPGZMyYMWlsbMz27dtrZ8o6OzvT2tqa008//aBrBzN27Ng0NzcfyaEBAAx7/f6QwH333Ze//du/zT/90z/ltNNOq22/7rrrsmLFiiTJunXrsm3bttpLlYdaAwCgt36dQdu6dWs+85nPZOrUqbnyyiuTHHj58V//9V+zfPny3HTTTWlra8vo0aOzatWqjBo1KkkOuQYAQG/9CrTJkyenqqq6axMnTsxjjz3W7zUAAHrzTQIAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhRFoAACFEWgAAIURaAAAhWkc6gMYKA0NDcmIxqTqSRpGJD37U1XVUB8WAEC/DdoZtA0bNuSSSy7J9OnTM2/evDz//PMDdtsNDQ0H/nDeguTaew789+e3AwAMI4MWaLfffnuWLFmSl156Ke3t7bn55psH7sZHNCYXfCj5w+8mV/3hgf+e/8Fk5KiBuw8AgEEyKIG2Y8eO/Nu//VtuvPHGJMk111yTl19+ORs3bhyYO6h6kgv/T/KzM2YNDcn7FiY97w7M7QMADKJBCbSXX345kyZNSmPjgbe8NTQ0pLW1NV1dXXX337NnT3bv3l372bt376HvoGFE8pNHk5+956yqkh9/78B2AIBhpsgPCSxYsKDX5aVLl6a9vf3gV+jZn7zwg+TPfuPAmbMffy9Z/3iSZOfOncfyUBkmdu3aNdSHQIHMBfWYC+o52rkYN25cv/YflEA7++yz093dnf3796exsTFVVaWrqyutra119+/o6Mjs2bNrl5uamtLU1HTQ2z9nynuzpXNzsmFt8uKTyYiRSZIpU6f2+y+E45dZoB5zQT3mgnoGcy4G5TXACRMmZM6cOVm1alWSZPXq1Zk8eXKmTZtWd/+xY8emubm59nOoOEuS+75074E/NJ1y4IMBTacc2H7vvQP3IAAABsmgvUnr/vvvz/3335/p06fnnnvuycqVKwfsthcvXpzVq1dn1nnnpmnUyMw679ysWbMmH/vYxwbsPgAABsugvQftvPPOyw9/+MNjdvuLFy/O4sWLs3PnTqemAYBhzcccAQAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACiPQAAAKI9AAAAoj0AAACnNcBdrevXuzfPny7N27d6gPhYKYC+oxF9RjLqhnKOaioaqqatDu7TCeeeaZzJ07N//+7/+eOXPm9Pv6u3fvTktLS1577bU0NzcfgyNkODIX1GMuqMdcUM9QzMVxdQYNAOB4INAAAArTONQH8PPeeuutJMkLL7xwRNffs2dPkuTZZ5/N2LFjB+y4GN7MBfWYC+oxF9QzUHNx/vnnZ8yYMX3at6j3oP3N3/xNbrzxxqE+DACAAdef99gXFWj/+7//m0cffTRTpkzJySefPNSHAwAwYIbtGTQAAHxIAACgOAINAKAwwybQNmzYkEsuuSTTp0/PvHnz8vzzz9fd78EHH0xbW1vOPffc3Hbbbdm3b1+f1hh++jITjz/+eH71V381M2bMyMyZM7N06dL09PQkSTo7OzNy5MjMnj279rNp06bBfhgMsL7MxZNPPpmTTz6517/9zz5FnniuOB71ZS5WrlzZaybGjx+fxYsXJ/F8cbz61Kc+lSlTpqShoSHPPvvsQfcbkraohokrr7yyWrlyZVVVVfXNb36zev/73/9L+/zXf/1XNWnSpKq7u7vq6empPvKRj1Rf/epXD7vG8NSXmXjmmWeqTZs2VVVVVW+99VZ16aWX1q6zefPmqqWlZbAOl0HSl7l44oknqlmzZtW9vueK41Nf5uIXzZw5s3r44YerqvJ8cbzq6OioXn755eqcc86pfvSjH9XdZ6jaYlgE2n//939Xp556arVv376qqqqqp6enmjhxYrVhw4Ze+33xi1+sbr/99trl73znO9Wll1562DWGn77OxC/6/d///epzn/tcVVWecI9HfZ2LQwWa54rjz5E8Xzz99NPVe97znuqdd96pqsrzxfHuUIE2VG0xLF7ifPnllzNp0qQ0Nh74vboNDQ1pbW1NV1dXr/26urpyzjnn1C5PmTKlts+h1hh++joTP2/79u15+OGH85u/+Zu1bW+88UbmzZuXOXPm5Atf+ELefffdY37sHDv9mYtNmzZlzpw5mTdvXv7yL/+ytt1zxfHnSJ4vHnzwwdx0000ZNWpUbZvnixPTULVFUd8kAMfK7t2785GPfCRLly7N+9///iTJpEmTsm3btkyYMCE7d+7Mb//2b+dLX/pSli5dOsRHy7E2Z86cbN26NS0tLdm6dWuuvvrqjB8/Ph//+MeH+tAowBtvvJG/+7u/y9NPP13b5vmCwTYszqCdffbZ6e7uzv79+5MkVVWlq6srra2tvfZrbW3Nli1bapc7Oztr+xxqjeGnrzORJK+//noWLlyYRYsW5dOf/nRte1NTUyZMmJAkGTduXD7xiU9k7dq1g/MAOCb6OhfNzc1paWlJkkyePDnXX3997d/ec8Xxpz/PF0nyzW9+MzNnzsyMGTNq2zxfnLiGqi2GRaBNmDAhc+bMyapVq5Ikq1evzuTJkzNt2rRe+11zzTX5h3/4h2zfvj1VVWXFihX5nd/5ncOuMfz0dSb27NmThQsXZuHChfnjP/7jXms7duyofdpm7969WbNmTS6++OLBeQAcE32di+7u7tqneV9//fX84z/+Y+3f3nPF8aevc/EzDz74YG655ZZe2zxfnLiGrC0G5J1sg2D9+vXV/Pnzq7a2tmru3LnVc889V1VVVd1yyy3Vt7/97dp+DzzwQDV16tRq6tSp1Sc+8YnaGzwPt8bw05eZuOuuu6rGxsZq1qxZtZ+77rqrqqqqWr16dTVz5szqoosuqmbMmFHdeeed1dtvvz1kj4eB0Ze5+MpXvlLNmDGj9m//uc99rurp6andhueK409f/x+yfv36auzYsdXu3bt7Xd/zxfFpyZIl1VlnnVWNHDmymjBhQnXuuedWVVVGW/iqJwCAwgyLlzgBAE4kAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDACDQCgMAINAKAwAg0AoDD/Dy+fBzN599O2AAAAAElFTkSuQmCC",
"text/plain": [
"PyPlot.Figure(PyObject <Figure size 600x400 with 1 Axes>)"
]
},
"metadata": {},
"output_type": "display_data"
},
{
"name": "stderr",
"output_type": "stream",
"text": [
"┌ Info: Saved animation to \n",
"│ fn = /mnt/juliabox/tmp.gif\n",
"└ @ Plots /home/jrun/.julia/packages/Plots/UQI78/src/animation.jl:90\n"
]
},
{
"data": {
"text/html": [
"<img src=\"tmp.gif\" />"
],
"text/plain": [
"Plots.AnimatedGif(\"/mnt/juliabox/tmp.gif\")"
]
},
"execution_count": 20,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using Statistics\n",
"using Plots \n",
"pyplot(reuse=true)\n",
"\n",
"# 質点個数, ばね定数, 質量, 重力加速度, 計算間隔\n",
"n=40; k=1.0; m=1.0; g=1.0; dt=0.2\n",
"\n",
"u = [(m*g/k)*i for i in 1:(n-1) ] # 質点間隔は下にある点の重さで決まる\n",
"y = [sum(u[1:j-1]) for j in 1:n] # 初期質点位置\n",
"ylimit = y[n]\n",
"v = [0.0 for i in 1:n] # 初期質点速度\n",
"zero = [0 for j in 1:n] # x座標\n",
"\n",
"function detect_crush() # 衝突検出\n",
" for i in 1:(n-1)\n",
" if y[i] > y[i+1]\n",
" mean_y = mean(y[i:n]) # ここより上は全部つぶれているので\n",
" mean_v = mean(v[i:n]) # それらの重心の運動に置き換える\n",
" for j in i:n\n",
" y[j] = mean_y\n",
" v[j] = mean_v\n",
" end\n",
" detect_crush() # さらにその下を追い越した可能性が\n",
" break # あるので再度チェックする\n",
" end\n",
" end\n",
"end\n",
"\n",
"@gif for t in 1:120\n",
" for i in 1:n\n",
" v[i] += dt * ( i == 1 ? -g + (k/m)*(y[i+1]-y[i]) :\n",
" i == n ? -g - (k/m)*(y[i]-y[i-1]) :\n",
" -g - (k/m)*(y[i]-y[i-1]) + (k/m)*(y[i+1]-y[i]) )\n",
" end\n",
" for i in 1:n\n",
" y[i] += dt * v[i]\n",
" end\n",
" detect_crush()\n",
"\n",
" plot(zero, y, seriestype=:scatter, ylims=(-ylimit*0.1,ylimit*1.1))\n",
"end"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
"この系で成り立つスケール則について\n",
"====\n",
"\n",
"全質量 $m$ のバネを、$N-1$ 個のバネで連結された $N$ 個の質点でモデル化する。\n",
"\n",
"- 質点１つあたりの質量は $m/N$\n",
"- バネ全体のばね定数を $k$ とすると、１区間当たりのばね定数は $(N-1)k$\n",
"- 個別の質点及びバネに下から $1,2,\\dots$ の番号を付ける\n",
"\n",
"初期状態で $i$ 番目のバネに $im/N$ の質量がぶら下がる。\n",
"ここではバネの自然長はゼロなので、この区間のバネの長さは、\n",
"\n",
"$\\Delta l_i=\\frac{i(m/N)g}{(N-1)k}=\\frac{img}{N(N-1)k}$\n",
"\n",
"となる。全長は、\n",
"\n",
"$\\begin{aligned}\n",
"l_0&=\\sum_{i=1}^{N-1}\\Delta l_i\\\\\n",
"&=\\sum_{i=1}^{N-1} \\frac{img}{N(N-1)k}\\\\\n",
"&=\\frac{\\cancel{N(N-1)}}{2}\\frac{mg}{\\cancel{N(N-1)}\\,k}\\\\\n",
"&=\\frac{mg}{2k}\n",
"\\end{aligned}$\n",
"\n",
"$i$ 番目の質点の $y$ 座標を $y_i$ とすると、初期位置は\n",
"\n",
"$\\begin{aligned}\n",
"y_i(t=0)&=\\sum_{j=1}^{i-1}\\Delta l_j=\\frac{i(i-1)}{N(N-1)}l_0\\\\\n",
"\\end{aligned}$\n",
"\n",
"にあるから、重心位置は\n",
"\n",
"$\\begin{aligned}\n",
"\\overline y&=\\frac{1}{N}\\sum_{i=1}^Ny_i\\\\\n",
"&=\\frac{l_0}N\\sum_{i=1}^N\\frac{i(i-1)}{N(N-1)}\\\\\n",
"&=\\frac{l_0}{3}\\frac{N+1}{N}\\sim\\frac{l_0}{3}\n",
"\\end{aligned}$\n",
"\n",
"運動方程式は、\n",
"\n",
"$$\n",
"\\begin{aligned}\n",
"\\frac mN\\frac{d^2y_i}{dt^2}\n",
"&=-\\frac mNg+\\frac k{N-1}\\Delta l_i-\\frac k{N-1}\\Delta l_{i-1}\\\\\n",
"\\end{aligned}\n",
"$$\n",
"\n",
"$N-1\\sim N$ とみなせるとき、\n",
"\n",
"$$\n",
"\\begin{aligned}\n",
"\\frac{d^2y_i}{dt^2}\n",
"&=-g+\\frac km(y_{i+1}-2y_i+y_{i-1})\\\\\n",
"\\end{aligned}\n",
"$$\n",
"\n",
"距離を $l_0$ を単位として測るなら、$y_i=l_0Y_i$ として、\n",
"\n",
"$$\n",
"\\begin{aligned}\n",
"\\frac{d^2Y_i}{dt^2}\n",
"&=-\\frac{2k}{m}+\\frac{k}{m}(Y_{i+1}-2Y_i+Y_{i-1})\\\\\n",
"\\end{aligned}\n",
"$$\n",
"\n",
"時間を $t_0=\\sqrt{m/k}$ を単位として測るなら $t=t_0T$ として、\n",
"\n",
"$$\n",
"\\begin{aligned}\n",
"\\frac{d^2Y_i}{dT^2}\n",
"&=-2+(Y_{i+1}-2Y_i+Y_{i-1})\\\\\n",
"\\end{aligned}\n",
"$$\n",
"\n",
"のように無次元化される。初期条件\n",
"\n",
"$$\\begin{aligned}\n",
"Y_i(t=0)&=\\frac{i(i-1)}{N(N-1)}\\\\\n",
"\\end{aligned}\n",
"$$\n",
"\n",
"にも物理パラメータは入らないので、\n",
"このスケールでは物理パラメータによらず同じ現象が観測されることになる。\n",
"\n",
"固体物理などの問題設定では $Y_{i+1}-2Y_i+Y_{i-1}$ \n",
"の部分を空間微分で書き直せるのだが、今の場合は $Y_{i+1}$ と $Y_i$ \n",
"との距離は近似的にも等間隔ではないため、空間微分への書き直しができない。\n"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": []
}
],
"metadata": {
"kernelspec": {
"display_name": "Julia 1.0.3",
"language": "julia",
"name": "julia-1.0"
},
"language_info": {
"file_extension": ".jl",
"mimetype": "application/julia",
"name": "julia",
"version": "1.0.3"
}
},
"nbformat": 4,
"nbformat_minor": 2
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment