anisotropy/data/analyze.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from anisotropy.database import Database, tables\n",
    "import pathlib\n",
    "import peewee as pw\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "db = Database(pathlib.Path(\"anisotropy.db\").resolve())\n",
    "savefig = False\n",
    "execution = 5"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "def load_data(execution: int, field: str):\n",
    "    if not db.getExecution(execution):\n",
    "        print(\"Execution not found\")\n",
    "\n",
    "    for model in db.tables:\n",
    "        try:\n",
    "            column = getattr(model, field)\n",
    "        \n",
    "        except AttributeError:\n",
    "            pass\n",
    "\n",
    "        else:\n",
    "            break\n",
    "\n",
    "    query = model.select(tables.Shape.alpha, column, tables.Shape.direction, tables.Shape.label)\n",
    "    idn = db.tables.index(model)\n",
    "\n",
    "    for table in reversed(db.tables[ :idn]):\n",
    "        query = query.join(table, pw.JOIN.LEFT_OUTER)\n",
    "        \n",
    "    query = query.switch(tables.Shape)\n",
    "    query = query.where(\n",
    "        tables.Shape.exec_id == execution,\n",
    "    #    tables.Shape.label == structure,\n",
    "    )\n",
    "    query = query.order_by(tables.Shape.label, tables.Shape.direction, tables.Shape.alpha)\n",
    "\n",
    "    with db:\n",
    "        if query.exists():\n",
    "            table = []\n",
    "            for row in query.dicts():\n",
    "                for k in row.keys():\n",
    "                    if type(row[k]) == list:\n",
    "                        row[k] = str(row[k])\n",
    "\n",
    "                table.append(row)\n",
    "        \n",
    "        else:\n",
    "            table = None\n",
    "\n",
    "    if table is None:\n",
    "        print(\"Results not found\")\n",
    "\n",
    "    else:\n",
    "        return pd.DataFrame(table)\n",
    "\n",
    "\n",
    "def nanmean(arr):\n",
    "    temp = arr.copy()\n",
    "\n",
    "    for n, item in enumerate(temp):\n",
    "        if np.all(np.isnan(item)):\n",
    "            \n",
    "            vals = temp[n - 1 : n + 2]\n",
    "\n",
    "            if np.sum(~np.isnan(vals)) <= 1:\n",
    "                vals = temp[n - 2 : n + 3]\n",
    "\n",
    "            temp[n] = vals[~np.isnan(vals)].mean()\n",
    "\n",
    "    return temp\n",
    "\n",
    "def filter_group(arr, nan = True, qhigh = True):\n",
    "    temp = arr.copy()\n",
    "    check = True\n",
    "    quan = np.quantile(temp[~np.isnan(temp)], 0.97)\n",
    "    limit = 1000\n",
    "\n",
    "    while check:\n",
    "        if nan and np.any(np.isnan(temp)):\n",
    "            temp = nanmean(temp)\n",
    "            check = True\n",
    "        \n",
    "        elif qhigh and np.any(quan < temp):\n",
    "            temp[quan < temp] = np.nan\n",
    "            check = True\n",
    "\n",
    "        else:\n",
    "            check = False \n",
    "        \n",
    "        if limit <= 0:\n",
    "            break\n",
    "\n",
    "        else:\n",
    "            limit -= 1\n",
    "\n",
    "    return temp"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Porosity"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/tmp/ipykernel_53258/3263217096.py:3: UserWarning: Boolean Series key will be reindexed to match DataFrame index.\n",
      "  simple = df[df.label == \"simple\"][df.direction == \"[0.0, 0.0, 1.0]\"][\"porosity\"].to_numpy()\n",
      "/tmp/ipykernel_53258/3263217096.py:4: UserWarning: Boolean Series key will be reindexed to match DataFrame index.\n",
      "  bodyCentered = df[df.label == \"bodyCentered\"][df.direction == \"[0.0, 0.0, 1.0]\"][\"porosity\"].to_numpy()\n",
      "/tmp/ipykernel_53258/3263217096.py:5: UserWarning: Boolean Series key will be reindexed to match DataFrame index.\n",
      "  faceCentered = df[df.label == \"faceCentered\"][df.direction == \"[0.0, 0.0, 1.0]\"][\"porosity\"].to_numpy()\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYgAAAEGCAYAAAB/+QKOAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAA6iklEQVR4nO3dd3xUVfr48c8zk0kChIRQDEKAIAEUkGIQpAqKioUigm1FdEVsrO6i7urXFV3Lb21rW1nFtmBbBBWNCNhIUFQUEJQOoUhHCDVA+vP74w4xYIBAcnMnk+f9es2LW+c+h4F55txzzzmiqhhjjDGH83kdgDHGmNBkCcIYY0yJLEEYY4wpkSUIY4wxJbIEYYwxpkQRXgdQXurWratJSUlF6/v27aNGjRreBVRBqko5wcoajqpKOSF0yzpv3rztqlqvpH1hkyCSkpKYO3du0Xp6ejq9evXyLqAKUlXKCVbWcFRVygmhW1YR+eVI++wWkzHGmBJZgjDGGFMiSxDGGGNKFDZtEMYYc6Ly8vLYsGED2dnZrl0jLi6OpUuXuvb+xxIdHU1iYiKBQKDU57iaIESkL/Ac4AdeVdXHjnDcZcB7wJmqOje47V7gBqAAuF1VP3UzVmNM1bVhwwZq1qxJUlISIuLKNfbu3UvNmjVdee9jUVUyMzPZsGEDTZs2LfV5rt1iEhE/MAa4EGgFXCUirUo4riZwB/B9sW2tgCuB1kBf4D/B9zPGmHKXnZ1NnTp1XEsOXhMR6tSpc9w1JDfbIDoBGaq6WlVzgQnAgBKOexh4HCge+QBggqrmqOoaICP4fsYY44pwTQ4HnUj53LzF1BBYX2x9A9C5+AEicgbQSFU/EZG7Dzt39mHnNjz8AiIyAhgBkJCQQHp6etG+rKysQ9aPJqdA2bZfSaxZ+drsj6eclZ2VNfyESjnj4uLYu3evq9coKChw/RrHkp2dfVx/3541UouID3gauO5E30NVXwZeBujYsaMW74RyPJ1SXpq5iie+XcaXd/aiad3Q6+l4NKHa+cYNVtbwEyrlXLp0qevtA8dqg4iJiSErKwuAzZs306NHD5555hn69etHUlISNWvWxO937rRnZWVx6qmnMmXKlOOKITo6mg4dOpT6eDcTxEagUbH1xOC2g2oCbYD0YNWnPpAqIv1LcW65uqJjI2pVCxQlh28ytnN6Yhyx0aVv7TfGmPKwd+9eLrroIv72t7/Rr1+/ou1paWnUrVsXcBLrU0895Xosbt5TmQM0F5GmIhKJ0+icenCnqu5W1bqqmqSqSTi3lPoHn2JKBa4UkSgRaQo0B35wK9D4GpFc2akxAFk5+dz05jz+kbrErcsZY0yJ8vLyGDRoEP379+fGG2/0Ohz3EoSq5gMjgU+BpcBEVV0sIg8FawlHO3cxMBFYAkwHblPVArdiLS4mKoIJI87ijnObA/Dr3mw+XbwFm5rVmKrjirHfMWmu04SaV1DIFWO/Y/L8DQAcyC3girHf8fFPmwDYk53HFWO/Y/qizQDs2JfLFWO/44slWwHnO6S0/vjHPzJz5kyuuuqq8izOCXO1VVZVp6pqC1VtpqqPBreNVtXUEo7tdbAPRHD90eB5LVV1mptxHq5Nwzga16kOwBvf/sLId35k0273OtAYY8y+ffvIzMxk3Lhx3HbbbV6HA1hP6mP6c5/m9GpZj4a1qgHw/rwN9GhRl5NqRnscmTHGLe/e1KVoOeD3HbJeLdJ/yHpsdOCQ9do1Ig9ZL+13RVRUFJMmTaJGjRqMHz+et99+mz/84Q9lKUaZVb7nOitYhN9Hx6TaAGzbm8O9kxfy2tdrPI7KGBNuIiIiiuaLGDNmDPfddx+7d+/2NCZLEMehXs0oPv1zT27tnQzAsi17mLpws7VPGGPKVXJyMtdffz3/93//52kcdovpOBXvJzH+21+Ytmgz3ZvXtUdijTFlcrAPxEEPPPBA0fLatWsP2derV68K6T9iNYgyeHhAaybe1IXY6ACqyitfrWbb3hyvwzLGmHJhCaIMIvw+WiQ4PSMzfs3i8enLmBZ81M0YYyo7u8VUTpon1OTzUWfTKN552unbjO0cyCvgnFNPCvtBwIwx4ckSRDkq3j7x2qw1rMncR88W9Qj4LUEYYyofSxAueWloCpt3ZRPw+8grKOTfMzK4rmsStWtEeh2aMcaUirVBuCTg9xX1xp6zdgdj0jJYsH6nx1EZY0zpWYKoAF2b1SX9rl70bnkSAB/O38j0RdZ/whjzm5iYmKLlzZs3k5yczMcffwzA7t27ufbaa0lOTqZZs2Zce+21RZ3o1q5dS7Vq1Wjfvn3Rq379+uUy2qsliArSqHZ1RARV5Z0f1vH29+u8DskYE4JKGu77hhtu4JRTTiEjI4NVq1bRtGlThg8fXnROs2bNWLBgQdHr5ptvLpdYrA2igokI7wzvzO4DeYgIu/fnMSY9g1t7NaNWdWufMKYqK2m474yMDObNm8e7775bdNzo0aNJTk5m1apVRZMIucEShAci/D7qxEQBMCtjO+O+WculHRpagjAmVPz34t9vaz0QOt0Iufvh7SG/39/+aujwB9iXCROvPXTf9Z+U6rIHh/v+97//XbRtyZIltG/f/pBE4Pf7ad++PYsXL6Zt27aleu8TYbeYPHZx25OZdU9vTjs5FoBx36xh3i87PI7KGFPRbLhvU6KDwwFn5xXw2jdr6NasLilNansclTFV2NF+8UdWP/r+GnVKXWMo7kjDfbdq1YoFCxZQWFiIz+f8pi8sLGTBggW0atXquK9zPKwGEUKiA34+/XNP7r3oNADWbN/HK1+tJq+g0OPIjDFuO9Jw38nJyXTo0IFHHnmk6NhHHnmEM844g+TkZFdjsgQRYqpHRhBXzRkZ9uOfNvHclyvZuT/X46iMMRXp8OG+X3vtNVasWEGzZs1o1qwZK1as4LXXXnM9DldvMYlIX+A5wA+8qqqPHbb/ZuA2oADIAkao6hIRScKZx3p58NDZqlo+z21VIref25xBZzQsugU1duYqLmxzclEHPGNM+DjacN/x8fG89dZbJZ6XlJTEokWLDtn24IMPlktMriUIEfEDY4DzgA3AHBFJVdUlxQ57R1VfCh7fH3ga6Bvct0pV27sVX2WRGO8kg027DvDclyspUOXWXu5WK40xBtytQXQCMlR1NYCITAAGAEUJQlX3FDu+BmBdi4+gQa1qzLizF7WqO7ef5qzdwYad+6llvbGNMS5xM0E0BNYXW98AdD78IBG5DRgFRALnFNvVVETmA3uAv6vq1yWcOwIYAZCQkEB6enrRvqysrEPWw80rP+ewbEcB93UoDOtyFhfun2lxVaWsoVLOuLg49u7d6+o1CgoKXL/GsWRnZx/f37equvICBuO0OxxcHwq8cJTjrwbGB5ejgDrB5RScRBN7tOulpKRocWlpaRrOCgoKdf2OfZqWlqb5BYX6woyVunNfjtdhuSrcP9PiqkpZQ6WcS5Yscf0ae/bscf0ax1JSOYG5eoTvVTefYtoINCq2nhjcdiQTgIEAqpqjqpnB5XnAKqCFO2FWTj6fFLVPzF+3k399tpyvV273OCpjTDhxM0HMAZqLSFMRiQSuBFKLHyAizYutXgysDG6vF2zkRkROAZoDq12MtVLrmFSbL0adzSVtTwYgbfmvLNq42+OojDGVnWttEKqaLyIjgU9xHnN9XVUXi8hDOFWaVGCkiPQB8oCdwLDg6T2Bh0QkDygEblZVG3/iKE6p5wwVrKo8Pm0ZMVERTLq5i013aow5Ya52lFPVqaraQlWbqeqjwW2jg8kBVb1DVVurantV7a2qi4Pb3y+2/QxV/djNOMOJiPDuTV341+XtEBH25eQz4Yd1FBTa007GhLKjzQeRlJTE6aefXjTfw+jRowHo1asXc+fOLTpv7dq1tGnTptxisrGYwlBctUBRb+wPF2zkvsmLaNMwjjYN4zyOzBhzLCXNBwGQlpZG3bp1KzQWSxBh7upOjWnd4LfkMH3RZjo0jichNtrjyIwJXddPv/6Yx5ydeDbXtbmu6PgByQMYmDyQndk7GZU+6pBj/9v3v6W6bknzQXjJxmIKcyJC+0a
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "df = load_data(execution, \"porosity\")\n",
    "alpha = df[\"alpha\"].unique()\n",
    "simple = df[df.label == \"simple\"][df.direction == \"[0.0, 0.0, 1.0]\"][\"porosity\"].to_numpy()\n",
    "bodyCentered = df[df.label == \"bodyCentered\"][df.direction == \"[0.0, 0.0, 1.0]\"][\"porosity\"].to_numpy()\n",
    "faceCentered = df[df.label == \"faceCentered\"][df.direction == \"[0.0, 0.0, 1.0]\"][\"porosity\"].to_numpy()\n",
    "\n",
    "fig, ax = plt.subplots(nrows = 1, ncols = 1)\n",
    "ax.plot(alpha, np.pad(simple, (0, alpha.size - simple.size), 'constant', constant_values = np.nan), \":\", label = \"КП\")\n",
    "ax.plot(alpha, np.pad(bodyCentered, (0, alpha.size - bodyCentered.size), 'constant', constant_values = np.nan), \"--\", label = \"КОП\")\n",
    "ax.plot(alpha, np.pad(faceCentered, (0, alpha.size - faceCentered.size), 'constant', constant_values = np.nan), \"-.\", label = \"КГЦ\")\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.xlabel(r\"$\\alpha$\")\n",
    "plt.ylabel(r\"$m$\")\n",
    "plt.show()\n",
    "\n",
    "if savefig:\n",
    "    fig.savefig(\"porosity-rounded.png\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Simple structure"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEWCAYAAABrDZDcAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAqCUlEQVR4nO3de5gcZZn+8e8NBKIkHCSQJZwCKCqyKiR4QH5LIiqngK4iC6ssIIpBLnddXBYFVhTxACKKIkZUQBQJLCAiCgiYhJPoEiCAxEgMpwDhECRkgiQheX5/VE1oZrp6qma6+jB1f66rr+muqq5+3u6e9+n3UFWKCMzMrLrWancAZmbWXk4EZmYV50RgZlZxTgRmZhXnRGBmVnFOBGZtImlEu2MwAycCs5aRtImkcyQ9IGkxcGW7YzIDWKfdAdjwIukhYCywqmbxOsAdEbF7W4LqAOmv/xnADcCuEfFceyMye5lbBFaG/SNiVO8NmNrugDrAR4DHIuJYJwHrNE4E1nKSHpL0eUn3S/qbpPMljaxZ/wlJ8yU9K+kqSeNq1r1J0vXpuiclnSDpnZJ60ttKSStqHm8taXtJv5O0WNIzki6StFHNPneXdLekpelzVkualBH7TEkfT+9vL+lRSQdIWjeN6R9rtt1M0guSNgXeBqyS9CdJS9L9vLHPe/Ke9P6otGy3pI8nSVpYs23fx+MkXS7paUkPSvr3mnVrp+/RX9PyzZa0VbouJL02vX+ApEckbZc+3k/SXZKeT8v4xT7vwyclPZy+X8sk+RQFXcyJwNrlI8BewPbADsBJAJLeDXwNOAjYHHgYmJ6uG03StXItMA54LXBjRPy+pvVxEXB6TYvkEUDpPscBbwS2Ar5YE8sZwC+ADdJ9PD5Q8JI2B64DToiIqyJiRRrnR2s2OySN72ng1cB7gc8AY4DfAL+StG6d3R8HrKx5vJqM/1VJawG/AuYAWwB7Ap+RtFe6ybFpHPsCGwAfA17os489gGnAfhGxIF28DPg3YCNgP+BoSR9It18fOAc4LH2/3lIvNuseTgTWLmdHxKMR8SzwFZLKCpIEcV5E3BkRy4HPA++UNB6YAiyKiG9GxIsRsTQi/jDQC0XE/Ii4PiKWp5XymcAefTZbmyRh5LExSRK4KCJ+WrP8J8Ahknr3cyhQu/6XaRwrSZLPq4Ddancs6R+AI9MYez0KbCapXoW7K7BpRJwSESvSivyHwMHp+o8DJ0XEvEjMiYjFNc/fGbgK+EhE3Nu7MCJmRsS9EbE6Iu4BLubl92wtkuTkMcZhwonA2uXRmvsPk/xaJ/37cO+KiOgBFpP82t0K+GvRF5I0VtJ0SY9Jeh74Gcmv8l6fBg4AXpT0XE0sWb4E9ADvrqn0SZPSC8AkSW8gabFcla5e3qdcq0negy367Ptk4LvAszXbPgicAlyfxnd1zfbbAOMkPdd7A04gGbCHgd+zHwEPkLRW1pD0dkkz0u6mJSTjPGPSeJaSJKsLJb0A3Nlg/9YFnAisXbaqub81L3fHPE5SuQFruiE2AR4jqTi3G8RrfRUI4B8jYgOS7pvaCvz/SCrpEyNiIwbuGroU2D3dx6f7rPtJuv9Dgcsi4sV0+SN9yiWS9+CxmufuQNJddlbfF0x/8W+WxjelZtWjwIMRsVHNbXRE7FuzfvsGZflMur8jJe1Ss/znJElsq4jYkKTrqLbFdCVJ99VkoPZ51oWcCKxdjpG0paTXACcCl6TLLwaOkPRWSeuRVOJ/iIiHSH4Jby7pM5LWkzRa0ttzvNZokl/wSyRtQdIHv4akg0iS0bdyxn5L+ov+Y8AXegdYUz8D/pkkGVxYs/xSYD9JeyqZSvpZklbCbTXbnAScUpM88vgjsFTS8ZJelQ4O7yRp13T9j4AvS3qdEm+WtEnN82+OiEXAfwHn6+WD3EYDz0bEi5LeBvxrn9f9OnBVnq4563xOBNYuPwd+Cywg6bo4FSAibgD+B7gceILk1+zB6bqlJF0Y+wOLSLo0Jud4rS+R/GpdAvwauKJ3haSNSRLAJyLipSIFiIi/kFSIP+rtIoqIR0m6SgK4uWbbB0haCd8FnknLsH86yNzrGV6ZPPLEsIrkF/1bgQfTffwI2DDd5EySJPRb4HngxyRjE33381OS1sMJ6aJPAadIWgp8Id0HAJLeRTKAfELf/Vh3ki9MY62m5KCzj6eV/rAj6Tzg8Yg4qd2xmOXhUX+zJkpnN32QZDaOWVdw15BZk0j6MnAf8I10po9ZV3DXkJlZxblFYGZWcV03RjBmzJgYP378msfLli1j/fXXb19ALVKVckJ1ylqVckJ1ytrJ5Zw9e/YzEbFpvXVdlwjGjx/PHXfcsebxzJkzmTRpUvsCapGqlBOqU9aqlBOqU9ZOLqekh7PWuWvIzKzinAjMzCrOicDMrOK6boygnpUrV7Jw4UJefLHIKVq6y4YbbsjcuXMH3G7kyJFsueWWjBjh66KbWT7DIhEsXLiQ0aNHM378eGrOCjysLF26lNGjRzfcJiJYvHgxCxcuZNttt21RZGbW7YZFInjxxReHdRLISxKbbLIJTz/9dLtDMbMmmXjq9TzTs6Lf8jGj1uWOk95b5xnFDZsxgqongV5+H8yGl3pJoNHywRg2icDMzAbHicDMrOKcCMzMKs6JoIluvPFGDj300MLPmzp1KrfeemvT9mdmVsSwmDVURJkj8HPmzGHnnYtfj+T222/ne9/7XtP2Z2bDx5hR62bWWc1SuURQ5gj8nDlzOPzww1m+fDmf/OQnGTduHF/5ylcazuSZO3cuO+ywA2uvvXbD/U2dOpVtttlmwP2Z2fDSrCmijbhrqInuueceNttsM/baay/e85738NWvfnXASvuaa65h7733HnB/kyZNyrU/M7OiKtciKMvKlStZsGABhxxyCD/4wQ945zvfyZVXXsmvf/1rnn/+eY488kje97739Xveddddx/nnnz/g/nbaaScAFixYwFe+8hWWLFnCZZddVnq5zGz4c4ugSebOncuuu+7KOuuss6ab5wMf+AA//OEPmTZtGpdcckm/57zwwgs899xzjBs3Ltf+ALbbbjt+/OMfl1cQM6uc0hKBpK0kzZB0v6Q/SfqPOttI0nckzZd0j6RdyoqnbHPmzGG33XZj+vTpHHHEETz55JNr1p166qkcc8wx/Z4zY8YMJk+enGt/Tz31VGmxm1m1ldkieAn4bETsCLwDOEbSjn222Qd4XXo7Cvh+ifEA2SPtQx2BnzNnDjvttBM77LADp512GgcddBArV67k+OOPZ5999mGXXfrnuEbjA333d9hhh7Fy5cohxWhmVk9pYwQR8QTwRHp/qaS5wBbA/TWbvR+4MCICuF3SRpI2T59birJG4M8444w196dMmcKUKVP4zne+ww033MCSJUuYP38+U6dOfcVzbrvtNr71rW/l2t8ee+zBiBEjWLx4MSeeeCJ33XUXX/va1/j85z9fSnnMrDqU1MElv4g0HrgJ2Ckinq9ZfjXw9Yi4JX18I3B8RNzR5/lHkbQYGDt27ITp06evWdfT08MWW2zBa1/72tLL0U6rVq2qO8W0nvnz57NkyZKSIypPT08Po0aNancYpatKOaE6Ze3kck6ePHl2REyst670WUOSRgGXA5+pTQJFRMS5wLkAEydOjNqLQ8+cOZORI0cOeK7+bpfnegS9Ro4c2dUHonXyBcCbqSrlhOqUtVvLWeqsIUkjSJLARRFxRZ1NHgO2qnm8ZbrMzMxapMxZQwJ+DMyNiDMzNrsK+Ld09tA7gCVljg+YmVl/ZXYNvQs4FLhX0t3pshOArQEiYhrwG2BfYD7wAnBEifGYmVkdZc4augVoeD6EdLZQ/wn2ZmbWMj6y2Mys4pwIzMwqzonAzKzinAjMzCrOiaCJfKlKM+tG1U0ESxfB+fvA0icH3janoVyq8h3veEfT9mdmVkR1E8Gs0+GR22HWaU3bZW/FvXz5cg4//HBOOOEEBjqX00CXquzd39SpU3Ptz8ysqGomgqWL4O6LIFYnf5vUKvClKs2sGw14QJmke4Han6EiORbszaVFVbZZpydJAJK/s06DKVlnwcin3qU
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "df = load_data(execution, \"flowRate\")\n",
    "simple = df[df.label == \"simple\"].groupby(df.direction)\n",
    "alpha = simple.get_group(\"[0.0, 0.0, 1.0]\")[\"alpha\"].to_numpy()\n",
    "\n",
    "anisotropy_21 = 2 * simple.get_group('[0.0, 0.0, 1.0]')[\"flowRate\"].to_numpy() / simple.get_group('[1.0, 0.0, 0.0]')[\"flowRate\"].to_numpy()\n",
    "anisotropy_31 = 2 * simple.get_group('[1.0, 1.0, 1.0]')[\"flowRate\"].to_numpy() / simple.get_group('[1.0, 0.0, 0.0]')[\"flowRate\"].to_numpy()\n",
    "#poly = np.polynomial.Polynomial.fit(alpha, anisotropy_21, 1)\n",
    "\n",
    "anisotropy_21 = filter_group(anisotropy_21, qhigh = False)\n",
    "anisotropy_31 = filter_group(anisotropy_31, qhigh = False)\n",
    "\n",
    "fig, ax = plt.subplots(nrows = 1, ncols = 1)\n",
    "ax.plot(alpha, anisotropy_21, \"s\", label = \"$k_2$ / $k_1$\")\n",
    "ax.plot(alpha, anisotropy_31, \"^\", label = \"$k_3$ / $k_1$\")\n",
    "#ax.plot(alpha, poly(alpha), \"-\")\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.xlabel(r\"$\\alpha$\")\n",
    "plt.ylabel(\"Анизотропия проницаемости\")\n",
    "plt.title(\"Простая кубическая\")\n",
    "plt.show()\n",
    "\n",
    "if savefig:\n",
    "    fig.savefig(\"anisotropy-simple.png\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Body-centered structure"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEWCAYAAABrDZDcAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAAvcElEQVR4nO3deZwU1bn/8c9XQFFH0ICiBg0axcR44wLu3htwiwtqfgkhJlevGBVJ9BpvNDGoMcaoRINZjFHUGFGjQYPLRdziwhAVlwsobohBggpuMMoyGPbn90fVYNP0TJ+e6eqq7n7er1e/pqvqVNXT1dP1VJ1TdUpmhnPOufq1QdoBOOecS5cnAuecq3OeCJxzrs55InDOuTrnicA55+qcJ4IaIKlL2jFkjW8T58J5IqhCkvpIuk3SbEkfA79LO6YskDRM0vOS3gWaJG2fdkzOVQNPBCmQNEfSoTnDn4136lcGzLsF8DTwMrCbmW1hZt9PMNyqIOknwBnAKWa2rZl1M7O3047LuWrQOe0A6p2kLYHHgIfM7McBs/wAGG9mRZNGvZC0KfBjYHczeyfteJyrNn5GkCJJmwN/A54HzozHbS3pE0k9csrtJWl+XO+9D9AjPoP4SNJ4SdvmlDVJO8Xvt5f0L0l/joeHSnoqp2z+8BckPRovd6akITnTNpZ0laS3JC2S9FQ8rk+8zs5xue9LerUlfkknS5ohaUkc8+l52+Bnkt6T1BzHOqeN7XVsvOyFkholfTGe9CVgCXCNpI8lzZJ0WjzPv8fLXpobZ4Flt/XZx8Tz7pkz7sp43KHx8MUt2zke7hxP7xMPd5d0a/w9viXpQkkbFPoe4nFzJQ3IGR4k6cX4s0+W9OWcaflnmKdKaozf35/3+Zvj1+iceUdIei3edjdL6hpPGyBpbs5yh8TLODUn7qclXRP/T7wu6ZCc8tvG/58f5X4nOdtrZRzLQkn3StosnvZ5SU9IapK0QNLtin4rRT9vPLz2NxAPXyppTM7wXyW9H8f8d0lfypm2taS/xTE1xzFeTI3zRJCeBuAhorOy71rc14eZvQ80AkNyyp4IjDWzlcAmwMHAYGAb4C1gbCvr+AXQlDO8hla+c0VH1Y8CdwBbAccD10raNS4yCugHHAB8hugIfE3eMo4HzgW+amYt6/0QGAR0A04GfiNpr7j8F4ARwEAzawCOaeVzIKkv8BfgbGBL4EHgfkkbxttke2BevE0GA5dLOtjMnoyX/aWCCw777ACvAy07wC5xrB+0tswCfg90B3YEvgL8V7w9iooT0J+A04EewPXAeEkbFZvXzI7J+/ybm1mDmQ3PKfafwFeBzwN9gQsLxNCF6P/pvbxJ+wJvAj2BnwH3SPpMPG0sMBfYlpzvJGfeO+PYtgd2AE5qWR0wMp7vi8B2wMXFPmsJHgJ2JvqupwG350w7G1gNbBPHdmcZ15tZngjScx3QDPQGDsybdgtwAoCkTsC3gdtypv/JzKaZ2XKiHen+LUeeLeIjxv3jZbV4G/iipN4F4hkEzDGzm81slZm9ANwNfDM+cv0u8AMzm2dmq81scrz+FkcANwFHmtnao0gze8DM3rTIJKIzoH/PW3dIFeW3gAfM7NE4IY4CNiZKTBAlpfPMbJmZvQj8kWhnG6LVz55TZjxwmKSNiZLAY8CykIXH3+HxwAgzW2Jmc4CriBJ8iGHA9Wb2XLztbwGWA/sFzl/MNWb2jpl9BFxG9P+W73TgOeCNvPEfAr81s5VmdicwEzha0nZE/9ch30knon1RE4CZzYq/5+VmNh/4NVHyLAsz+1P8PSwnSjC7S+qeU2QD6mzfWFcfNmNeJ9qh/Bj4Y7yDafG/wK6SdgAOAxaZ2fPxtOVEZwEAmFkz0Q/os3nLvwL4KbAyZ9wkoiOc6ZIWAtfmTPscsG98Srwwnv6fwNZER3tdiY78WvNHYA55P1hJR0p6Nq4eWAgcFS8PM3ud6OjzSUnNwH1tLH/bvM+9Bngn/tzLgY/MbElO+bdYf5sALIirQGZIOiHgs7dYSZQMBgOnxJ8335Cc+RfkjO8JdMmNv434CvkccE5efNsRbZMW9+VMuzpwuS1y21XeylsucZXNj4n+n/LNazmbzZt/W4p/J0PieOcDS4H74/X1kjRW0jxJi4E/E//P5Cj2eaflTD8357N0kvRLSW/Gy54TT2pZ/lXAJ8CSeN7cM/Oa5YkgPZfFR0o3Ev0Qf9EywcyWAXcRnRWcyLpnA28T7RiAtdUaPYiqRVocHI+7K3eF8VH5cDPrYWabA7lXG70DTDKzzXNeDWb2PaKd2jKiqoPWfJvoqP2yljOOuOribqKj917xOh8kOvVvcRfRjnw34GttLP/dvM8top3hvHibfKaljjnWUlWUr6eZbQF8DxgjqaHIZ8/1R6IdYg8zm15g2Xe1zM+6O64FRInkcznjWouvkHeI/l9y49vEzP6SU+ZrOes+K3C5LbbLi+vdvOk/Ivpsb7G+z8bfRf7871L8O7krjncToqvgrorHXw4Y8G9m1o3od5C7Dij+effKmT4qZ/x3gOOAQ4mq6vrE4wUQn4E8SXTxxubk/YZqlSeCbDgNGCZpn5xxtwJDgWNZNxH8BThZ0h7xjvZy4Lm4uqHFxcCP847UipkA9JV0oqQu8WtvSV+Mj77/BPw6bgDsJGn/vDrqJ83sFaKjsxvicRsCGxEd8a2SdCRweN56rwN+lRd/IXcRVTkcEtdXn0OUQCab2TzgKWCkpK5xtdgpREeSrVlF9ONXW589d4b4DOZhovrrYGa2Oo7/MkmbSfoc8MMi8eW6ERguaV9FNpV0dN5OtiPOkNQ7rtu/gHXrxTcjasu4rJV5twLOirfZN4nq9B+Mr96aTNh3soZox79lzjqbgUWSPkuUiMplM6L/myaiBHR57sS4ivU81j1IqnmeCDLAzGYDFwE3x42fmNnTRD+QablHYmb2RFz2bqKGu88T1T/nesHMGkuMYQnRTvp4oqO594mql1p29ucSHbX9H/BRPK3Q/88vgW0knRQv8yyineDHREdj41sKSvoOUfXLbwPim0l0ZPh7oiPsY4BjzGxFXOQ/iRoc3wPuAS4ys8cKLGqhpCVEyXVYXFdc7LPnxvEjM/vfYvEW8N9E1R+ziZLWHUTJtcXeiq4UmqvoSp2tgb9K6m1mU4gOFq4h2o6ziA4SyuUOorab2UTVf5fmTOsGXG1mH7cy73NEDa8LiJLFYPv0QoFvEx1xvwvcC/ws7zv5Vlwl2ATsCpwfj/85sBewCHiA6Pssl1uJqqjmAa8Bz+ZNvx74ZStnPzVLpR00ukqS9ARwh5kVqo92NU7RJZFDA86WOrKOOcCprSTNYvMOjec9qNxxucryG8oyStLeREdFx6Udi0vNVOBfaQfhap8nggySdAtRw+kP8q66cHXEzM5JOwZXH7xqyDnn6pw3FjvnXJ2ruqqhnj17Wp8+fdIOY62lS5ey6aabph1Gq7IeH2Q/xqzHBx5jOWQ9PuhYjFOnTl1gZlsWnGhmVfXq16+fZcnEiRPTDqFNWY/PLPsxZj0+M48x39NvP21Pv/10SfPU+jYEplgr+9WqOyNwzrlizn88uiWhcWhjuoFUCW8jcM65OueJwDnn6pwnAuecq3OJthHEt68vIXrQwyoz6583fQBRl8v/jEfdY2aXJBmTc865dVWisXigmS1oY/qTZjYoqZX3v/RRFjSvWG98z4YNmXLhYUmt1jnnqkbNXzVUKAm0Nt6ThnOuHiXaxYSkfxJ1m2tEj9q7IW/6AKLulOcSdVV7rpm9WmA5w4ge10evXr36jR3b2iN61zf04aWtThtzxLo3ZpRStkVzczMNDQ3B8VRa1uOD7MdY6fjOemIpiwscv3TbEK4+uDr/D6GyMZ794tkA/HaP3wbPU+vbcODAgVPzq+dbJH1GcJCZzZO0FfCopNfN7O8506cBnzOzZklHET2qcOf8hcQJ5AaA/v3724ABA8IjePiBViett5zAsuueOYiom/n1zxyycIbR2Ni4/ufMmKzHWOn4Frfyf7h4RYH/2VjWtyFUNsbN52wOtL69CqnnbZhoIrDoyVGY2YeS7gX
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "df = load_data(execution, \"flowRate\")\n",
    "bodyCentered = df[df.label == \"bodyCentered\"].groupby(df.direction)\n",
    "alpha = bodyCentered.get_group(\"[0.0, 0.0, 1.0]\")[\"alpha\"].to_numpy()\n",
    "\n",
    "anisotropy_21 = 2 * bodyCentered.get_group('[0.0, 0.0, 1.0]')[\"flowRate\"].to_numpy() / bodyCentered.get_group('[1.0, 0.0, 0.0]')[\"flowRate\"].to_numpy()\n",
    "anisotropy_31 = 2 * bodyCentered.get_group('[1.0, 1.0, 1.0]')[\"flowRate\"].to_numpy() / bodyCentered.get_group('[1.0, 0.0, 0.0]')[\"flowRate\"].to_numpy()\n",
    "\n",
    "anisotropy_21 = filter_group(anisotropy_21)\n",
    "anisotropy_31 = filter_group(anisotropy_31)\n",
    "#poly = np.polynomial.Polynomial.fit(alpha, anisotropy_21, 10)\n",
    "\n",
    "fig, ax = plt.subplots(nrows = 1, ncols = 1)\n",
    "ax.plot(alpha, anisotropy_21, \"s\", label = r\"$k_2$ / $k_1$\")\n",
    "ax.plot(alpha, anisotropy_31, \"^\", label = r\"$k_3$ / $k_1$\")\n",
    "ax.axvline(0.134, linestyle = \"-.\", color = \"green\")\n",
    "#ax.plot(alpha, poly(alpha), \"-\")\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.xlabel(r\"$\\alpha$\")\n",
    "plt.ylabel(\"Анизотропия проницаемости\")\n",
    "plt.title(\"Кубическая объемноцентрированная\")\n",
    "\n",
    "if savefig:\n",
    "    fig.savefig(\"anisotropy-bodycentered.png\")"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Face-centered structure"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEWCAYAAABrDZDcAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjUuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/YYfK9AAAACXBIWXMAAAsTAAALEwEAmpwYAAApYklEQVR4nO3de5gcZZn+8e9NSIyQBDCRCARMUFCRVYFwEncNiBAUwdWICQsryMGssuhPQUQREBGEBRVWJIJyENGAUdi4gJHTRAU5JRCUhLgxIISDQIAwEyTH5/dH1YSiM4fqnq7q6Z77c119Tdf5eat76qn3rberFBGYmdnAtUGjAzAzs8ZyIjAzG+CcCMzMBjgnAjOzAc6JwMxsgHMiGCAkDW50DGbWPzkRtChJYyVdJWmxpBeACxodk5n1T04E/ZSkRyXtmxneKj2on5tj2c2AO4A/ATtGxGYR8dkCwzWzJrZhowOw3kl6I3ALcFNEfDnHIp8HZkZEr0nDzMw1gn5O0qbAb4F7gOPScW+S9LKkkZn5dpb0bHotYDdgZFqDeF7STElbZuYNSW9N328j6R+SfpoOHyHpD5l5K4ffLunmdL0LJR2SmfZ6SedL+pukZZL+kI4bm25zw3S+z0p6qDN+SUdKWiCpPY35MxX74DRJT0nqSGN9tJt91bmdjszrnHRam6SzJd0j6SVJ/yPpDZllfyHp6TTu30l6Z2baFZLOzAy/VVJkhjeR9OM0xicknSlpUFf7Lx23RNKEzPCBkh6Q9KKkOyW9KzOtsmZ4tKS29P2v0zIuryj3tMyyJ0uaL+kFSZdLGppOmyBpSWa9h6TrODoT9x2Svp/uk4clfSAz/5bp9+p5SYskHZOZdrqkVWksL0q6TtLwdNpbJN0maamk5yRdreQ73mt50+F13910+ExJV+T8HN8k6bdpTB1pjKdjTgT93DDgJpKa26cjvR9IRDwNtAGHZOY9HJgeEauAjYB9gEnAFsDfgOndbOObwNLM8Fq6+V5I2hi4GfgZsDkwGfiBpB3SWc4DdgHeC7wB+HK6vuw6JgMnAPtHROd2nwEOBEYARwLflbRzOv/bgZOBvSNiGPCRbsqRtWlEDEtfJ2XG/zvwaZJ9shq4MDPtJmC7tFxzgatzbKfTFen63grsBOwHHJ1nQUk7AZcBnwFGAj8EZkp6XW/LRsRH0n3SebDrLPfUzGz/BuwPvAXYHjilixgGk3wPnqqYtDvwV2AUcBrwq0zynA4sAbYk+Z6dJWmfzLLXpLFtA4wDPtW5OeDsdLl3AFsDp/dW1ir09Dl+AVgDbJHGdk0dt9vUnAj6t4uBDmAMsFfFtCuBwwDSs88pwFWZ6ZdFxNyIWEFyIN1T0tjsCtIzzz3TdXV6DHiHpDFdxHMg8GhEXB4RqyPifuCXwCckbUBykP18RDwREWsi4s50+50mAj8GDoiIdWejEXFDRPw1ErNJakD/XLHtejRjXhURf46I5cDXgUM6z9wj4rKIaE/jPR14t6RNeluhpNHAh4AvRMTyiHgG+C5JkszjWOCHEXF3us+uBFYAe1Rduq59PyIej4jngW+RfE8qfQa4G/hLxfhngO9FxKqIuAZYCHxY0tYk38eTIuKViHgA+BFJoq00iOQ4sxQgIhZFxM0RsSIingW+A7y/z6VM5fgcN8DHvfV4h/RvD5OcAX8Z+JGk12em/Q+wg6RxwAeBZRFxTzptBUktAICI6CD5R9yqYv3nkBwQV2XGzSY5U5on6UXgB5lpbwZ2T6vWL6bT/w14E8lZ41CSM8ju/Ah4lIp/fEkHSLorbWZ4keTAOiqN/WGSs9jfS+oAru9h/b15PPP+b8BgYJSkQZK+Lemvkl5KY6QzhtQJmTLPzYx/c7qepzLTf0hyRprHm4EvVezTrUnOmDtdn5l2YRfr6EllmbPrJW2y+TLJ96DSE5210IrltwSej4j2imnZ79chabzPAsuBX6fbGy1petqE9hLwU167n6H38s7NTD8hU5bePsfzgZeB9nTZbI16QHMi6N++lZ5xXUryD/3NzgkR8QpwLUmt4HBeWxt4jOQAA6xr0hkJPJGZZ5903LXZDaZn5VMjYmREbApkexs9DsyOiE0zr2ER8R/Ac8ArJE0Q3ZkCfBL4VmeNI20C+SVJs9LodJs3kjQhdLqWJLntCHy0h/X3ZuvM+21IEuBzwKHAwcC+wCbA2HSebAzndZYZ2Dkz/vE0tlGZfTIiIt5JPo+TfM7ZfbpRRPw8M89HM9s+Pud6O1WW+cmK6ScC10bE31jfVpKy+6Bz+SeBN3S2+2emZb9f16bxbkTSe+38dPxZQAD/FBEjSL6/2W1A7+XdOTP9vMz4Hj/HtAbye5JOF5tS8d0fyJwImscxwLGSdsuM+wlwBHAQr00EPweOlPSe9EB7FnB3RDyamed04MsVZ3y9+V9ge0mHSxqcvnaV9I6IWEvS1v2d9ELiIEl7VrR1/z4i/kxylndJOm4I8DqSM8fVkg4gaWPPuhj4r4r4a3GYpB0kbQScAcyIiDXAcJKD+VKSA9dZeVcYEU+RNGWdL2mEpA3SC6J5mzsuBaZK2l2JjSV9uOIg2xefkzQmbdv/Gq9tFx9Ock3mW90suzlwfPo5f4KkTf/GiHgcuBM4W9LQtInxKJKz+0prSQ78b8xsswNYJmkrkkRULz1+jmnT6Em89uTGcCJoGhGxGDgVuFzSkHTcHST/aHOzZ3QRcVs67y9JLgC+hfXbrO+PiLYqY2gnOUhPJjkrfJqkeanzYH8CydnfvcDz6bSuvmPfBraQ9Kl0nceTnJ29QHJWN7NzRkmHkjQ9fa+aWLtxFcmF3adJmrE6zzZ/QtK08QQwH7iryvX+O0lCm09ShhkkF6Q77aqkp9ASJT113gT8QtKYiLiPJMl/P112EUlyr5efkSSqxSTNdmdmpo0ALoyIF7pZ9m6SC6/PkSSLSZkL/FNIzrifBK4DTouIWzLLfjJtylsK7AB8NR3/DZIa1TLgBuBXfSlchd4+xx8C3+6m9jOgqboTQutvJN0G/CwiftToWPozJV0Qf9pf9lMazxF1qOX0tI1HgaMrDtB5lz0iXfZ99Y7L+h//oKyJSdqV5Ozq4EbHYlWbA/yj0UGYgRNB05J0JcmF089X9N6wJhARX2p0DGad3DRkZjbA+WKxmdkA13RNQ6NGjYqxY8c2Ooxcli9fzsYbb9zoMArRymWD1i6fy9a8+lK+OXPmPBcRb+xqWtMlgrFjx3Lfffc1Ooxc2tramDBhQqPDKEQrlw1au3wuW/PqS/kkddtt1k1DZmYDnBOBmdkA50RgZjbAORGYmQ1wTgRmZgNc0/UaMjMbSMafeTPPdax8dcRvbgBg1LAh3HfKB+uyDdcIzMz6sdckgRzja+FEYGY2wLlpyMysBOs18aTq2cRTq0JrBJImSlooaZGkr3Qx/c2SbpX0oKS2bh6YbmbW9Mpo4qlVYTUCSYOAi0gerL4EuFfSzIiYn5ntPOAnEXGlpH2As0mev2tm1i/15zP7WhXZNLQbsCh9xCKSppM8QCWbCHYAvpi+vx24vsB4zKwktRwsyz7A1rq9ss/sRw0b0m2c9VJkItgKeDwzvATYvWKeecDHgAuAfwWGSxqZeS4qAJKOBY4FGD16NG1tbUXFXFcdHR2Fx3r8bct5qYvv34ghcOE+xd2FsYyyNVIrl6+astX6/erpYNndtmtZpssY0+6VRcTYmyKWO+99g4HBQPLZDRs2rM/bq9Toi8UnAN9Pn4/6O5KHTq+pnCkiLgEuARg/fnw0y90Fq7lTYK1nJy+lX/r1xq+k223X48zLd3lsLq/9zAUsB4r5fgHrDsZd6Xa5Wpah5BgbsVxGUd/LIhPBE8DWmeEx6bh1IuJJkhoBkoYBH4+IFwuMqXC1/vijzOpmrdsq44ctVoz+fKFyoCijiadWRSaCe4HtJI0jSQCTgUOzM0gaBTwfEWuBk4HLCoynFK38D1d22VrxolxfeZ80r/78+RSWCCJitaTjgFnAIOCyiHhI0hnAfRExE5gAnC0pSJqGPldUPNZ86lZzSfV2sGyGGk8rn2g0i/58Zl+rQq8RRMS
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "df = load_data(execution, \"flowRate\")\n",
    "faceCentered = df[df.label == \"faceCentered\"].groupby(df.direction)\n",
    "alpha = faceCentered.get_group(\"[0.0, 0.0, 1.0]\")[\"alpha\"].to_numpy()\n",
    "\n",
    "anisotropy_21 = 2 * faceCentered.get_group('[0.0, 0.0, 1.0]')[\"flowRate\"].to_numpy() / faceCentered.get_group('[1.0, 0.0, 0.0]')[\"flowRate\"].to_numpy()\n",
    "anisotropy_31 = 2 * faceCentered.get_group('[1.0, 1.0, 1.0]')[\"flowRate\"].to_numpy() / faceCentered.get_group('[1.0, 0.0, 0.0]')[\"flowRate\"].to_numpy()\n",
    "\n",
    "anisotropy_21 = filter_group(anisotropy_21)\n",
    "anisotropy_31 = filter_group(anisotropy_31)\n",
    "#poly = np.polynomial.Polynomial.fit(alpha, anisotropy_21, 10)\n",
    "\n",
    "fig, ax = plt.subplots(nrows = 1, ncols = 1)\n",
    "ax.plot(alpha, anisotropy_21, \"s\", label = r\"$k_2$ / $k_1$\")\n",
    "ax.plot(alpha, anisotropy_31, \"^\", label = r\"$k_3$ / $k_1$\")\n",
    "#ax.plot(alpha, poly(alpha), \"-\")\n",
    "plt.legend()\n",
    "plt.grid(True)\n",
    "plt.xlabel(r\"$\\alpha$\")\n",
    "plt.ylabel(\"Анизотропия проницаемости\")\n",
    "plt.title(\"Кубическая гранецентрированная\")\n",
    "plt.show()\n",
    "\n",
    "if savefig:\n",
    "    fig.savefig(\"anisotropy-facecentered.png\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "19649669bd52b0be75e091dcf60d2128e4a347083ff474cfec5ff9275df3ceed"
  },
  "kernelspec": {
   "display_name": "Python 3.9.9 64-bit ('anisotropy': conda)",
   "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.10"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}