bluemira.radiation_transport.neutronics.make_pre_cell ===================================================== .. py:module:: bluemira.radiation_transport.neutronics.make_pre_cell .. autoapi-nested-parse:: Make pre-cells using bluemira wires. Attributes ---------- .. autoapisummary:: bluemira.radiation_transport.neutronics.make_pre_cell.CCW_90 bluemira.radiation_transport.neutronics.make_pre_cell.CW_90 Classes ------- .. autoapisummary:: bluemira.radiation_transport.neutronics.make_pre_cell.PreCell bluemira.radiation_transport.neutronics.make_pre_cell.PreCellArray bluemira.radiation_transport.neutronics.make_pre_cell.DivertorPreCell bluemira.radiation_transport.neutronics.make_pre_cell.DivertorPreCellArray Functions --------- .. autoapisummary:: bluemira.radiation_transport.neutronics.make_pre_cell.ratio_of_distances bluemira.radiation_transport.neutronics.make_pre_cell.find_equidistant_point bluemira.radiation_transport.neutronics.make_pre_cell.pick_higher_point bluemira.radiation_transport.neutronics.make_pre_cell.calculate_new_circle Module Contents --------------- .. py:data:: CCW_90 .. py:data:: CW_90 .. py:function:: ratio_of_distances(point_of_interest: numpy.typing.NDArray[numpy.float64], anchor1: numpy.typing.NDArray[numpy.float64], normal1: numpy.typing.NDArray[numpy.float64], anchor2: numpy.typing.NDArray[numpy.float64], normal2: numpy.typing.NDArray[numpy.float64]) -> numpy.typing.NDArray[numpy.float64] Find how close a point is to line 1 and line 2, and then express that ratio as a tuple of floats that sums up to unit. Each line is given by the user by specifying any point on that line, and a direction NORMAL to that line. The point_of_interest must lie on the positive side of the line. :param point_of_interest: point to which we want to calculate the ratio of distances. :param anchor1: Any point on line 1 and line 2 respectively. :param anchor2: Any point on line 1 and line 2 respectively. :param normal1: The positive distance direction of line 1 and line 2 respectively. :param normal2: The positive distance direction of line 1 and line 2 respectively. :returns: ratio of distances. Sum of these two numbers should yield unity (1.0). :rtype: dist_to_1, dist_to_2 :raises GeometryError: distance not on positive side of line .. py:function:: find_equidistant_point(point1: numpy.typing.NDArray[numpy.float64], point2: numpy.typing.NDArray[numpy.float64], distance: float) -> tuple[numpy.typing.NDArray[numpy.float64], numpy.typing.NDArray[numpy.float64]] Find the two (or 0) points on a 2D plane that are equidistant to each other. :param point1: 2D points, each with shape (2,) :param point2: 2D points, each with shape (2,) :param distance: the distance that both points must obey by. :returns: The two intersection points of circle1 and circle2. :rtype: intersection1, intersection2 :raises GeometryError: Points too close .. py:function:: pick_higher_point(point1: numpy.typing.NDArray[numpy.float64], point2: numpy.typing.NDArray[numpy.float64], vector: numpy.typing.NDArray[numpy.float64]) -> numpy.typing.NDArray[numpy.float64] Pick the point that, when projected onto `vector`, gives a higher value. :param point1: Point to choose from. :param point2: Point to choose from. :param vector: The co-vector direction that we want to project onto. :returns: The point further among the two choices [point1, point2]. :rtype: point .. py:function:: calculate_new_circle(old_circle_info: bluemira.radiation_transport.neutronics.wires.CircleInfo, new_points: numpy.typing.NDArray) -> bluemira.radiation_transport.neutronics.wires.CircleInfo Calculate how far does the new circle get shifted. :param old_circle_info: an object accessed by WireInfoList[i].key_points. info on circle where the start_point and end_point are each of shape (3,). :param new_points: array of shape (2, 3) :returns: An instance of CircleInfo representing the new (scaled) arc of circle. :rtype: new_circle_info .. py:class:: PreCell(interior_wire: bluemira.geometry.wire.BluemiraWire | bluemira.geometry.coordinates.Coordinates, vv_wire: bluemira.geometry.wire.BluemiraWire, exterior_wire: bluemira.geometry.wire.BluemiraWire) A pre-cell is the BluemiraWire outlining the reactor cross-section BEFORE they have been simplified into straight-lines. Unlike a Cell, its outline may be constructed from curved lines. .. py:attribute:: interior_wire .. py:attribute:: vv_wire .. py:attribute:: exterior_wire .. py:attribute:: vertex .. py:attribute:: vv_point .. py:attribute:: outline .. py:property:: half_solid :type: bluemira.geometry.solid.BluemiraSolid Create the 180° revolved shape on demand only. Revolved 180° instead of 360° for easier viewing .. py:property:: blanket_outline :type: bluemira.geometry.wire.BluemiraWire Create the outline of the blanket, i.e. the part excluding the vacuum vessel. .. py:property:: blanket_half_solid :type: bluemira.geometry.solid.BluemiraSolid Get the volume of the blanket .. py:method:: plot_2d(*args, **kwargs) -> matplotlib.axes.Axes Plot the outline in 2D :returns: Axes on which the pre-cell is plotted. .. py:method:: show_cad(*args, **kwargs) -> None Plot the outline in 3D .. py:property:: cell_walls :type: bluemira.radiation_transport.neutronics.radial_wall.CellWalls The side (clockwise side and counter-clockwise) walls of this cell. Only create it when called, because some instances of PreCell will never use it. it is of type :class:`~bluemira.radiation_transport.neutronics.radial_wall.CellWalls`. .. py:property:: normal_to_interior :type: numpy.typing.NDArray The vector pointing from the interior_wire direction towards the exterior_wire, specifically, it's perpendicular to the interior_wire. Also only created when called, because it's not needed .. py:method:: get_cell_wall_cut_points_by_fraction(fraction: float) -> numpy.typing.NDArray Find the cut points on the cell's side walls by multiplying the original lengths by a fraction. When fraction=0, this returns the interior_start and interior_end. :param fraction: A scalar value :returns: The position of the pre-cell wall end points at the required fraction, array of shape (2, 2) [[cw_wall x, cw_wall z], [ccw_wall x, ccw_wall z]]. :rtype: new end points .. py:method:: get_cell_wall_cut_points_by_thickness(thickness: float) Offset a line parallel to the interior_wire towards the exterior direction. Then, find where this line intersect the cell's side walls. :param thickness: A scalar value :returns: The position of the pre-cell wall end points at the required thickness, array of shape (2, 2) [[cw_wall x, cw_wall z], [ccw_wall x, ccw_wall z]]. :rtype: new end points .. py:class:: PreCellArray(list_of_pre_cells: list[PreCell]) A list of pre-cells materials :param list_of_pre_cells: An adjacent list of pre-cells :raises GeometryError: Precells must share corners .. rubric:: Notes The list of pre-cells must be adjacent to each other. .. py:attribute:: pre_cells .. py:attribute:: cell_walls .. py:method:: check_convexity() Check that the outermost points of self forms a convex hull. :raises GeometryError: PreCellArray is not convex .. py:property:: volumes :type: tuple[float] Create the iterable of volumes on demand. .. py:method:: straighten_exterior(*, preserve_volume: bool = False) -> PreCellArray Turn the exterior curves of each cell into a straight edge. This is done at the PreCellArray level instead of the PreCell level to allow volume preservation, see Parameters below for more details. :param preserve_volume: Whether to preserve the volume of each cell during the transformation from pre-cell with curved-edge to pre-cell with straight edges. If True, increase the length of the cut lines appropriately to compensate for the volume loss due to the straight line approximation. :returns: An array of the copies of the pre-cells, where the exterior wall and interior walls (as well as the dividing walls between adjacent pre-cells) were made of straight BluemiraWire's, as opposed to the curved lines initally provided. .. py:method:: plot_2d(*args, **kwargs) Plot pre cells in 2d :returns: Axes on which the pre-cell array is plotted. .. py:method:: show_cad(*args, **kwargs) -> None Show pre cell CAD .. py:method:: exterior_vertices() -> numpy.typing.NDArray Returns all of the vertices on the exterior side of the pre-cell array. :returns: array of shape (N+1, 3) arranged clockwise (inboard to outboard). :rtype: exterior_vertices .. py:method:: interior_vertices() -> numpy.typing.NDArray Returns all of the vertices on the interior side of the pre-cell array. :returns: array of shape (N+1, 3) arranged clockwise (inboard to outboard). :rtype: interior_vertices .. py:method:: __len__() -> int Number of pre cells .. py:method:: __getitem__(index_or_slice) -> list[PreCell] | PreCell Get pre cell .. py:method:: __setitem__(index_or_slice, new_pre_cell: PreCell | PreCellArray) Set an element to be a new Precell, or a slice to be a new PreCellarray. .. py:method:: __add__(other_array) -> PreCellArray Adding two list together to create a new one. :returns: Both pre-cell arrays (self and other_array) concatenated together to form a new PreCellArray. :raises TypeError: Operation not possible between types .. py:method:: __repr__() -> str String representation :returns: A string that includes the number of pre-cells in the array. .. py:method:: copy() Each element of the new_copy.pre_cells list points to the same items as the self.pre_cells :returns: copy of itself :rtype: new_copy .. py:class:: DivertorPreCell(interior_wire: bluemira.radiation_transport.neutronics.wires.WireInfoList, vv_wire: bluemira.radiation_transport.neutronics.wires.WireInfoList, exterior_wire: bluemira.radiation_transport.neutronics.wires.WireInfoList) An intermediate class between the bluemira wire and the final csg product. A divertor pre-cell is the equivalent of a blanket's pre-cell, but for the divertor. .. py:attribute:: interior_wire .. py:attribute:: vv_wire .. py:attribute:: exterior_wire .. py:attribute:: cw_wall .. py:attribute:: ccw_wall .. py:attribute:: vertex .. py:method:: plot_2d(*args, **kwargs) -> matplotlib.axes.Axes Plot 2d precell :returns: The Axes object on which the divertor pre-cell was plotted. .. py:method:: show_cad(*args, **kwargs) -> None Show precell CAD .. py:property:: outline :type: bluemira.geometry.wire.BluemiraWire We don't need the volume value, so we're only going to generate the outline when the user wants to plot it. .. py:property:: half_solid :type: bluemira.geometry.solid.BluemiraSolid Create the 180° revolved shape on demand only. Revolved 180° instead of 360° for easier viewing .. py:method:: offset_interior_wire(thickness: float) -> bluemira.radiation_transport.neutronics.wires.WireInfoList Offset the interior wire towards the exterior_wire. The true problem of expanding/shrinking a wire is a much more difficult one, so I've only opted for a simpler (but incorrect) approach of pushing the wire to a desired direction determined by how close it is to the wall. TODO: Re-write this method, as it currently do weird things to circles. New method should do the following: 1. Find a point to be our new origin. (Name that point "p") This is likely to be a point near exterior_vertices. 2. Scale down the interior_wire by x%. We should give more thought of how to derive/search for the optimal vector (p[0], p[-1], x), such that all lines are displaced by This proposed method has several new benefits: 1. The circles will be scaled correctly (center moves towards the new origin by x%, radius scaled down by x%). 2. All tangents are preserved, so no need to change them. :returns: The interior wire offsetted by thickness, represented as a WireInfoList. .. py:class:: DivertorPreCellArray(list_of_div_pc: list[DivertorPreCell]) An array of Divertor pre-cells :raises GeometryError: Divertor precells must share corners .. py:attribute:: pre_cells .. py:method:: check_convexity() Check that the outermost points of self forms a convex hull. :raises GeometryError: Precell array is not convex .. py:method:: exterior_vertices() -> numpy.typing.NDArray Returns all of the tokamak's poloidal cross-section's outside corners' coordinates, in 3D. :returns: aray of shape (N+1, 3) arranged counter-clockwise (inboard to outboard). :rtype: exterior_vertices .. py:method:: interior_vertices() -> numpy.typing.NDArray Returns all of the tokamak's poloidal cross-section's inside corners' coordinates, in 3D. :returns: A numpy array of every vertex in each pre-cell's vertices. .. py:method:: __len__() -> int Number of pre cells .. py:method:: __getitem__(index_or_slice) -> list[DivertorPreCell] | DivertorPreCell Get pre cell .. py:method:: __repr__() -> str String representation .. py:method:: plot_2d(*args, **kwargs) Plot precell array cad in 2d :returns: The Axes on which the divertor pre-cell array was plotted. .. py:method:: show_cad(*args, **kwargs) -> None Show precell array CAD .. py:method:: copy() NOT a deepcopy, each element of the new_copy.pre_cells list points to the same items as the self.pre_cells :returns: copy of itself