Skip to content

Instantly share code, notes, and snippets.

@pshriwise
Last active August 10, 2023 20:36
Show Gist options
  • Save pshriwise/85f670f9e29f677d4a863c91093bf726 to your computer and use it in GitHub Desktop.
Save pshriwise/85f670f9e29f677d4a863c91093bf726 to your computer and use it in GitHub Desktop.
SurfaceFilter and MuFitler combination
Display the source blob
Display the rendered blob
Raw
{
"cells": [
{
"cell_type": "code",
"execution_count": 1,
"metadata": {},
"outputs": [],
"source": [
"%matplotlib inline\n",
"import numpy as np\n",
"import matplotlib.pyplot as plt\n",
"import openmc\n",
"import os"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {},
"outputs": [],
"source": [
"# openmc_data = \"/home/zoeprieto/Doctorado/openmc_ncrystal_conda/endfb-viii.0-hdf5/cross_sections.xml\"\n",
"openmc_data = '/Users/pshriwise/data/xs/openmc/nndc_hdf5/cross_sections.xml'\n",
"os.environ[\"OPENMC_CROSS_SECTIONS\"] = openmc_data\n",
"openmc.config['cross_sections'] = openmc_data"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {},
"outputs": [],
"source": [
"mat=openmc.Material(material_id=1) \n",
"mat.add_nuclide(\"H1\",2)\n",
"mat.add_nuclide(\"O16\",1)\n",
"mat.add_s_alpha_beta('c_H_in_H2O')\n",
"mat.set_density('g/cm3', 1)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"<matplotlib.image.AxesImage at 0x13bd955d0>"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
},
{
"data": {
"image/png": "iVBORw0KGgoAAAANSUhEUgAAARoAAAEMCAYAAAAbPHk8AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjcuMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8pXeV/AAAACXBIWXMAAA9hAAAPYQGoP6dpAAAcLUlEQVR4nO3de1BU5/0G8GdBWfDC4gVZ1iwXxfsFvIEQEzDQIioJTmoZjEK8pVpJNTu2jYmKkTQ0URtqpaiTKJ14iZpEdBzDFDHRWFELiNVMdYCCqLiIRnaB1IWw5/dHfm6zwiKrvOwuPJ+ZM+N5z/ue8z3M7uO57O6RSZIkgYhIICdbF0BEXR+DhoiEY9AQkXAMGiISjkFDRMIxaIhIOAYNEQnHoCEi4Rg0RCQcg4aIhLOboDl9+jRiY2OhUqkgk8mQnZ1tWtbU1ITf//73GDduHHr37g2VSoXExERUVVW1uc4NGzZAJpOZTSNHjhS8J0T0KLsJmoaGBgQGBiIjI6PFsu+//x5FRUVYt24dioqK8MUXX+DatWt48cUXH7veMWPG4Pbt26bpzJkzIsonojb0sHUBD8XExCAmJqbVZQqFArm5uWZt27ZtQ3BwMCorK+Hj42NxvT169IBSqezQWonIOnYTNNbS6XSQyWTw8PBos19JSQlUKhVcXV0RGhqKtLS0NoPJYDDAYDCY5o1GI7777jsMGDAAMpmso8on6hIkSUJdXR1UKhWcnNo4QZLsEADp8OHDFpf/97//lSZOnCjNmzevzfUcP35cOnjwoHTp0iUpJydHCg0NlXx8fCS9Xm9xTEpKigSAEydOVkw3btxo870o+/83tl2RyWQ4fPgw4uLiWixramrCyy+/jJs3b+Lrr7+Gu7t7u9dbW1sLX19f/OlPf8LixYtb7fPoEY1Op4OPjw9u3Lhh1baIugO9Xg+1Wo3a2looFAqL/Rzq1KmpqQm//OUvcf36dZw8edLqN76HhweGDx+O0tJSi33kcjnkcnmLdnd3dwYNkQWPu6xgN3edHudhyJSUlODEiRMYMGCA1euor69HWVkZvL29BVRIRJbYTdDU19ejuLgYxcXFAIDy8nIUFxejsrISTU1N+MUvfoGCggLs3bsXzc3N0Gq10Gq1aGxsNK0jMjIS27ZtM82vXr0ap06dQkVFBc6ePYs5c+bA2dkZCQkJnb17RN2a3Zw6FRQUYPr06aZ5jUYDAEhKSsKGDRtw9OhRAEBQUJDZuK+++goREREAgLKyMty9e9e07ObNm0hISMC9e/fg6emJadOm4dy5c/D09BS7M0Rkxi4vBtsTvV4PhUIBnU7HazREj2jv+8NuTp2IqOti0BCRcAwaIhKOQUNEwjFoiEg4Bg0RCcegISLhGDREJByDhoiEY9AQkXAMGiISjkFDRMIxaIhIOAYNEQnHoCEi4Rg0RCQcg4aIhGPQEJFwDBoiEo5BQ0TCMWiISDgGDREJZzdBc/r0acTGxkKlUkEmkyE7O9tsuSRJWL9+Pby9veHm5oaoqCiUlJQ8dr0ZGRnw8/ODq6srQkJCcOHCBUF7QESW2E3QNDQ0IDAwEBkZGa0u/+CDD7B161Zs374d58+fR+/evREdHY0HDx5YXOeBAweg0WiQkpKCoqIiBAYGIjo6Gnfu3BG1G0TUGskOAZAOHz5smjcajZJSqZQ2bdpkaqutrZXkcrm0f/9+i+sJDg6WVqxYYZpvbm6WVCqVlJaW1u5adDqdBEDS6XTW7QRRN9De94fdHNG0pby8HFqtFlFRUaY2hUKBkJAQ5OfntzqmsbERhYWFZmOcnJwQFRVlcQwAGAwG6PV6s4mIno7dPHu7LVqtFgDg5eVl1u7l5WVa9qi7d++iubm51TFXr161uK20tDS88847T1lxxzt+aq6tSyA7NDP8kK1LaBeHOKLpTGvWrIFOpzNNN27csHVJRA7PIYJGqVQCAKqrq83aq6urTcseNXDgQDg7O1s1BgDkcjnc3d3NJiJ6Og4RNP7+/lAqlcjLyzO16fV6nD9/HqGhoa2OcXFxwaRJk8zGGI1G5OXlWRxDRGLYzTWa+vp6lJaWmubLy8tRXFyM/v37w8fHB6tWrcK7776LYcOGwd/fH+vWrYNKpUJcXJxpTGRkJObMmYPk5GQAgEajQVJSEiZPnozg4GCkp6ejoaEBCxcu7OzdI+rW7CZoCgoKMH36dNO8RqMBACQlJSErKwu/+93v0NDQgNdeew21tbWYNm0acnJy4OrqahpTVlaGu3fvmubj4+NRU1OD9evXQ6vVIigoCDk5OS0uEBORWDJJkiRbF2HP9Ho9FAoFdDqdTa/X8K4TtcbWd53a+/5wiGs0ROTYGDREJByDhoiEY9AQkXAMGiISjkFDRMIxaIhIOAYNEQnHoCEi4Rg0RCQcg4aIhGPQEJFwDBoiEo5BQ0TCMWiISDgGDREJx6AhIuEYNEQkHIOGiIRj0BCRcAwaIhKOQUNEwjFoiEg4hwkaPz8/yGSyFtOKFSta7Z+VldWi708fNkdEncdunlT5OP/85z/R3Nxsmr9y5Qp+9rOfYe5cyw9Wc3d3x7Vr10zzMplMaI1E1DqHCRpPT0+z+T/+8Y8YOnQowsPDLY6RyWRQKpVWbcdgMMBgMJjm9Xq9dYUSUQsOc+r0U42NjdizZw8WLVrU5lFKfX09fH19oVar8dJLL+Hbb7997LrT0tKgUChMk1qt7sjSibolhwya7Oxs1NbW4tVXX7XYZ8SIEdi1axeOHDmCPXv2wGg0IiwsDDdv3mxz3WvWrIFOpzNNN27c6ODqibofhzl1+qmPP/4YMTExUKlUFvuEhoYiNDTUNB8WFoZRo0Zhx44dSE1NtThOLpdDLpd3aL1E3Z3DBc3169dx4sQJfPHFF1aN69mzJyZMmIDS0lJBlRGRJQ536rR7924MGjQIs2bNsmpcc3MzLl++DG9vb0GVEZElDhU0RqMRu3fvRlJSEnr0MD8YS0xMxJo1a0zzGzduxN///nf85z//QVFREebPn4/r169jyZIlnV02UbfnUKdOJ06cQGVlJRYtWtRiWWVlJZyc/peb9+/fx9KlS6HVatGvXz9MmjQJZ8+exejRozuzZCICIJMkSbJ1EfZMr9dDoVBAp9PB3d3dZnUcP2X5g4nUfc0MP2TT7bf3/eFQp05E5JgYNEQkHIOGiIRj0BCRcAwaIhKOQUNEwrXrczQajcbqFa9duxb9+/e3ehwRdT3tCpr09HSEhobCxcWlXSs9c+YMkpOTGTREBMCKTwYfPnwYgwYNalffvn37PnFBRNT1tOsaze7du6FQKNq90h07dsDLy+uJiyKirqVdRzRJSUlWrXTevHlPVAwRdU1P9aXK+vp6GI1GszZbfh+IiOyT1be3y8vLMWvWLPTu3RsKhQL9+vVDv3794OHhgX79+omokYgcnNVHNPPnz4ckSdi1axe8vLz4CBMieiyrg+bSpUsoLCzEiBEjRNRDRF2Q1adOU6ZM4ZMBiMgqVh/RfPTRR1i2bBlu3bqFsWPHomfPnmbLx48f32HFEVHXYHXQ1NTUoKysDAsXLjS1yWQySJIEmUxm9thaIiLgCYJm0aJFmDBhAvbv38+LwUTULlYHzfXr13H06FEEBASIqIeIuiCrLwa/8MILuHTpkohaiKiLsvqIJjY2Fm+88QYuX76McePGtbgY/OKLL3ZYcUTUNVgdNMuWLQPw4wPaHiXyYvCGDRvwzjvvmLWNGDECV69etTjm0KFDWLduHSoqKjBs2DC8//77mDlzppD6iMgyq0+djEajxUn0HacxY8bg9u3bpunMmTMW+549exYJCQlYvHgxLl68iLi4OMTFxeHKlStCaySilhzqpzx79OgBpVJpmgYOHGix75///GfMmDEDv/3tbzFq1CikpqZi4sSJ2LZtWydWTETAEwTNb37zG2zdurVF+7Zt27Bq1aqOqMmikpISqFQqDBkyBK+88goqKyst9s3Pz0dUVJRZW3R0NPLz89vchsFggF6vN5uI6OlYHTSff/45nn322RbtYWFh+OyzzzqkqNaEhIQgKysLOTk5yMzMRHl5OZ577jnU1dW12l+r1bb48S0vLy9otdo2t5OWlgaFQmGa1Gp1h+0DUXdlddDcu3ev1V/bc3d3x927dzukqNbExMRg7ty5GD9+PKKjo3H8+HHU1tbi4MGDHbqdNWvWQKfTmSZ+r4vo6VkdNAEBAcjJyWnR/uWXX2LIkCEdUlR7eHh4YPjw4SgtLW11uVKpRHV1tVlbdXU1lEplm+uVy+Vwd3c3m4jo6Vh9e1uj0SA5ORk1NTV44YUXAAB5eXnYsmUL0tPTO7o+i+rr61FWVoYFCxa0ujw0NBR5eXlm141yc3MRGhraSRUS0UNP9F0ng8GAP/zhD0hNTQUA+Pn5ITMzE4mJiR1e4EOrV69GbGwsfH19UVVVhZSUFDg7OyMhIQEAkJiYiMGDByMtLQ0AsHLlSoSHh2PLli2YNWsWPv30UxQUFGDnzp3CaiSi1j3RbwYvX74cy5cvR01NDdzc3NCnT5+OrquFmzdvIiEhAffu3YOnpyemTZuGc+fOwdPTEwBQWVkJJ6f/nQmGhYVh3759WLt2Ld566y0MGzYM2dnZGDt2rPBaicicTJIkydZF2DO9Xg+FQgGdTmfT6zXHT8212bbJfs0MP2TT7bf3/dGui8ETJ07E/fv3273xadOm4datW+3uT0RdW7tOnYqLi3Hp0qV2P+K2uLgYBoPhqQojoq6j3ddoIiMj0d6zLP4YFhH9VLuCpry83OoVP/PMM1aPIaKuqV1B4+vrK7oOIurCHOrb20TkmBg0RCQcg4aIhGPQEJFwVgdNUlISTp8+LaIWIuqirA4anU6HqKgoDBs2DO+99x4/AUxEj2V10GRnZ+PWrVtYvnw5Dhw4AD8/P8TExOCzzz5DU1OTiBqJyME90TUaT09PaDQaXLp0CefPn0dAQAAWLFgAlUqFN954AyUlJR1dJxE5sKe6GHz79m3k5uYiNzcXzs7OmDlzJi5fvozRo0fjww8/7KgaicjBWR00TU1N+PzzzzF79mz4+vri0KFDWLVqFaqqqvC3v/0NJ06cwMGDB1t9wBwRdU9W//CVt7c3jEYjEhIScOHCBQQFBbXoM336dHh4eHRAeUTUFVgdNB9++CHmzp0LV1dXi308PDye6IuYRNQ1WR00ln4MnIjIEn4ymIiEY9AQkXAMGiISjkFDRMI5TNCkpaVhypQp6Nu3LwYNGoS4uDhcu3atzTFZWVmQyWRmU1t3y4hIDIcJmlOnTmHFihU4d+4ccnNz0dTUhJ///OdoaGhoc5y7uztu375tmq5fv95JFRPRQ0/0pEpbyMnJMZvPysrCoEGDUFhYiOeff97iOJlMBqVSKbo8ImqDwxzRPEqn0wHAY581VV9fD19fX6jVarz00kv49ttv2+xvMBig1+vNJiJ6Og4ZNEajEatWrcKzzz7b5rO0R4wYgV27duHIkSPYs2cPjEYjwsLCcPPmTYtj0tLSoFAoTJNarRaxC0TdikM+e3v58uX48ssvcebMGaueH9XU1IRRo0YhISEBqamprfYxGAxmT9nU6/VQq9V89jbZJUd59rbDXKN5KDk5GceOHcPp06etfkhdz549MWHCBJSWllrsI5fLIZfLn7ZMIvoJhzl1kiQJycnJOHz4ME6ePAl/f3+r19Hc3IzLly/D29tbQIVEZInDHNGsWLEC+/btw5EjR9C3b19otVoAgEKhgJubGwAgMTERgwcPRlpaGgBg48aNmDp1KgICAlBbW4tNmzbh+vXrWLJkic32g6g7cpigyczMBABERESYte/evRuvvvoqAKCyshJOTv87SLt//z6WLl0KrVaLfv36YdKkSTh79ixGjx7dWWUTERz0YnBnau/FLtF4MZha4ygXgx3mGg0ROS4GDREJx6AhIuEYNEQkHIOGiIRj0BCRcAwaIhKOQUNEwjFoiEg4Bg0RCcegISLhGDREJByDhoiEY9AQkXAMGiISjkFDRMIxaIhIOAYNEQnHoCEi4Rg0RCQcg4aIhGPQEJFwDBoiEs7hgiYjIwN+fn5wdXVFSEgILly40Gb/Q4cOYeTIkXB1dcW4ceNw/PjxTqqUiB5yqKA5cOAANBoNUlJSUFRUhMDAQERHR+POnTut9j979iwSEhKwePFiXLx4EXFxcYiLi8OVK1c6uXKi7s2hnlQZEhKCKVOmYNu2bQAAo9EItVqN119/HW+++WaL/vHx8WhoaMCxY8dMbVOnTkVQUBC2b9/e6jYMBgMMBoNpXq/XQ61W80mVZJf4pMoO1tjYiMLCQkRFRZnanJycEBUVhfz8/FbH5Ofnm/UHgOjoaIv9ASAtLQ0KhcI0qdXqjtkBom6sh60LaK+7d++iubkZXl5eZu1eXl64evVqq2O0Wm2r/bVarcXtrFmzBhqNxjT/8IjG1mz9PxfR03CYoOkscrkccrnc1mUQdSkOc+o0cOBAODs7o7q62qy9uroaSqWy1TFKpdKq/kQkhsMEjYuLCyZNmoS8vDxTm9FoRF5eHkJDQ1sdExoaatYfAHJzcy32JyIxHOrUSaPRICkpCZMnT0ZwcDDS09PR0NCAhQsXAgASExMxePBgpKWlAQBWrlyJ8PBwbNmyBbNmzcKnn36KgoIC7Ny505a7QdTtOFTQxMfHo6amBuvXr4dWq0VQUBBycnJMF3wrKyvh5PS/g7SwsDDs27cPa9euxVtvvYVhw4YhOzsbY8eOtdUuEHVLDvU5Glto7+cEiLqjLvc5GiJyXAwaIhKOQUNEwjFoiEg4Bg0RCcegISLhGDREJByDhoiEY9AQkXAMGiISjkFDRMIxaIhIOAYNEQnHoCEi4Rg0RCQcg4aIhGPQEJFwDBoiEo5BQ0TCMWiISDgGDREJx6AhIuEcImgqKiqwePFi+Pv7w83NDUOHDkVKSgoaGxvbHBcREQGZTGY2LVu2rJOqJqKHHOIBclevXoXRaMSOHTsQEBCAK1euYOnSpWhoaMDmzZvbHLt06VJs3LjRNN+rVy/R5RLRIxwiaGbMmIEZM2aY5ocMGYJr164hMzPzsUHTq1cvKJVK0SUSURscImhao9Pp0L9//8f227t3L/bs2QOlUonY2FisW7euzaMag8EAg8Fgth3gxyfyEZG5h++Lxz7wVnJAJSUlkru7u7Rz5842++3YsUPKycmR/vWvf0l79uyRBg8eLM2ZM6fNMSkpKRIATpw4WTHduHGjzfeVTZ+9/eabb+L9999vs8+///1vjBw50jR/69YthIeHIyIiAh999JFV2zt58iQiIyNRWlqKoUOHttrn0SMao9GI7777DgMGDIBMJrNqex1Fr9dDrVbjxo0b3f753/xb/Mhe/g6SJKGurg4qlQpOTpbvLdk0aGpqanDv3r02+wwZMgQuLi4AgKqqKkRERGDq1KnIyspqc8da09DQgD59+iAnJwfR0dFPXHdna++D1LsD/i1+5Gh/B5teo/H09ISnp2e7+t66dQvTp0/HpEmTsHv3bqtDBgCKi4sBAN7e3laPJaIn5xCfo7l16xYiIiLg4+ODzZs3o6amBlqtFlqt1qzPyJEjceHCBQBAWVkZUlNTUVhYiIqKChw9ehSJiYl4/vnnMX78eFvtClG35BB3nXJzc1FaWorS0lI888wzZssenvk1NTXh2rVr+P777wEALi4uOHHiBNLT09HQ0AC1Wo2XX34Za9eu7fT6n5ZcLkdKSgrkcrmtS7E5/i1+5Gh/B5teoyGi7sEhTp2IyLExaIhIOAYNEQnHoCEi4Rg0DiAjIwN+fn5wdXVFSEiI6RZ+d7Fhw4YWP/fx00+Ld2WnT59GbGwsVCoVZDIZsrOzzZZLkoT169fD29sbbm5uiIqKQklJiW2KbQODxs4dOHAAGo0GKSkpKCoqQmBgIKKjo3Hnzh1bl9apxowZg9u3b5umM2fO2LqkTtHQ0IDAwEBkZGS0uvyDDz7A1q1bsX37dpw/fx69e/dGdHQ0Hjx40MmVPkZ7v8hIthEcHCytWLHCNN/c3CypVCopLS3NhlV1rpSUFCkwMNDWZdgcAOnw4cOmeaPRKCmVSmnTpk2mttraWkkul0v79++3QYWW8YjGjjU2NqKwsBBRUVGmNicnJ0RFRSE/P9+GlXW+kpISqFQqDBkyBK+88goqKyttXZLNlZeXQ6vVmr0+FAoFQkJC7O71waCxY3fv3kVzczO8vLzM2r28vMy+ftHVhYSEICsrCzk5OcjMzER5eTmee+451NXV2bo0m3r4GnCE14dDfAWBureYmBjTv8ePH4+QkBD4+vri4MGDWLx4sQ0ro/biEY0dGzhwIJydnVFdXW3WXl1d3a1/ntTDwwPDhw9HaWmprUuxqYevAUd4fTBo7JiLiwsmTZqEvLw8U5vRaEReXh5CQ0NtWJlt1dfXo6ysrNv/3Ie/vz+USqXZ60Ov1+P8+fN29/rgqZOd02g0SEpKwuTJkxEcHGz6NvrChQttXVqnWb16NWJjY+Hr64uqqiqkpKTA2dkZCQkJti5NuPr6erMjt/LychQXF6N///7w8fHBqlWr8O6772LYsGHw9/fHunXroFKpEBcXZ7uiW2Pr2170eH/5y18kHx8fycXFRQoODpbOnTtn65I6VXx8vOTt7S25uLhIgwcPluLj46XS0lJbl9Upvvrqq1Z/ozcpKUmSpB9vca9bt07y8vKS5HK5FBkZKV27ds22RbeCPxNBRMLxGg0RCcegISLhGDREJByDhoiEY9AQkXAMGiISjkFDRMIxaIhIOAYN2Y2KigrTT3UGBQUJ3VZWVpZpW6tWrRK6LWLQkB06ceKE2RcFRYiPj8ft27ft7suHXRW/VEl2Z8CAARgwYIDQbbi5ucHNzQ0uLi5Ct0M/4hENCVFTUwOlUon33nvP1Hb27Fm4uLg80dHKrl27MGbMGMjlcnh7eyM5Odm0TCaTYceOHZg9ezZ69eqFUaNGIT8/H6WlpYiIiEDv3r0RFhaGsrKyDtk3sh6DhoTw9PTErl27sGHDBhQUFKCurg4LFixAcnIyIiMjrVpXZmYmVqxYgddeew2XL1/G0aNHERAQYNYnNTUViYmJKC4uxsiRIzFv3jz86le/wpo1a1BQUABJkszCiTqZjb89Tl3cr3/9a2n48OHSvHnzpHHjxkkPHjyw2Le8vFwCIF28eNGsXaVSSW+//bbFcQCktWvXmubz8/MlANLHH39satu/f7/k6uraYmx4eLi0cuXK9u8QPREe0ZBQmzdvxg8//IBDhw5h7969kMvlVo2/c+cOqqqqHnsUNH78eNO/H/5Y97hx48zaHjx4AL1eb9X2qWMwaEiosrIyVFVVwWg0oqKiwurxbm5u7erXs2dP079lMpnFNqPRaHUN9PQYNCRMY2Mj5s+fj/j4eKSmpmLJkiVWP2Gzb9++8PPzE367m8Ti7W0S5u2334ZOp8PWrVvRp08fHD9+HIsWLcKxY8esWs+GDRuwbNkyDBo0CDExMairq8M//vEPvP7664Iqp47GIxoS4uuvv0Z6ejo++eQTuLu7w8nJCZ988gm++eYbZGZmWrWupKQkpKen469//SvGjBmD2bNn2+WD7Mky/mYw2Y2Kigr4+/vj4sWLwr+C8FBERASCgoKQnp7eKdvrrnhEQ3YnLCwMYWFhQrexd+9e9OnTB998843Q7dCPeERDduOHH34w3ZmSy+VQq9XCtlVXV2d6wqOHhwcGDhwobFvEoCGiTsBTJyISjkFDRMIxaIhIOAYNEQnHoCEi4Rg0RCQcg4aIhGPQEJFw/wf3w2Ev40RBSAAAAABJRU5ErkJggg==",
"text/plain": [
"<Figure size 258.065x259.74 with 1 Axes>"
]
},
"metadata": {},
"output_type": "display_data"
}
],
"source": [
"Lx,Ly,Lz=10.0,10.0,10.0\n",
"\n",
"surf01=openmc.XPlane(boundary_type = \"vacuum\", surface_id = 1)\n",
"surf02=openmc.YPlane(boundary_type = \"vacuum\", surface_id = 2)\n",
"surf03=openmc.ZPlane(boundary_type = \"vacuum\", surface_id = 3) \n",
"surf04=openmc.XPlane(Lx,boundary_type = \"vacuum\", surface_id = 4)\n",
"surf05=openmc.YPlane(Ly,boundary_type = \"vacuum\", surface_id = 5)\n",
"surf06=openmc.ZPlane(Lz,boundary_type = \"vacuum\", surface_id = 6)\n",
"\n",
"cell=openmc.Cell(region = +surf01 & -surf04 & +surf02 & -surf05 & +surf03 & -surf06, fill=mat , cell_id = 1 )\n",
"\n",
"univ=openmc.Universe(cells = [cell], universe_id = 1)\n",
"\n",
"univ.plot(origin=(Lx/2, Ly/2, Lz/2), width=(Lx*1.5, Lz*1.5), pixels=(200, 200), basis='xy', color_by='cell')"
]
},
{
"cell_type": "code",
"execution_count": 5,
"metadata": {},
"outputs": [],
"source": [
"sets = openmc.Settings()\n",
"sets.run_mode = \"fixed source\"\n",
"sets.particles = 10000\n",
"sets.batches = 110\n",
"sets.inactive = 10\n",
"sets.photon_transport = True\n",
"\n",
"sets.surf_source_write = {\n",
" 'surface_ids': [6],\n",
" 'max_particles': 10000\n",
"}"
]
},
{
"cell_type": "code",
"execution_count": 6,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/Users/pshriwise/soft/openmc/openmc/source.py:388: FutureWarning: This class is deprecated in favor of 'IndependentSource'\n",
" warnings.warn(\"This class is deprecated in favor of 'IndependentSource'\", FutureWarning)\n"
]
}
],
"source": [
"E_f=1.0e6\n",
"\n",
"X=np.linspace(0,Lx,50)\n",
"Y=np.linspace(0,Ly,50)\n",
"PX=np.sin(np.pi*X/(Lx))\n",
"PY=np.sin(np.pi*Y/Lx)\n",
"Z=[1e-6,]\n",
"PZ=[1.0,]\n",
"\n",
"source=openmc.Source()\n",
"source.space = openmc.stats.CartesianIndependent(openmc.stats.Tabular(X,PX),openmc.stats.Tabular(Y,PY),openmc.stats.Discrete(Z,PZ))\n",
"\n",
"mu=openmc.stats.Uniform(1.0,1.0)\n",
"phi=openmc.stats.Uniform(-np.pi,+np.pi)\n",
"source.angle= openmc.stats.PolarAzimuthal(mu, phi, reference_uvw=(0.0, 0.0, 1.0))\n",
"\n",
"source.energy = openmc.stats.Discrete([E_f], [1.0]) \n",
"sets.source=source"
]
},
{
"cell_type": "code",
"execution_count": 7,
"metadata": {},
"outputs": [],
"source": [
"tallies = openmc.Tallies()\n",
"\n",
"# Grilla\n",
"mesh_xy = openmc.RegularMesh(mesh_id = 1)\n",
"mesh_xy.lower_left = [0, 0, 0]\n",
"mesh_xy.upper_right = [Lx, Ly, Lz]\n",
"mesh_xy.dimension = [100, 100, 1]\n",
"\n",
"mesh_xz = openmc.RegularMesh(mesh_id = 2)\n",
"mesh_xz.lower_left = [0, 0, 0]\n",
"mesh_xz.upper_right = [Lx, Ly, Lz]\n",
"mesh_xz.dimension = [100, 1, 100]\n",
"\n",
"\n",
"# Filtro de grilla\n",
"filter_mesh_xy = openmc.filter.MeshFilter(mesh_xy, filter_id = 1)\n",
"filter_mesh_xz = openmc.filter.MeshFilter(mesh_xz, filter_id = 2)\n",
"\n",
"# Filtro de particulas\n",
"filter_neutron = openmc.ParticleFilter(['neutron'], filter_id = 6)\n",
"filter_photon = openmc.ParticleFilter(['photon'], filter_id = 7)\n",
"\n",
"#Filtro de energías\n",
"filter_energy_neutron = openmc.filter.EnergyFilter(np.logspace(np.log10(1e-5), np.log10(2e7), 201),filter_id = 8)\n",
"\n",
"#Filtro de ángulos\n",
"filter_mu = openmc.MuFilter(np.linspace(-1.0,1.0,50),filter_id=10)\n",
"filter_theta = openmc.PolarFilter(np.linspace(0,np.pi,50),filter_id=11)\n",
"filter_phi= openmc.AzimuthalFilter(np.linspace(-np.pi,+np.pi,50),filter_id=12)\n",
"\n",
"# Tallies\n",
"flux_neutron_xy = openmc.Tally(name = 'flux_neutron_xy', tally_id = 1)\n",
"flux_neutron_xy.scores = [\"flux\"]\n",
"flux_neutron_xy.filters = [filter_mesh_xy, filter_neutron]\n",
"tallies.append(flux_neutron_xy)\n",
"\n",
"\n",
"flux_neutron_xz = openmc.Tally(name = 'flux_neutron_xz', tally_id = 3)\n",
"flux_neutron_xz.scores = [\"flux\"]\n",
"flux_neutron_xz.filters = [filter_mesh_xz, filter_neutron]\n",
"tallies.append(flux_neutron_xz)\n",
"\n",
"spectra_neutron = openmc.Tally(name = 'spectra_neutron', tally_id = 9)\n",
"spectra_neutron.scores = ['flux']\n",
"spectra_neutron.filters = [filter_energy_neutron, filter_neutron]\n",
"tallies.append(spectra_neutron)\n",
"\n",
"leak_salida_mu = openmc.Tally(name='leakage_salida_mu')\n",
"leak_salida_mu.filters = [openmc.SurfaceFilter(surf06), filter_mu]\n",
"leak_salida_mu.scores = ['current']\n",
"tallies.append(leak_salida_mu)\n",
"\n",
"leak_salida_phi = openmc.Tally(name='leakage_salida_phi')\n",
"leak_salida_phi.filters = [openmc.SurfaceFilter(surf06), filter_phi]\n",
"leak_salida_phi.scores = ['current']\n",
"tallies.append(leak_salida_phi)\n",
"\n",
"leak_salida_theta = openmc.Tally(name='leakage_salida_theta')\n",
"leak_salida_theta.filters = [openmc.SurfaceFilter(surf06), filter_theta]\n",
"leak_salida_theta.scores = ['current']\n",
"tallies.append(leak_salida_theta)"
]
},
{
"cell_type": "code",
"execution_count": 8,
"metadata": {},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"/Users/pshriwise/soft/openmc/openmc/mixin.py:70: IDWarning: Another Filter instance already exists with id=3.\n",
" warn(msg, IDWarning)\n"
]
}
],
"source": [
"geom = openmc.Geometry(univ)\n",
"mats = openmc.Materials(geom.get_all_materials().values()) \n",
"mats.cross_sections = openmc_data\n",
"\n",
"\n",
"geom.export_to_xml()\n",
"mats.export_to_xml()\n",
"sets.export_to_xml()\n",
"tallies.export_to_xml()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
"rm: statepoint.110.h5: No such file or directory\n",
"rm: summary.h5: No such file or directory\n",
" %%%%%%%%%%%%%%%\n",
" %%%%%%%%%%%%%%%%%%%%%%%%\n",
" %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\n",
" %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\n",
" %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\n",
" %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%\n",
" %%%%%%%%%%%%%%%%%%%%%%%%\n",
" %%%%%%%%%%%%%%%%%%%%%%%%\n",
" ############### %%%%%%%%%%%%%%%%%%%%%%%%\n",
" ################## %%%%%%%%%%%%%%%%%%%%%%%\n",
" ################### %%%%%%%%%%%%%%%%%%%%%%%\n",
" #################### %%%%%%%%%%%%%%%%%%%%%%\n",
" ##################### %%%%%%%%%%%%%%%%%%%%%\n",
" ###################### %%%%%%%%%%%%%%%%%%%%\n",
" ####################### %%%%%%%%%%%%%%%%%%\n",
" ####################### %%%%%%%%%%%%%%%%%\n",
" ###################### %%%%%%%%%%%%%%%%%\n",
" #################### %%%%%%%%%%%%%%%%%\n",
" ################# %%%%%%%%%%%%%%%%%\n",
" ############### %%%%%%%%%%%%%%%%\n",
" ############ %%%%%%%%%%%%%%%\n",
" ######## %%%%%%%%%%%%%%\n",
" %%%%%%%%%%%\n",
"\n",
" | The OpenMC Monte Carlo Code\n",
" Copyright | 2011-2023 MIT, UChicago Argonne LLC, and contributors\n",
" License | https://docs.openmc.org/en/latest/license.html\n",
" Version | 0.13.4-dev\n",
" Git SHA1 | 922b3f9de6ccf3eb1e61171c3b48d51b62875884\n",
" Date/Time | 2023-08-10 15:33:03\n",
"\n",
" Reading settings XML file...\n",
" Reading cross sections XML file...\n",
" Reading materials XML file...\n",
" Reading geometry XML file...\n",
" Reading H1 from /Users/pshriwise/data/xs/openmc/nndc_hdf5/H1.h5\n",
" Reading H from /Users/pshriwise/data/xs/openmc/nndc_hdf5/photon/H.h5 \n",
" Reading O16 from /Users/pshriwise/data/xs/openmc/nndc_hdf5/O16.h5\n",
" Reading O from /Users/pshriwise/data/xs/openmc/nndc_hdf5/photon/O.h5 \n",
" Reading c_H_in_H2O from /Users/pshriwise/data/xs/openmc/nndc_hdf5/c_H_in_H2O.h5\n",
" Minimum neutron data temperature: 294 K\n",
" Maximum neutron data temperature: 294 K\n",
" Reading tallies XML file...\n",
" WARNING: You are tallying the 'current' score and haven't used a particle\n",
" filter. This score will include contributions from all particles.\n",
" WARNING: You are tallying the 'current' score and haven't used a particle\n",
" filter. This score will include contributions from all particles.\n",
" WARNING: You are tallying the 'current' score and haven't used a particle\n",
" filter. This score will include contributions from all particles.\n",
" Preparing distributed cell instances...\n",
" Reading plot XML file...\n",
" Writing summary.h5 file...\n",
" Maximum neutron transport energy: 20000000 eV for H1\n",
"\n",
" ===============> FIXED SOURCE TRANSPORT SIMULATION <===============\n",
"\n",
" Simulating batch 1\n",
" Simulating batch 2\n",
" Simulating batch 3\n",
" Simulating batch 4\n",
" Simulating batch 5\n",
" Simulating batch 6\n",
" Simulating batch 7\n",
" Simulating batch 8\n",
" Simulating batch 9\n",
" Simulating batch 10\n",
" Simulating batch 11\n",
" Simulating batch 12\n",
" Simulating batch 13\n",
" Simulating batch 14\n",
" Simulating batch 15\n",
" Simulating batch 16\n",
" Simulating batch 17\n",
" Simulating batch 18\n",
" Simulating batch 19\n",
" Simulating batch 20\n",
" Simulating batch 21\n",
" Simulating batch 22\n",
" Simulating batch 23\n",
" Simulating batch 24\n",
" Simulating batch 25\n",
" Simulating batch 26\n",
" Simulating batch 27\n",
" Simulating batch 28\n",
" Simulating batch 29\n",
" Simulating batch 30\n",
" Simulating batch 31\n",
" Simulating batch 32\n",
" Simulating batch 33\n",
" Simulating batch 34\n",
" Simulating batch 35\n",
" Simulating batch 36\n",
" Simulating batch 37\n",
" Simulating batch 38\n",
" Simulating batch 39\n",
" Simulating batch 40\n",
" Simulating batch 41\n",
" Simulating batch 42\n",
" Simulating batch 43\n",
" Simulating batch 44\n",
" Simulating batch 45\n",
" Simulating batch 46\n",
" Simulating batch 47\n",
" Simulating batch 48\n",
" Simulating batch 49\n",
" Simulating batch 50\n",
" Simulating batch 51\n",
" Simulating batch 52\n",
" Simulating batch 53\n",
" Simulating batch 54\n",
" Simulating batch 55\n",
" Simulating batch 56\n",
" Simulating batch 57\n",
" Simulating batch 58\n",
" Simulating batch 59\n",
" Simulating batch 60\n",
" Simulating batch 61\n",
" Simulating batch 62\n",
" Simulating batch 63\n",
" Simulating batch 64\n",
" Simulating batch 65\n",
" Simulating batch 66\n",
" Simulating batch 67\n",
" Simulating batch 68\n",
" Simulating batch 69\n",
" Simulating batch 70\n",
" Simulating batch 71\n"
]
}
],
"source": [
"!rm statepoint.110.h5\n",
"!rm summary.h5\n",
"openmc.run()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"sp=openmc.StatePoint('statepoint.110.h5')\n",
"S0=1.0"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"tally=sp.get_tally(name='flux_neutron_xy') \n",
"mesh_filter=tally.find_filter(openmc.filter.MeshFilter)\n",
"data=tally.get_slice(scores=['flux'])\n",
"\n",
"Nx=mesh_filter.mesh.dimension[0]\n",
"xmin=mesh_filter.mesh.lower_left[0]\n",
"xmax=mesh_filter.mesh.upper_right[0]\n",
"dx=xmax-xmin\n",
"\n",
"Ny=mesh_filter.mesh.dimension[1]\n",
"ymin=mesh_filter.mesh.lower_left[1]\n",
"ymax=mesh_filter.mesh.upper_right[1]\n",
"dy=ymax-ymin\n",
"\n",
"Nz=mesh_filter.mesh.dimension[2]\n",
"zmin=mesh_filter.mesh.lower_left[2]\n",
"zmax=mesh_filter.mesh.upper_right[2]\n",
"dz=zmax-zmin\n",
"\n",
"data.mean.shape=(Nx,Ny)\n",
"data.std_dev.shape=(Nx,Ny)\n",
"\n",
"data_mean = data.mean*S0/(dx/Nx*dy/Ny*dz/Nz)\n",
"data_stdv = data.std_dev*S0/(dx/Nx*dy/Ny*dz/Nz)\n",
"\n",
"plt.imshow(data_mean, origin='lower', interpolation='none', extent=(xmin,xmax,ymin,ymax), cmap='viridis')\n",
"plt.colorbar()\n",
"plt.xlabel('x [cm]')\n",
"plt.ylabel('y [cm]')\n",
"plt.title('flux XY')\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"tally=sp.get_tally(name='flux_neutron_xz') \n",
"mesh_filter=tally.find_filter(openmc.filter.MeshFilter)\n",
"data=tally.get_slice(scores=['flux'])\n",
"\n",
"Nx=mesh_filter.mesh.dimension[0]\n",
"xmin=mesh_filter.mesh.lower_left[0]\n",
"xmax=mesh_filter.mesh.upper_right[0]\n",
"dx=xmax-xmin\n",
"\n",
"Ny=mesh_filter.mesh.dimension[1]\n",
"ymin=mesh_filter.mesh.lower_left[1]\n",
"ymax=mesh_filter.mesh.upper_right[1]\n",
"dy=ymax-ymin\n",
"\n",
"Nz=mesh_filter.mesh.dimension[2]\n",
"zmin=mesh_filter.mesh.lower_left[2]\n",
"zmax=mesh_filter.mesh.upper_right[2]\n",
"dz=zmax-zmin\n",
"\n",
"data.mean.shape=(Nx,Nz)\n",
"data.std_dev.shape=(Nx,Nz)\n",
"\n",
"data_mean = data.mean*S0/(dx/Nx*dy/Ny*dz/Nz)\n",
"data_stdv = data.std_dev*S0/(dx/Nx*dy/Ny*dz/Nz)\n",
"\n",
"plt.imshow(data_mean, origin='lower', interpolation='none', extent=(xmin,xmax,zmin,zmax), cmap='viridis')\n",
"plt.colorbar()\n",
"plt.xlabel('x [cm]')\n",
"plt.ylabel('z [cm]')\n",
"plt.title('flux XZ')\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"dV=Lx*Ly*Lz\n",
"data=sp.get_tally(name='spectra_neutron').get_pandas_dataframe(nuclides=False)\n",
"data.columns=['Emin','Emax','particle','score','mean','stdv']\n",
"Emin=data['Emin'].values\n",
"Emax=data['Emax'].values\n",
"E=(Emin+Emax)/2.0\n",
"dE=Emax-Emin\n",
"\n",
"data_mean=data['mean']*S0/(dE*dV)\n",
"data_stdv=data['stdv']*S0/(dE*dV)\n",
"\n",
"plt.loglog(E,data_mean)\n",
"plt.xlabel('E [eV]')\n",
"plt.ylabel('$\\phi(E)$ [n/(cm² s eV)]')\n",
"plt.grid()\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"dS=Lx*Ly\n",
"data=sp.get_tally(name='leakage_salida_mu').get_pandas_dataframe(nuclides=False)\n",
"data.columns=['surface','mu_min','mu_max','score','mean','stdv']\n",
"mu_min=data['mu_min'].values\n",
"mu_max=data['mu_max'].values\n",
"mu=(mu_min+mu_max)/2.0\n",
"dmu=mu_max-mu_min\n",
"data_mean_mu=data['mean']*S0/(dmu*dS)\n",
"data_stdv_mu=data['stdv']*S0/(dmu*dS)\n",
"\n",
"plt.plot(mu,data_mean_mu)\n",
"plt.xlabel('$\\mu$')\n",
"plt.ylabel('$J(\\mu)$')\n",
"plt.grid()\n",
"plt.show()"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {},
"outputs": [],
"source": [
"dS=Lx*Ly\n",
"data=sp.get_tally(name='leakage_salida_theta').get_pandas_dataframe(nuclides=False)\n",
"\n",
"data.columns=['surface','theta_min','theta_max','score','mean','stdv']\n",
"theta_min=data['theta_min'].values\n",
"theta_max=data['theta_max'].values\n",
"cos_theta=(np.cos(theta_min)+np.cos(theta_max))/2.0\n",
"dcos_theta=np.cos(theta_min)-np.cos(theta_max)\n",
"data_mean_theta=data['mean']*S0/(dcos_theta*dS)\n",
"data_stdv_theta=data['stdv']*S0/(dcos_theta*dS)\n",
"\n",
"plt.plot(cos_theta,data_mean_theta)\n",
"plt.xlabel('$\\mu$')\n",
"plt.ylabel('$J(\\mu)$')\n",
"plt.grid()\n",
"plt.show()"
]
}
],
"metadata": {
"kernelspec": {
"display_name": "Python 3 (ipykernel)",
"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.10.1"
}
},
"nbformat": 4,
"nbformat_minor": 4
}
Display the source blob
Display the rendered blob
Raw
Loading
Sorry, something went wrong. Reload?
Sorry, we cannot display this file.
Sorry, this file is invalid so it cannot be displayed.
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment