diff --git a/docs/source/examples/hexagonal.rst b/docs/source/examples/hexagonal.rst new file mode 100644 index 000000000..8955fdb3f --- /dev/null +++ b/docs/source/examples/hexagonal.rst @@ -0,0 +1,13 @@ +.. _notebook_hexagonal: + +=========================== +Modeling Hexagonal Lattices +=========================== + +.. only:: html + + .. notebook:: ../../../examples/jupyter/hexagonal-lattice.ipynb + +.. only:: latex + + IPython notebooks must be viewed in the online HTML documentation. diff --git a/docs/source/examples/index.rst b/docs/source/examples/index.rst index 7adbe13ec..1633f9950 100644 --- a/docs/source/examples/index.rst +++ b/docs/source/examples/index.rst @@ -8,9 +8,9 @@ The following series of `Jupyter `_ Notebooks provide examples for how to use various features of OpenMC by leveraging the :ref:`pythonapi`. ------------ -Basic Usage ------------ +------------- +General Usage +------------- .. toctree:: :maxdepth: 1 @@ -21,12 +21,21 @@ Basic Usage tally-arithmetic expansion-filters search - triso - candu nuclear-data nuclear-data-resonance-covariance cad-geom +-------- +Geometry +-------- + +.. toctree:: + :maxdepth: 1 + + hexagonal + triso + candu + ------------------------------------ Multi-Group Cross Section Generation ------------------------------------ diff --git a/docs/source/pythonapi/model.rst b/docs/source/pythonapi/model.rst index 14f4a4beb..ee038987b 100644 --- a/docs/source/pythonapi/model.rst +++ b/docs/source/pythonapi/model.rst @@ -12,8 +12,8 @@ Convenience Functions openmc.model.borated_water openmc.model.cylinder_from_points - openmc.model.get_hexagonal_prism - openmc.model.get_rectangular_prism + openmc.model.hexagonal_prism + openmc.model.rectangular_prism openmc.model.subdivide TRISO Fuel Modeling diff --git a/examples/jupyter/expansion-filters.ipynb b/examples/jupyter/expansion-filters.ipynb index 7870faffe..69d315e24 100644 --- a/examples/jupyter/expansion-filters.ipynb +++ b/examples/jupyter/expansion-filters.ipynb @@ -62,7 +62,7 @@ "source": [ "# Define surfaces used to construct regions\n", "zmin, zmax = -10., 10.\n", - "box = openmc.model.get_rectangular_prism(10., 10., boundary_type='reflective')\n", + "box = openmc.model.rectangular_prism(10., 10., boundary_type='reflective')\n", "bottom = openmc.ZPlane(z0=zmin, boundary_type='vacuum')\n", "boron_lower = openmc.ZPlane(z0=-0.5)\n", "boron_upper = openmc.ZPlane(z0=0.5)\n", diff --git a/examples/jupyter/hexagonal-lattice.ipynb b/examples/jupyter/hexagonal-lattice.ipynb new file mode 100644 index 000000000..e0a48e236 --- /dev/null +++ b/examples/jupyter/hexagonal-lattice.ipynb @@ -0,0 +1,372 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "In this example, we will create a hexagonal lattice and show how the orientation can be changed via the cell rotation property. Let's first just set up some materials and universes that we will use to fill the lattice." + ] + }, + { + "cell_type": "code", + "execution_count": 1, + "metadata": {}, + "outputs": [], + "source": [ + "%matplotlib inline\n", + "import openmc" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "fuel = openmc.Material(name='fuel')\n", + "fuel.add_nuclide('U235', 1.0)\n", + "fuel.set_density('g/cm3', 10.0)\n", + "\n", + "fuel2 = openmc.Material(name='fuel2')\n", + "fuel2.add_nuclide('U238', 1.0)\n", + "fuel2.set_density('g/cm3', 10.0)\n", + "\n", + "water = openmc.Material(name='water')\n", + "water.add_nuclide('H1', 2.0)\n", + "water.add_nuclide('O16', 1.0)\n", + "water.set_density('g/cm3', 1.0)\n", + "\n", + "mats = openmc.Materials((fuel, fuel2, water))\n", + "mats.export_to_xml()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "With our three materials, we will set up two universes that represent pin-cells: one with a small pin and one with a big pin. Since we will be using these universes in a lattice, it's always a good idea to have an \"outer\" universe as well that is applied outside the defined lattice." + ] + }, + { + "cell_type": "code", + "execution_count": 3, + "metadata": {}, + "outputs": [], + "source": [ + "r_pin = openmc.ZCylinder(r=0.25)\n", + "fuel_cell = openmc.Cell(fill=fuel, region=-r_pin)\n", + "water_cell = openmc.Cell(fill=water, region=+r_pin)\n", + "pin_universe = openmc.Universe(cells=(fuel_cell, water_cell))\n", + "\n", + "r_big_pin = openmc.ZCylinder(r=0.5)\n", + "fuel2_cell = openmc.Cell(fill=fuel2, region=-r_big_pin)\n", + "water2_cell = openmc.Cell(fill=water, region=+r_big_pin)\n", + "big_pin_universe = openmc.Universe(cells=(fuel2_cell, water2_cell))\n", + "\n", + "all_water_cell = openmc.Cell(fill=water)\n", + "outer_universe = openmc.Universe(cells=(all_water_cell,))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now let's create a hexagonal lattice using the `HexLattice` class:" + ] + }, + { + "cell_type": "code", + "execution_count": 4, + "metadata": {}, + "outputs": [], + "source": [ + "lat = openmc.HexLattice()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "We need to set the `center` of the lattice, the `pitch`, an `outer` universe (which is applied to all lattice elements outside of those that are defined), and a list of `universes`. Let's start with the easy ones first. Note that for a 2D lattice, we only need to specify a single number for the pitch." + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [], + "source": [ + "lat.center = (0., 0.)\n", + "lat.pitch = (1.25,)\n", + "lat.outer = outer_universe" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now we need to set the `universes` property on our lattice. It needs to be set to a list of lists of Universes, where each list of Universes corresponds to a ring of the lattice. The rings are ordered from outermost to innermost, and within each ring the indexing starts at the \"top\". To help visualize the proper indices, we can use the `show_indices()` helper method." + ] + }, + { + "cell_type": "code", + "execution_count": 6, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + " (0, 0)\n", + " (0,11) (0, 1)\n", + "(0,10) (1, 0) (0, 2)\n", + " (1, 5) (1, 1)\n", + "(0, 9) (2, 0) (0, 3)\n", + " (1, 4) (1, 2)\n", + "(0, 8) (1, 3) (0, 4)\n", + " (0, 7) (0, 5)\n", + " (0, 6)\n" + ] + } + ], + "source": [ + "print(lat.show_indices(num_rings=3))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Let's set up a lattice where the first element in each ring is the big pin universe and all other elements are regular pin universes. From the diagram above, we see that the outer ring has 12 elements, the middle ring has 6, and the innermost degenerate ring has a single element." + ] + }, + { + "cell_type": "code", + "execution_count": 7, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + "HexLattice\n", + "\tID =\t4\n", + "\tName =\t\n", + "\tOrientation =\ty\n", + "\t# Rings =\t3\n", + "\t# Axial =\tNone\n", + "\tCenter =\t(0.0, 0.0)\n", + "\tPitch =\t(1.25,)\n", + "\tOuter =\t3\n", + "\tUniverses \n", + " 2\n", + " 1 1\n", + "1 2 1\n", + " 1 1\n", + "1 2 1\n", + " 1 1\n", + "1 1 1\n", + " 1 1\n", + " 1\n" + ] + } + ], + "source": [ + "outer_ring = [big_pin_universe] + [pin_universe]*11\n", + "middle_ring = [big_pin_universe] + [pin_universe]*5\n", + "inner_ring = [big_pin_universe]\n", + "lat.universes = [outer_ring, middle_ring, inner_ring]\n", + "print(lat)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now let's put our lattice inside a circular cell that will serve as the top-level cell for our geometry." + ] + }, + { + "cell_type": "code", + "execution_count": 8, + "metadata": {}, + "outputs": [], + "source": [ + "outer_surface = openmc.ZCylinder(r=4.0, boundary_type='vacuum')\n", + "main_cell = openmc.Cell(fill=lat, region=-outer_surface)\n", + "geom = openmc.Geometry([main_cell])\n", + "geom.export_to_xml()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Now let's create a plot to see what our geometry looks like." + ] + }, + { + "cell_type": "code", + "execution_count": 9, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZAAAAGQAgMAAAD90d5fAAAABGdBTUEAALGPC/xhBQAAACBjSFJNAAB6JgAAgIQAAPoAAACA6AAAdTAAAOpgAAA6mAAAF3CculE8AAAADFBMVEX///8AAP///wCAgACerKf2AAAAAWJLR0QAiAUdSAAAAAd0SU1FB+MGCAIxFif9eE4AAAhVSURBVHja7Z3LdeM4EEXlBUJgPgyBC9E8xxvtlQSj6BBmMc5HoWjZR2OLlkVK9Xn1IeRRd626DdDXVQ8ECiAIbjYmK/3ZOttVJutnVgGxEualv7M2m1F6wrpcRtOTtq3ASKWUnrW0iAmMNMpLL1qbAukVy2A0GiRB/KIxMmTRGfGAqcFKCNgLwoi2MIwRCxigelx7lBFxBXYk4grO8LticMTvioXhdcXkiNcVG8PnitERnyuNFeLowcBea27t+o44XHE4YnfFLLtHeg/D2opd0bLGyyH7p9mk9zFs8SpeiEV6Z7RM8XLK/mm49O5oWeLljpYhXoFo4fEKRAuPVyBacLxC0ULjFYoWGq9QtNB4xRhY/xWUBBMlGC0sXlEGEi8tWq/v77/C8SoqQ6XojViRZP8J+TcqCuCI6kpQkrcJ8k9MFEWS9y+LiSJL8nqByPHSRIGipcYrJMn+AlHalywKKElMlCYHIosiX/t6hfzyi4LdJcE7pWRBJFEUSfZXSKD76rMggijaWPL+jjYvQZSSB+FFafIgvCh9HqT3Ql7nkF9OCJBD4BBO+ZIJ4ZRvMiGc8n0mhBFFTYJtkPZxkJILoZVvciG08tpVRkjvg5i6FRoCzLBsEEr5kg2hlG+yIVsfxDD8MhD9KiOEUB6ZWeMp0dnax0AKcBWepp7tvnk1wFVwwj3Z9jEQ5Cp4EvRlHt3x6dyXtS4IODHlIEWouzsdL/+8lWRWRNht82r4qsPpdLqJ13e0PooO/JVbHLL7+E2Xv/dm2WNeBEAEp8eP3/T78p/9om0tigjDIafTLF7LpahFkQoRGtew/E2LRbVzkSBKa4McjEUUpPAVdxrkyF+7bMONBjkaiybb1ofw9c7NVITAbbgGROoeQ5C2NqQI9SLCL9pwFUizFmQLQiLdygIiVItBehBy05/v56mK1tXPIXIWsWimb4vRd1QhLQgZ5+Le5hHybTKHFLHeMIv7TSIxaJLM2rAM6W+H+OUgL1/bQS34/Pd+//MuFVYcmbVhBXI1Ywa5gKBXWHPhT7NDjFm9D2Kbnywg8NNe43TubG1NSLkvou9jbfY7EjdNx0J29CihQMireAjTJSkrEiPVW3bcvTgwnasCIQeXLQfZMWORvEo0kGPxlrtNRmboliE7Oj8SIb+tkNEGOfkhJwbyQkPuqyurqcxV7QMhQwRyoCFlXUj3ZJBmXciWSSNSW1dFSE9Xz7rjewaS2neJkKMVsrNBMseTCfLCVL//sWtknPqVF7p61hjPQxKzlQlSiJ8n5l1Tv0JCaPNkkD8T4snqJ0hjqL+XGhdnWyPkTdQ9CWKfM8KQyOwXhkTm8V8QwJHAisSnIZDQ2goKCa0SWSDu9S4QElu5myDq/DoOaXVIaDX1z4OE1uq/IGV9SPcnQeLC/yzIgSpSnzPCkHi3AkO+/2t69nuBNFqVyFPsybYAZDy5n8fjkN3Ju7PAAPHvkTBA/Ls9LJCruTJIK8SVC1shrqzeCnHNT2pBLLV907mPxMtSWZv92iADfY8pkIG7M2nI6F2RoDsyGuJfW8Ehg3+ViIwXCdn517uO1O8jIbmLahIkbXmQgTA5m7KaOrKDMQ/JWrJlILkr3H8hfyFEdR+E+n2PuxmrdCtVOsgqXX2VQavK8FsnkaiSEjGWm9wx5k5TG0ttqXHx9hOnDnupcWVB3ly6/5yJaY0pdpXFgirLHlUWcBKWoopWJWNRDYQcjEVGSJWFzudZF36etfoMSJWHND/jcdPzQDIeZlZ5LGt9wIwWOSBMMx3TICMv7g4aGa3bF/CiJaQB4sX+tSO4EUOHDMIvSttSErTngpS1Idb9Xf9DiHAfj0IRu7HPsEXRV2TcbDn9tVKRYbPlwHeugzCCMEWtbQNsL603iRtgDVt5ZQi3qGbclOwrMm6vBiDc9uqGhtDK24uMu9FTIUMEcqAhZV2I8V0HiJ8AcXsCv+Tigxhf1/EVVYTgr1D1njve+jKYq8j6Wlvv6YWtL+j1jvGEfQswc2S8QApVnRvIhWkPOcabX//0FPEQKe8SZiTSi6zhU54lM7/3G4EkHPPM2+b5IM16DPhN/8js9wopcsXIPL4DIaEViSvEcFjF0vS1lRaEhFaJrhDLASJo0ZdtMEjaUSjNWhD05JjzgJtxPE1ZC9LVhjzPCUvPc+rVE51EVjTIgSoyHtxW5Qg6/DC9hdkg8LGAcNHZNjBklActwwGHDV9zJw+/R/5Kw3mQgxASWZJbSBGcFrZ02I7PXCWzbx8BWSXpvmWskXQ7z7KNQko+pLuDrKB864YMyJNFDoI2rxF5RjrZPQNVXk2xv815DncvL+AAkAJduVMG3Jl1BARTftSSoKsRuoPKGyAbN+SkjbgKpMmF0Mf7l1xIR0Ig5XFI+zgIpDwOoRmQ8jCE+6BHyYR0DAQRBYYwkkCiwJBNAAJ3KyykyYPwn/Ap+sVoV9+xEEB5dNBidYeVB6rxDFCUkCSQKFhK1AmQtAxPkCQvt5cYWWm3/Em4kgPpREiSKKIkWaLIjBxRtE8nlgxIp0BSRFEkyRFFY2SIon9itMQhmiQpoqiS1PmAbZ1P8Vb5qHCVzyPX+dBziUH0BpwgCiRJnc+I1/kgepVPu1f5SH0kXmi0QvGCoxWIFx6tQLzwaG38/ZeF4Y2XJVpu6Q2y++NlY/ikN8nujZcxWi7pbbI7XTE74nDF7ohDeqvsZ7NCPAyrKy5HrK74GDZXnI7YXPEyLK64HbG44mfgrgQcwV2JMNAerA1BsB7M02uZAxZlINp3YYgesHCwkIBlMLQW1qZAZFm6HIZISWMI4qeIrlBSGUzEEmM1GdHG2mzG5v5+WQFxi1kJ8WnFJcZ/bkPEpRmBEaUAAAAldEVYdGRhdGU6Y3JlYXRlADIwMTktMDYtMDhUMDc6NDk6MjItMDU6MDA/X9I4AAAAJXRFWHRkYXRlOm1vZGlmeQAyMDE5LTA2LTA4VDA3OjQ5OjIyLTA1OjAwTgJqhAAAAABJRU5ErkJggg==\n", + "text/plain": [ + "" + ] + }, + "execution_count": 9, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "p = openmc.Plot.from_geometry(geom)\n", + "p.color_by = 'material'\n", + "p.colors = colors = {\n", + " water: 'blue',\n", + " fuel: 'olive',\n", + " fuel2: 'yellow'\n", + "}\n", + "p.to_ipython_image()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "At this point, if we wanted to simulate the model, we would need to create an instance of `openmc.Settings`, export it to XML, and run." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Lattice orientation\n", + "\n", + "Now let's say we want our hexagonal lattice orientated such that two sides of the lattice are parallel to the x-axis. This can be achieved by two means: either we can rotate the cell that contains the lattice, or we can can change the `HexLattice.orientation` attribute. By default, the `orientation` is set to \"y\", indicating that two sides of the lattice are parallel to the y-axis, but we can also change it to \"x\" to make them parallel to the x-axis." + ] + }, + { + "cell_type": "code", + "execution_count": 10, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZAAAAGQAgMAAAD90d5fAAAABGdBTUEAALGPC/xhBQAAACBjSFJNAAB6JgAAgIQAAPoAAACA6AAAdTAAAOpgAAA6mAAAF3CculE8AAAADFBMVEX///8AAP+AgAD//wDoUCWoAAAAAWJLR0QAiAUdSAAAAAd0SU1FB+MGCAIxFif9eE4AAAeGSURBVHja7Z1NkuMoEIWrFhyB++gILCxNhDZezUaX8CnqCF607tNH6WWFpmxXWfDEX0KSdmucO4dxfmQ+hARG8PZGMtVfzdB+RbLeMgFEI8x7v7GOm6F6jxlehu69dhBgsFJUHzS2jEUYbJT3PmodC6RPGAdDpyAM4qsUg0OWNKM+YclkMSTsPYdR28LyGHUJy1C9XvtcRk0o2YHUhJLPKA+FEEh5KBRGaSikQEpDoTHKQiEGUhaKpkIKerDMXsu2rn0gBaEUBEIPhSx7ifQlDGorLsoWNV8Fsl+MJn0Zg5YvVQqhSF+YLVK+CmW/WL70xdmi5Ks4W4R8VWQrP18V2crPV0W2svNVla3cfFVlKzdfVdnKzVcdI6//qpQkT5TKbOXlq5aRk6/qbOXkS9VD0o24WpIcUeoZaVEYJEmLojggKVE0ByQlCgcjJQqLJClRFA8kLormgcRF4WHERWGSJC6K4oLERNFckJgoXIyYKGySxERRfJCwKJoPEhaFjxERRQLCqHtYecUJCSmvOSEh5TkZIVFYJQmJIgJRvBC/8poX4lc+UHiY3c/TL/fz/BH4IQUyz47XYXa9TvNMgAR0H2fX6+RCv5jz2f9Ln/LKX3QCL18f7apf6vDL/0uf8jqULcfLpeZ2ZBNALTtkQ65OLS+XmtuRzQBNQPosyASRRSAe5SO6214AcqvD2f/brhQyu5ERIcpfcAIvALnVIbt5aX7I4TGQnh+yaV7vLSDdIyCqBQSblw6Uq7lONsqLQPpMCKVb2TSvULGaDhIh4cchcIpdfeR+0mPzCkPgzog3rcidcQNRwXJ4j4eaDzFJoA3rYLkB0jGC05gk0LzCkK98uTUFCcZItgDSNzJxCPMD/WqdNES1ghhpiG4FOUhDWjGcNiwBadaC7TYsAlHtIEYWottBDrKQdgyrDe8G0vAyWS8UEYhqCTG7g+iWkIMkpCXjfqHsBhK6FmEANcDT+/QRLX63LgaZYEwA45IRhhBTaAgRgwwwuhlhMAQjrCE4GOoi1yIOBmGsOMCwbgoO60wEAsPaAUa9I0DDsx8xyOyBWF4mD2SOQLTvqwFmAUbwAjMRWNyyQwqyepl8kDWy8Skg3kyOXsjdy+CFnL2uUpDVy+yDrJFN/xvIe2tI9wyQM3gBSF7r6qL9IxPEBCF8F2MS8gFeAHIOFs+D9D4I9sII8Xsyb9F/5Vi6+mvnFYBMnjvjWvNbZJHieZARWsu0vf3+ihTPg+A06Qg1nbYPEgFHMci4eSSymRevkN5AIFdI39xekGeENB1f36x7QZ4QotpDzAvygjwR5J9l+W1/Pi6fzvfL8mcL0WSGQzl+fbQpl6+RcqBCTovrdXGhx+vnSsg1EMvrzekKvX0NoVAhR/CyuNDvOrgqkSHfThfX6R16AmgR5Mfpj5fvwO5VXxZfvoiQI3g5uZHd6/BZAzmBl8WN7F6HpQayuF7uNf+O7ARQFsha8884hMJYa37zcgLoApH9GAlyBC8AserwmYCM4TEBQhY3MgJkDo9uTuAFIFYdljgEx2mTBQWIVfNrZPkQHHHao97F9YIQqw5u89pAYOw8NIHMkVkAgNjp+ayBOJMmXBD/9M8NYqfn4gUh9td/KJBJDnJ+QV6Qv7UJN7kYRboVmQ5SpKsXuWmJ3H5FHiRkHok29m8fgpzcmtdAVit+TNUUCjg9Qs1Pbh1+7BmHDie35m0GQUe35m2GcyIDU5Ehtsxkgci0h8gETuFUlKJByibVqJACe0FekJaQ/fwT9IKQIPv57/cFoUF04CvWhRghiMSSEpHFMSLLfHgXLKkIZPVSufTKDxFZRCayHG70Qu5eqAv7RJYoPnyxJRtkP6tsHw45gxeAEJZXa983f+Vq9AjkA7wA5BwsvoEobybZuvqHv1Axzqyvhoi85PLA13WYXzzaz3taMhDdkiH6FqBqCdnf65/7ee93P+9iy0B0O4bwm/6qHcTIQvazI8Z+dinZ0c4xqhXESEP2s8PSfna9ktkkLAxh3O5MBcuN4AWGcZSN20S2oBPZTE9kW0CZDQ51JgSmoJ5wP0gVKIezZQDBuTfXDED2s9uoyOasO9rLVvFDzAYisomxyHbMIhtLy2yR7Ycwb/at/EVxmhScjpTGJbQBu8hW8jKb4uue1fzb+yteiPFC9nOug8gxGDIHeihOiAlARA5Z2c+ZNDJH+Cg+iAlCRI5VEjkgSuaoK8UFMRGIyPFjIgepyRwJp3ggJgoROaZP5MBBmaMTFQfEJCAix1mKHMwpc8SoqoekJBE69lXkAFuZo3hFDhUWOR5Z5qBnVQdJN2AGUbIkkTlGXOZAdJGj3UUOqa/JV262qvKVna2KfOVnqyJf+dl6K++/KIzSfFGyVSw9QfbyfNEYZdKTZC/NFzFbRdLTZC8MhRxIQSj0QAqkp8p+NSqkhEENpSgQaihlDFoohYHQQillUEIpDoQSSjkjP5SKQPJDqWHk9mBdFSSvByvptcgJq2XkaG+qIemEVScrJ2EcjFQL61ggcVkMDyNKYWNExGcRPUFhZQQyxpirm3naWMfNeNteLw0QiGmEuJgqEuM/rS4lqiAX++MAAAAldEVYdGRhdGU6Y3JlYXRlADIwMTktMDYtMDhUMDc6NDk6MjItMDU6MDA/X9I4AAAAJXRFWHRkYXRlOm1vZGlmeQAyMDE5LTA2LTA4VDA3OjQ5OjIyLTA1OjAwTgJqhAAAAABJRU5ErkJggg==\n", + "text/plain": [ + "" + ] + }, + "execution_count": 10, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "# Change the orientation of the lattice and re-export the geometry\n", + "lat.orientation = 'x'\n", + "geom.export_to_xml()\n", + "\n", + "# Run OpenMC in plotting mode\n", + "p.to_ipython_image()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "When we change the orientation to 'x', you can see that the first universe in each ring starts to the right along the x-axis. As before, the universes are defined in a clockwise fashion around each ring. To see the proper indices for a hexagonal lattice in this orientation, we can again call `show_indices` but pass an extra orientation argument:" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [ + { + "name": "stdout", + "output_type": "stream", + "text": [ + " (0, 8) (0, 9) (0,10)\n", + "\n", + " (0, 7) (1, 4) (1, 5) (0,11)\n", + "\n", + "(0, 6) (1, 3) (2, 0) (1, 0) (0, 0)\n", + "\n", + " (0, 5) (1, 2) (1, 1) (0, 1)\n", + "\n", + " (0, 4) (0, 3) (0, 2)\n" + ] + } + ], + "source": [ + "print(lat.show_indices(3, orientation='x'))" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Hexagonal prisms\n", + "\n", + "OpenMC also contains a convenience function that can create a hexagonal prism representing the interior region of six surfaces defining a hexagon. This can be useful as a bounding surface of a hexagonal lattice. For example, if we wanted the outer boundary of our geometry to be hexagonal, we could change the `region` of the main cell:" + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [ + { + "data": { + "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZAAAAGQBAMAAABykSv/AAAABGdBTUEAALGPC/xhBQAAACBjSFJNAAB6JgAAgIQAAPoAAACA6AAAdTAAAOpgAAA6mAAAF3CculE8AAAAElBMVEX///+TUVD//wCAgAA4vPIAAP/pte6jAAAAAWJLR0QAiAUdSAAAAAd0SU1FB+MGCAIxFif9eE4AAAn2SURBVHja7d3tceM4DAbgsANJHVh0BVIFl7E7WPffymV346wlEh8UCVDGAP9u3jkvHhP+JCN/fHh5eXl5eXl5eXl5eXl5eXl5eXl5eXl5eXl5eXl5WatRtPQcQRYyOKS0ZB16syW8IHpLYgYi7VCbLSuQIA8ZHFJS8g6d2VJYEJ0lMQPRcKjMlhVI0IEMDuGWjkN+trILMl/QniYinnssSQ4yz6hkIuOcRBoCdDLjjWISKFaH/GkUlnzHoASMZR0BahSUMONUMuguyI8DaJUdpxJdyIy1solzkgmLJR0BazQj2cQz6kjjQRGybTSRTEXxrAjBHftWp8J4L1GDpJ3MRLyRELEkJFCOV0k2vmCOvWTQWZBsoy+SQ/FWogOZ8VaAeMbvhlkFEhiOpwSMZ9SxlQwKELjRP5LpcDwrQHiO361OFfGrRByCdzIT8YWIpSGBuSDVJT1b3MlqCJFZEisQ5nNWW4jEbJmBwK+GcV2xvr7iBcuTWPbFHX5ZX38X3GpkxNs7QvbFHYSsKyqJrHjVg4wA5NkJNF3MeAEg7WcLgqwruiSRGa9akABA/nWSv8/Z8QJAWs/WQUhc8RXTh4wAZMU7jex4BSCNZyschKwrumIcSNslgSCR6LQgXlQgYz9I29myAgk9IS1nywxk7AlpOFuhL6TdksAQhdeRlpCxL6TdbB2HRHasAQkIJKL3eIM3jS1nC4NQn0dwJ+PzSEvIiEEi2ij1wYr+hNhwtgIKoT6z407yM3vLJSEgEW30KSHizd0gBRlxyF/JPKOtEvF2OTP/oAqkeQlBwhkgLWZrPAOkxZJYgYRzQOpnywwkd6uTLOQiMlv6DpmzmvpPvr9LYLa6OCReE/s42r8DDp0czT+UJJBJC3JpC+nmyEhaQhQdqaThZKk6EknNbHV6oD+r3ZL0dbTbhAudHc32RTeQqQfk0gbS3bGTtIB0cmwlDSarm2MjOTpbnR/oz6pfklMsSIPjs+eYrAazZQZyhuesHeTQbIUzQo4sCQLBT5XOET9VSsTJrddCwJf1iO3rzN87Hofj9NZrX9whSMR3qGZ8h4qKM7deCQHfMB7YFEwbJbccIUj5bEGQ7D4sO2ZvAr/ElZARgLD3m5cDcfbWN5Di2QoAhH8CYD0QZ299CyldEgaEf7gkC1m4cR1kPCukdLYgyMqHLMVx/tarIIEDWfBOCMjKjXeQstkyAxnPCymarXBmSMmSmIGMLMiMt7IUx6xnraLZsgIJMCTyITPe6cKNEwh/tniQhej0QBwbQ0YYMrMhC9EpO04g7NkKPMiMt7IciBlv40uWBIVEvJMZd1Jx7taPQ0YMwj02uhyKc7eeQrizhUNmvBMqjujdkIsPQwIBmdFGyThijlycgfBma6QgypWB8JbECiS8A4QzW2Yg4ztAOLNlBZKbrMq90OsnHt/w//2W64ieLQHH7YZKyDgr6fDk+9XJjYgRydfdcPvvACSIOBDJ3xiUXP/EGQk1W0HGAUquvPi/YoiUA2j1yo1TSSmkzvGvk+z0vMRZyUtcCAlijpxkE2ckm7hstgQdGUlRfCtaEklH0mphfCuA7CZrqoMkndyIeCNJ7ob9A35gQ9o47o+v+pVKnu29xi+S60v8yEowSHvHH8ZLr/tGmfEjJ2EuSFvHvtV9o1Sckww8SJ3junc8W/0EGt1K0viRPuBhSPMnrEemlU+g0VfoLb0bnlDWbMk6XiS5Rsk4lTAma6qD5Dv5mZ5b3vkdfwJ3Q/owGWiIyIL8PI0CTipOlgSCSD7SN61ATiJOH+/0c1YLyANsBXJS8a89ZDgD5AF2Wgtp/OR7xzvNxg/c+eA9ATeGPPBO70fiXxwI/LK+lm4EcCAPuNMCyFACWYu3Zq5gJ/86PRQ/OJARgDD30tYyyB3ptAAyEguSvThQHpLbpr3CnTx4ECj+tYcMbMiR/eaukBGA8M8qFEEeSKclkLEV5CVeTgEJACRmO6Uhd6LTg/EjgQxtIPF0kLEesnaBjOiCvBNkcIhDNCBmHux2nn7NvCCaeYti501jACARnawTvo2HIG/3wcrMR107Xz6AkHf7OqjiC7q4+28O5H4kTiAfbSFzDnIHO+FAar4yDW0h/b7EVobcy2MmxMxGj52ttxfIVCcBWvnuRHwz1Mz2tJ0DA6G15A40Kn2EYwOZ6iTJnbppVPhQjZljTnYOnm0hb3wU0MzhTDvHZYOg5BOP2x5gNnOk3M4h/yAlIZaMiFPHUAx50z+EMfOnSXb+WCwISD7xWObP98z8QaUdSHgHCGOy7EDM/GG+HUigIMi2DyOOpZsUGQhrsuxcTsTOBV4CBvnZM1zwTs5wyR0Ugm+7n+siSGYuS2XnQmE8yBtcus3MxfTsQAILsuCdEBC/cmYJxMxFWe1cJtcMxMylpO1c3DucFVI4WXYugG/mJwns/EgEBHm7n+0w80Mqdn7aJkCQd/uxIQSiW7UQMz/IZecn0sxAzPyMoJ0fdjzHbNVP1jmet6qfs3aQ9/452tBdUvuGMQd5559sNvMj2nZ+1nw7W2/8Q/PwV6f6jpoFwfYX1B11kNBNkjiqJiuzwasFGdtCMocHOjnqJit7wqaPoxYSukgy/2jlZPU59HTJ/KO1jh5HHHOOekhutoSHK/cvVk+WHYjA0ewjkHqHHYj+E7DMZHV4TZRZEBoS8X2CmREvKpCAQyK+czOz4o1EaLIoCL7Vxj2UqgHB3wFzNwUBaC4WmiziQ8mB/eacc8UgjRYEhUS8U/75gkUBgv2V5Yp3GtnxikBaOSog64qumDYk9IU0mywEEolOC+JFAfLRF9LOwbuwphSk4YLYgbCuECoFaemwAwn9IE0nC4QovI60hXCuELrgnS4HIW0djCuECr1pbLwgRyHU23h9CHiFUPwepz5YxZxTcrLoS51Cn2WZ8aIGCQCE+sweWTH4mb35ZFGXOoUaJQ+lxnS9ZCHIBeoifuo0rvix0iQWnayWF9akSnZB7ED0TgbKTpYhSNCHiExWwfFZIr4QsfCC8GfrgkOJWHyy2M9bF3zJiFj8OeuDe6nTCwElZlP4Zf1vjQzJhYASj7LNPyHlYF3qlIAS8UUHEmgJEb80SjsGMQj4h6/ZyZiIRimH3ILQx2cvOJSIRb/Q2lbAJbt7dN8qEYt+VUpARrTRPTSJJzSWhKCXOk0dW0kmnrBY0oFejGfMFTtO7wZRSDJb/yTjiLZKxKljEIXA194CGn3GFygHY1kHeFZzBGtCHX8luVgYEvKtII3+lhCx2AmtQsiINvrVKhFLHWLEa1QqaUd+SdqX+ILYgSjNlrzDDiRoOAaHFJQGRMOhsSQqC2IHojBbOg47kCDtGBxSWNIQLYf0kqgtiB2I8GzpOby8vLy8vLy8vLy8vLy8vLy8vLy8vLy8vLy8vLy8lOp/nb6TXPV6EWgAAAAldEVYdGRhdGU6Y3JlYXRlADIwMTktMDYtMDhUMDc6NDk6MjItMDU6MDA/X9I4AAAAJXRFWHRkYXRlOm1vZGlmeQAyMDE5LTA2LTA4VDA3OjQ5OjIyLTA1OjAwTgJqhAAAAABJRU5ErkJggg==\n", + "text/plain": [ + "" + ] + }, + "execution_count": 12, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "main_cell.region = openmc.model.hexagonal_prism(\n", + " edge_length=3*lat.pitch[0],\n", + " orientation='x',\n", + " boundary_type='vacuum'\n", + ")\n", + "geom.export_to_xml()\n", + "\n", + "# Run OpenMC in plotting mode\n", + "p.color_by = 'cell'\n", + "p.to_ipython_image()" + ] + } + ], + "metadata": { + "anaconda-cloud": {}, + "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.7.0" + } + }, + "nbformat": 4, + "nbformat_minor": 1 +} diff --git a/examples/jupyter/pincell.ipynb b/examples/jupyter/pincell.ipynb index 72b65a265..f8641db13 100644 --- a/examples/jupyter/pincell.ipynb +++ b/examples/jupyter/pincell.ipynb @@ -903,8 +903,8 @@ } ], "source": [ - "box = openmc.get_rectangular_prism(width=pitch, height=pitch,\n", - " boundary_type='reflective')\n", + "box = openmc.rectangular_prism(width=pitch, height=pitch,\n", + " boundary_type='reflective')\n", "type(box)" ] }, diff --git a/include/openmc/particle.h b/include/openmc/particle.h index 62923c11f..a2b7ce557 100644 --- a/include/openmc/particle.h +++ b/include/openmc/particle.h @@ -9,6 +9,7 @@ #include // for unique_ptr #include #include +#include #include "openmc/constants.h" #include "openmc/position.h" @@ -38,7 +39,13 @@ constexpr double CACHE_INVALID {-1.0}; // Class declarations //============================================================================== -struct LocalCoord { +class LocalCoord { +public: + void rotate(const std::vector& rotation); + + //! clear data from a single coordinate level + void reset(); + Position r; //!< particle position Direction u; //!< particle direction int cell {-1}; @@ -48,9 +55,6 @@ struct LocalCoord { int lattice_y {-1}; int lattice_z {-1}; bool rotated {false}; //!< Is the level rotated? - - //! clear data from a single coordinate level - void reset(); }; //============================================================================== diff --git a/include/openmc/position.h b/include/openmc/position.h index 728763f25..2be72b27c 100644 --- a/include/openmc/position.h +++ b/include/openmc/position.h @@ -3,6 +3,7 @@ #include #include // for sqrt +#include #include // for out_of_range #include @@ -55,13 +56,16 @@ struct Position { //! Dot product of two vectors //! \param[in] other Vector to take dot product with //! \result Resulting dot product - inline double dot(Position other) { + inline double dot(Position other) const { return x*other.x + y*other.y + z*other.z; } - inline double norm() { + inline double norm() const { return std::sqrt(x*x + y*y + z*z); } + //! Rotate the position based on a rotation matrix + Position rotate(const std::vector& rotation) const; + // Data members double x = 0.; double y = 0.; @@ -91,6 +95,8 @@ inline bool operator==(Position a, Position b) inline bool operator!=(Position a, Position b) {return a.x != b.x || a.y != b.y || a.z != b.z;} +std::ostream& operator<<(std::ostream& os, Position a); + //============================================================================== //! Type representing a vector direction in Cartesian coordinates //============================================================================== diff --git a/openmc/__init__.py b/openmc/__init__.py index 715550302..1661c4ac4 100644 --- a/openmc/__init__.py +++ b/openmc/__init__.py @@ -32,6 +32,6 @@ from openmc.polynomial import * from . import examples # Import a few convencience functions that used to be here -from openmc.model import get_rectangular_prism, get_hexagonal_prism +from openmc.model import rectangular_prism, hexagonal_prism __version__ = '0.10.0' diff --git a/openmc/cell.py b/openmc/cell.py index 8381fa8d8..3fd70b345 100644 --- a/openmc/cell.py +++ b/openmc/cell.py @@ -227,10 +227,6 @@ class Cell(IDManagerMixin): @rotation.setter def rotation(self, rotation): - if not isinstance(self.fill, openmc.Universe): - raise TypeError('Cell rotation can only be applied if the cell ' - 'is filled with a Universe.') - cv.check_type('cell rotation', rotation, Iterable, Real) cv.check_length('cell rotation', rotation, 3) self._rotation = np.asarray(rotation) diff --git a/openmc/mgxs/__init__.py b/openmc/mgxs/__init__.py index e505c7ebc..886cc904a 100644 --- a/openmc/mgxs/__init__.py +++ b/openmc/mgxs/__init__.py @@ -7,6 +7,7 @@ from openmc.mgxs.mdgxs import * GROUP_STRUCTURES = {} """Dictionary of commonly used energy group structures: + - "CASMO-X" (where X is 2, 4, 8, 16, 25, 40 or 70) from the CASMO_ lattice physics code - "XMAS-172_" designed for LWR analysis ([SAR1990]_, [SAN2004]_) @@ -19,7 +20,7 @@ GROUP_STRUCTURES = {} .. _XMAS-172: https://www-nds.iaea.org/wimsd/energy.htm .. _SHEM-361: https://www.polymtl.ca/merlin/libraries.htm .. _activation: https://fispact.ukaea.uk/wiki/Keyword:GETXS -.. _CCFE-709: https://fispact.ukaea.uk/wiki/CCFE-709_group_structure +.. _CCFE-709: https://fispact.ukaea.uk/wiki/CCFE-709_group_structure .. _UKAEA-1102: https://fispact.ukaea.uk/wiki/UKAEA-1102_group_structure .. [SAR1990] Sartori, E., OECD/NEA Data Bank: Standard Energy Group Structures of Cross Section Libraries for Reactor Shielding, Reactor Cell and Fusion diff --git a/openmc/model/funcs.py b/openmc/model/funcs.py index 5f2db4d30..143033eca 100644 --- a/openmc/model/funcs.py +++ b/openmc/model/funcs.py @@ -3,6 +3,7 @@ from collections.abc import Iterable from math import sqrt from numbers import Real from functools import partial +from warnings import warn from openmc import XPlane, YPlane, Plane, ZCylinder, Quadric from openmc.checkvalue import check_type, check_value @@ -101,10 +102,14 @@ def borated_water(boron_ppm, temperature=293., pressure=0.1013, temp_unit='K', return out -def get_rectangular_prism(width, height, axis='z', origin=(0., 0.), - boundary_type='transmission', corner_radius=0.): +def rectangular_prism(width, height, axis='z', origin=(0., 0.), + boundary_type='transmission', corner_radius=0.): """Get an infinite rectangular prism from four planar surfaces. + .. versionchanged:: 0.11 + This function was renamed from `get_rectangular_prism` to + `rectangular_prism`. + Parameters ---------- width: float @@ -202,10 +207,21 @@ def get_rectangular_prism(width, height, axis='z', origin=(0., 0.), return prism -def get_hexagonal_prism(edge_length=1., orientation='y', origin=(0., 0.), - boundary_type='transmission', corner_radius=0.): +def get_rectangular_prism(*args, **kwargs): + warn("get_rectangular_prism(...) has been renamed rectangular_prism(...). " + "Future versions of OpenMC will not accept get_rectangular_prism.", + FutureWarning) + return rectangular_prism(*args, **kwargs) + + +def hexagonal_prism(edge_length=1., orientation='y', origin=(0., 0.), + boundary_type='transmission', corner_radius=0.): """Create a hexagon region from six surface planes. + .. versionchanged:: 0.11 + This function was renamed from `get_hexagonal_prism` to + `hexagonal_prism`. + Parameters ---------- edge_length : float @@ -346,6 +362,13 @@ def get_hexagonal_prism(edge_length=1., orientation='y', origin=(0., 0.), return prism +def get_hexagonal_prism(*args, **kwargs): + warn("get_hexagonal_prism(...) has been renamed hexagonal_prism(...). " + "Future versions of OpenMC will not accept get_hexagonal_prism.", + FutureWarning) + return hexagonal_prism(*args, **kwargs) + + def cylinder_from_points(p1, p2, r, **kwargs): """Return cylinder defined by two points passing through its center. diff --git a/src/geometry.cpp b/src/geometry.cpp index a92a4d34e..0f89a0ecd 100644 --- a/src/geometry.cpp +++ b/src/geometry.cpp @@ -173,32 +173,19 @@ find_cell_inner(Particle* p, const NeighborList* neighbor_list) //! Found a lower universe, update this coord level then search the next. // Set the lower coordinate level universe. - p->coord_[p->n_coord_].universe = c.fill_; + auto& coord {p->coord_[p->n_coord_]}; + coord.universe = c.fill_; // Set the position and direction. - p->coord_[p->n_coord_].r = p->coord_[p->n_coord_-1].r; - p->coord_[p->n_coord_].u = p->coord_[p->n_coord_-1].u; + coord.r = p->r_local(); + coord.u = p->u_local(); // Apply translation. - p->coord_[p->n_coord_].r -= c.translation_; + coord.r -= c.translation_; // Apply rotation. if (!c.rotation_.empty()) { - Position r = p->coord_[p->n_coord_].r; - p->coord_[p->n_coord_].r.x = r.x*c.rotation_[3] + r.y*c.rotation_[4] - + r.z*c.rotation_[5]; - p->coord_[p->n_coord_].r.y = r.x*c.rotation_[6] + r.y*c.rotation_[7] - + r.z*c.rotation_[8]; - p->coord_[p->n_coord_].r.z = r.x*c.rotation_[9] + r.y*c.rotation_[10] - + r.z*c.rotation_[11]; - Direction u = p->coord_[p->n_coord_].u; - p->coord_[p->n_coord_].u.x = u.x*c.rotation_[3] + u.y*c.rotation_[4] - + u.z*c.rotation_[5]; - p->coord_[p->n_coord_].u.y = u.x*c.rotation_[6] + u.y*c.rotation_[7] - + u.z*c.rotation_[8]; - p->coord_[p->n_coord_].u.z = u.x*c.rotation_[9] + u.y*c.rotation_[10] - + u.z*c.rotation_[11]; - p->coord_[p->n_coord_].rotated = true; + coord.rotate(c.rotation_); } // Update the coordinate level and recurse. @@ -211,26 +198,37 @@ find_cell_inner(Particle* p, const NeighborList* neighbor_list) Lattice& lat {*model::lattices[c.fill_]}; - // Determine lattice indices. - auto i_xyz = lat.get_indices(p->r_local(), p->u_local()); + // Set the position and direction. + auto& coord {p->coord_[p->n_coord_]}; + coord.r = p->r_local(); + coord.u = p->u_local(); - // Store lower level coordinates. - Position r = lat.get_local_position(p->r_local(), i_xyz); - p->coord_[p->n_coord_].r = r; - p->coord_[p->n_coord_].u = p->u_local(); + // Apply translation. + coord.r -= c.translation_; + + // Apply rotation. + if (!c.rotation_.empty()) { + coord.rotate(c.rotation_); + } + + // Determine lattice indices. + auto i_xyz = lat.get_indices(coord.r, coord.u); + + // Get local position in appropriate lattice cell + coord.r = lat.get_local_position(coord.r, i_xyz); // Set lattice indices. - p->coord_[p->n_coord_].lattice = c.fill_; - p->coord_[p->n_coord_].lattice_x = i_xyz[0]; - p->coord_[p->n_coord_].lattice_y = i_xyz[1]; - p->coord_[p->n_coord_].lattice_z = i_xyz[2]; + coord.lattice = c.fill_; + coord.lattice_x = i_xyz[0]; + coord.lattice_y = i_xyz[1]; + coord.lattice_z = i_xyz[2]; // Set the lower coordinate level universe. if (lat.are_valid_indices(i_xyz)) { - p->coord_[p->n_coord_].universe = lat[i_xyz]; + coord.universe = lat[i_xyz]; } else { if (lat.outer_ != NO_OUTER_UNIVERSE) { - p->coord_[p->n_coord_].universe = lat.outer_; + coord.universe = lat.outer_; } else { std::stringstream err_msg; err_msg << "Particle " << p->id_ << " is outside lattice " @@ -297,27 +295,32 @@ find_cell(Particle* p, bool use_neighbor_lists) void cross_lattice(Particle* p, const BoundaryInfo& boundary) { - auto& lat {*model::lattices[p->coord_[p->n_coord_-1].lattice]}; + auto& coord {p->coord_[p->n_coord_ - 1]}; + auto& lat {*model::lattices[coord.lattice]}; if (settings::verbosity >= 10 || simulation::trace) { std::stringstream msg; msg << " Crossing lattice " << lat.id_ << ". Current position (" - << p->coord_[p->n_coord_-1].lattice_x << "," - << p->coord_[p->n_coord_-1].lattice_y << "," - << p->coord_[p->n_coord_-1].lattice_z << ")"; + << coord.lattice_x << "," << coord.lattice_y << "," + << coord.lattice_z << "). r=" << p->r(); write_message(msg, 1); } // Set the lattice indices. - p->coord_[p->n_coord_-1].lattice_x += boundary.lattice_translation[0]; - p->coord_[p->n_coord_-1].lattice_y += boundary.lattice_translation[1]; - p->coord_[p->n_coord_-1].lattice_z += boundary.lattice_translation[2]; - std::array i_xyz {p->coord_[p->n_coord_-1].lattice_x, - p->coord_[p->n_coord_-1].lattice_y, - p->coord_[p->n_coord_-1].lattice_z}; + coord.lattice_x += boundary.lattice_translation[0]; + coord.lattice_y += boundary.lattice_translation[1]; + coord.lattice_z += boundary.lattice_translation[2]; + std::array i_xyz {coord.lattice_x, coord.lattice_y, coord.lattice_z}; // Set the new coordinate position. - p->r_local() = lat.get_local_position(p->coord_[p->n_coord_-2].r, i_xyz); + const auto& upper_coord {p->coord_[p->n_coord_ - 2]}; + const auto& cell {model::cells[upper_coord.cell]}; + Position r = upper_coord.r; + r -= cell->translation_; + if (!cell->rotation_.empty()) { + r = r.rotate(cell->rotation_); + } + p->r_local() = lat.get_local_position(r, i_xyz); if (!lat.are_valid_indices(i_xyz)) { // The particle is outside the lattice. Search for it from the base coords. @@ -362,9 +365,10 @@ BoundaryInfo distance_to_boundary(Particle* p) // Loop over each coordinate level. for (int i = 0; i < p->n_coord_; i++) { - Position r {p->coord_[i].r}; - Direction u {p->coord_[i].u}; - Cell& c {*model::cells[p->coord_[i].cell]}; + const auto& coord {p->coord_[i]}; + Position r {coord.r}; + Direction u {coord.u}; + Cell& c {*model::cells[coord.cell]}; // Find the oncoming surface in this cell and the distance to it. auto surface_distance = c.distance(r, u, p->surface_); @@ -372,10 +376,9 @@ BoundaryInfo distance_to_boundary(Particle* p) level_surf_cross = surface_distance.second; // Find the distance to the next lattice tile crossing. - if (p->coord_[i].lattice != C_NONE) { - auto& lat {*model::lattices[p->coord_[i].lattice]}; - std::array i_xyz {p->coord_[i].lattice_x, p->coord_[i].lattice_y, - p->coord_[i].lattice_z}; + if (coord.lattice != C_NONE) { + auto& lat {*model::lattices[coord.lattice]}; + std::array i_xyz {coord.lattice_x, coord.lattice_y, coord.lattice_z}; //TODO: refactor so both lattice use the same position argument (which //also means the lat.type attribute can be removed) std::pair> lattice_distance; @@ -384,8 +387,13 @@ BoundaryInfo distance_to_boundary(Particle* p) lattice_distance = lat.distance(r, u, i_xyz); break; case LatticeType::hex: - Position r_hex {p->coord_[i-1].r.x, p->coord_[i-1].r.y, - p->coord_[i].r.z}; + auto& cell_above {model::cells[p->coord_[i-1].cell]}; + Position r_hex {p->coord_[i-1].r}; + r_hex -= cell_above->translation_; + if (coord.rotated) { + r_hex = r_hex.rotate(cell_above->rotation_); + } + r_hex.z = coord.r.z; lattice_distance = lat.distance(r_hex, u, i_xyz); break; } diff --git a/src/particle.cpp b/src/particle.cpp index afcc56317..1fb1eadb0 100644 --- a/src/particle.cpp +++ b/src/particle.cpp @@ -33,6 +33,14 @@ namespace openmc { // LocalCoord implementation //============================================================================== +void +LocalCoord::rotate(const std::vector& rotation) +{ + this->r = this->r.rotate(rotation); + this->u = this->u.rotate(rotation); + this->rotated = true; +} + void LocalCoord::reset() { @@ -41,6 +49,7 @@ LocalCoord::reset() lattice = C_NONE; lattice_x = 0; lattice_y = 0; + lattice_z = 0; rotated = false; } @@ -337,9 +346,7 @@ Particle::transport() // If next level is rotated, apply rotation matrix const auto& m {model::cells[coord_[j].cell]->rotation_}; const auto& u {coord_[j].u}; - coord_[j + 1].u.x = m[3]*u.x + m[4]*u.y + m[5]*u.z; - coord_[j + 1].u.y = m[6]*u.x + m[7]*u.y + m[8]*u.z; - coord_[j + 1].u.z = m[9]*u.x + m[10]*u.y + m[11]*u.z; + coord_[j + 1].u = u.rotate(m); } else { // Otherwise, copy this level's direction coord_[j+1].u = coord_[j].u; diff --git a/src/position.cpp b/src/position.cpp index 8ab97c96d..1215bfb12 100644 --- a/src/position.cpp +++ b/src/position.cpp @@ -84,4 +84,21 @@ Position::operator-() const return {-x, -y, -z}; } +Position +Position::rotate(const std::vector& rotation) const +{ + return { + x*rotation[3] + y*rotation[4] + z*rotation[5], + x*rotation[6] + y*rotation[7] + z*rotation[8], + x*rotation[9] + y*rotation[10] + z*rotation[11] + }; +} + +std::ostream& +operator<<(std::ostream& os, Position r) +{ + os << "(" << r.x << ", " << r.y << ", " << r.z << ")"; + return os; +} + } // namespace openmc diff --git a/tests/regression_tests/lattice_hex_coincident/test.py b/tests/regression_tests/lattice_hex_coincident/test.py index eaea74b19..30bc470f2 100644 --- a/tests/regression_tests/lattice_hex_coincident/test.py +++ b/tests/regression_tests/lattice_hex_coincident/test.py @@ -103,10 +103,10 @@ class HexLatticeCoincidentTestHarness(PyAPITestHarness): inf_mat_univ = openmc.Universe(cells=[inf_mat,]) # a hex surface for the core to go inside of - hexprism = openmc.model.get_hexagonal_prism(edge_length=edge_length, - origin=(0.0, 0.0), - boundary_type = 'reflective', - orientation='x') + hexprism = openmc.model.hexagonal_prism(edge_length=edge_length, + origin=(0.0, 0.0), + boundary_type = 'reflective', + orientation='x') pincell_only_lattice = openmc.HexLattice(name="regular fuel assembly") pincell_only_lattice.center = (0., 0.) diff --git a/tests/regression_tests/lattice_hex_x/test.py b/tests/regression_tests/lattice_hex_x/test.py index 538f790da..ecb27c9c2 100644 --- a/tests/regression_tests/lattice_hex_x/test.py +++ b/tests/regression_tests/lattice_hex_x/test.py @@ -139,10 +139,10 @@ class HexLatticeOXTestHarness(PyAPITestHarness): # a hex surface for the core to go inside of - hexprism = openmc.model.get_hexagonal_prism(edge_length=edge_length, - origin=(0.0, 0.0), - boundary_type='reflective', - orientation='x') + hexprism = openmc.model.hexagonal_prism(edge_length=edge_length, + origin=(0.0, 0.0), + boundary_type='reflective', + orientation='x') region = hexprism & +fuel_bottom & -fuel_top inf_mat = openmc.Cell(cell_id=12) diff --git a/tests/regression_tests/lattice_mixed/geometry.xml b/tests/regression_tests/lattice_mixed/geometry.xml deleted file mode 100644 index 41aa05d41..000000000 --- a/tests/regression_tests/lattice_mixed/geometry.xml +++ /dev/null @@ -1,156 +0,0 @@ - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - 011 - 011 011 - 011 011 011 - 011 011 011 011 - 011 011 011 011 011 - 011 011 011 011 011 011 - 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 055 011 011 011 011 011 - 011 011 011 055 011 011 055 011 011 011 - 011 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 011 - 011 011 055 011 011 055 011 055 011 011 - 011 011 011 011 055 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 055 011 011 011 - 011 011 055 011 011 051 011 011 055 011 011 - 011 011 011 055 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 055 011 011 011 011 - 011 011 055 011 055 011 011 055 011 011 - 011 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 011 - 011 011 011 055 011 011 055 011 011 011 - 011 011 011 011 011 055 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 011 - 011 011 011 011 011 011 011 - 011 011 011 011 011 011 - 011 011 011 011 011 - 011 011 011 011 - 011 011 011 - 011 011 - 011 - - - - - - - - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 102 102 102 102 102 102 - 102 102 102 102 - - - - - - - - - - - - - - - - - - - diff --git a/tests/regression_tests/lattice_mixed/materials.xml b/tests/regression_tests/lattice_mixed/materials.xml deleted file mode 100644 index c7649fcf9..000000000 --- a/tests/regression_tests/lattice_mixed/materials.xml +++ /dev/null @@ -1,42 +0,0 @@ - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - diff --git a/tests/regression_tests/lattice_mixed/plots.xml b/tests/regression_tests/lattice_mixed/plots.xml deleted file mode 100644 index 903d69aef..000000000 --- a/tests/regression_tests/lattice_mixed/plots.xml +++ /dev/null @@ -1,32 +0,0 @@ - - - - - xy_cell - 0 0 0 - 30 30 - 500 500 - - - - xy_material - 0 0 0 - 30 30 - 500 500 - - - - yz_cell - 0 0 0 - 50 400 - 500 4000 - - - - yz_material - 0 0 0 - 5 5 - 500 500 - - - diff --git a/tests/regression_tests/lattice_mixed/results_true.dat b/tests/regression_tests/lattice_mixed/results_true.dat deleted file mode 100644 index 09af26768..000000000 --- a/tests/regression_tests/lattice_mixed/results_true.dat +++ /dev/null @@ -1,2 +0,0 @@ -k-combined: -9.877395E-01 1.744626E-02 diff --git a/tests/regression_tests/lattice_mixed/settings.xml b/tests/regression_tests/lattice_mixed/settings.xml deleted file mode 100644 index 69d50204e..000000000 --- a/tests/regression_tests/lattice_mixed/settings.xml +++ /dev/null @@ -1,14 +0,0 @@ - - - eigenvalue - 10 - 5 - 500 - - - - -8.0 -8.0 -176 - 8.0 8.0 176 - - - diff --git a/tests/regression_tests/lattice_mixed/test.py b/tests/regression_tests/lattice_mixed/test.py deleted file mode 100644 index 3df1d7a4f..000000000 --- a/tests/regression_tests/lattice_mixed/test.py +++ /dev/null @@ -1,6 +0,0 @@ -from tests.testing_harness import TestHarness - - -def test_lattice_mixed(): - harness = TestHarness('statepoint.10.h5') - harness.main() diff --git a/tests/regression_tests/lattice_mixed/__init__.py b/tests/regression_tests/lattice_rotated/__init__.py similarity index 100% rename from tests/regression_tests/lattice_mixed/__init__.py rename to tests/regression_tests/lattice_rotated/__init__.py diff --git a/tests/regression_tests/lattice_rotated/inputs_true.dat b/tests/regression_tests/lattice_rotated/inputs_true.dat new file mode 100644 index 000000000..8ee34b4f3 --- /dev/null +++ b/tests/regression_tests/lattice_rotated/inputs_true.dat @@ -0,0 +1,71 @@ + + + + + + + + + + + + 1.25 + 30 +
0.0 0.0
+ + 2 + 1 1 +1 2 1 + 1 1 +1 2 1 + 1 1 +1 1 1 + 1 1 + 1 + + + 1.25 1.25 + 30 + 4 4 + -2.5 -2.5 + +2 2 2 2 +1 1 1 1 +1 1 1 1 +1 1 1 1 + + + + + + + + + + + + + + + + + + + + + + + + + + + eigenvalue + 1000 + 5 + 0 + + + 0.0 0.0 0.0 + + + diff --git a/tests/regression_tests/lattice_rotated/results_true.dat b/tests/regression_tests/lattice_rotated/results_true.dat new file mode 100644 index 000000000..b92298cc2 --- /dev/null +++ b/tests/regression_tests/lattice_rotated/results_true.dat @@ -0,0 +1,2 @@ +k-combined: +4.538090E-01 2.452457E-02 diff --git a/tests/regression_tests/lattice_rotated/test.py b/tests/regression_tests/lattice_rotated/test.py new file mode 100644 index 000000000..9cab92e1e --- /dev/null +++ b/tests/regression_tests/lattice_rotated/test.py @@ -0,0 +1,85 @@ +import numpy as np +import openmc + +from tests.testing_harness import PyAPITestHarness + + +def rotated_lattice_model(): + model = openmc.model.Model() + + # Create some materials + fuel1 = openmc.Material() + fuel1.set_density('g/cm3', 10.0) + fuel1.add_nuclide('U235', 1.0) + fuel2 = openmc.Material() + fuel2.set_density('g/cm3', 10.0) + fuel2.add_nuclide('U238', 1.0) + water = openmc.Material() + water.set_density('g/cm3', 1.0) + water.add_nuclide('H1', 2.0) + water.add_nuclide('O16', 1.0) + water.add_s_alpha_beta('c_H_in_H2O') + model.materials.extend([fuel1, fuel2, water]) + + # Create universes for lattices + r_pin = openmc.ZCylinder(r=0.25) + fuel_cell = openmc.Cell(fill=fuel1, region=-r_pin) + water_cell = openmc.Cell(fill=water, region=+r_pin) + pin_universe = openmc.Universe(cells=(fuel_cell, water_cell)) + r_big_pin = openmc.ZCylinder(r=0.5) + fuel2_cell = openmc.Cell(fill=fuel2, region=-r_big_pin) + water2_cell = openmc.Cell(fill=water, region=+r_big_pin) + big_pin_universe = openmc.Universe(cells=(fuel2_cell, water2_cell)) + all_water_cell = openmc.Cell(fill=water) + outer_universe = openmc.Universe(30, cells=(all_water_cell,)) + + # Create hexagonal lattice + pitch = 1.25 + hexlat = openmc.HexLattice() + hexlat.center = (0., 0.) + hexlat.pitch = [pitch] + hexlat.outer = outer_universe + outer_ring = [big_pin_universe] + [pin_universe]*11 + middle_ring = [big_pin_universe] + [pin_universe]*5 + inner_ring = [big_pin_universe] + hexlat.universes = [outer_ring, middle_ring, inner_ring] + + # Create rectangular lattice + rectlat = openmc.RectLattice() + rectlat.center = (0., 0.) + rectlat.pitch = (pitch, pitch) + rectlat.lower_left = (-2*pitch, -2*pitch) + rectlat.outer = outer_universe + rectlat.universes = np.full((4, 4), pin_universe) + rectlat.universes[0] = big_pin_universe + + # Create cell filled with translated/rotated rectangular lattice on left + left_cyl = openmc.ZCylinder(x0=-4.0, r=4.0) + left_cell = openmc.Cell(fill=rectlat, region=-left_cyl) + left_cell.translation = (-4.0, 0.0, 0.0) + left_cell.rotation = (0.0, 0.0, 45.0) + + # Create cell filled with translated/rotated hexagonal lattice on right + right_cyl = openmc.ZCylinder(x0=4.0, r=4.0) + right_cell = openmc.Cell(fill=hexlat, region=-right_cyl) + right_cell.translation = (4.0, 0.0, 0.0) + right_cell.rotation = (0.0, 0.0, 30.0) + + # Finish up with the geometry + outer_cyl = openmc.ZCylinder(r=8.0, boundary_type='vacuum') + main_cell = openmc.Cell(fill=water, region=-outer_cyl & +left_cyl & +right_cyl) + model.geometry = openmc.Geometry([main_cell, left_cell, right_cell]) + + model.settings.batches = 5 + model.settings.inactive = 0 + model.settings.particles = 1000 + model.settings.source = openmc.Source(space=openmc.stats.Point()) + model.settings.export_to_xml() + + return model + + +def test(): + model = rotated_lattice_model() + harness = PyAPITestHarness('statepoint.5.h5', model) + harness.main() diff --git a/tests/unit_tests/test_filters.py b/tests/unit_tests/test_filters.py index 5817c4c3d..676ee452e 100644 --- a/tests/unit_tests/test_filters.py +++ b/tests/unit_tests/test_filters.py @@ -11,7 +11,7 @@ def box_model(): m.add_nuclide('U235', 1.0) m.set_density('g/cm3', 1.0) - box = openmc.model.get_rectangular_prism(10., 10., boundary_type='vacuum') + box = openmc.model.rectangular_prism(10., 10., boundary_type='vacuum') c = openmc.Cell(fill=m, region=box) model.geometry.root_universe = openmc.Universe(cells=[c]) diff --git a/tests/unit_tests/test_geometry.py b/tests/unit_tests/test_geometry.py index f3718630a..69211f72e 100644 --- a/tests/unit_tests/test_geometry.py +++ b/tests/unit_tests/test_geometry.py @@ -204,9 +204,9 @@ def test_get_by_name(): def test_hex_prism(): - hex_prism = openmc.model.get_hexagonal_prism(edge_length=5.0, - origin=(0.0, 0.0), - orientation='y') + hex_prism = openmc.model.hexagonal_prism(edge_length=5.0, + origin=(0.0, 0.0), + orientation='y') # clear checks assert (0.0, 0.0, 0.0) in hex_prism assert (10.0, 10.0, 10.0) not in hex_prism @@ -214,10 +214,10 @@ def test_hex_prism(): assert (0.0, 5.01, 0.0) not in hex_prism assert (0.0, 4.99, 0.0) in hex_prism - rounded_hex_prism = openmc.model.get_hexagonal_prism(edge_length=5.0, - origin=(0.0, 0.0), - orientation='y', - corner_radius=1.0) + rounded_hex_prism = openmc.model.hexagonal_prism(edge_length=5.0, + origin=(0.0, 0.0), + orientation='y', + corner_radius=1.0) # clear checks assert (0.0, 0.0, 0.0) in rounded_hex_prism