Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://mooseframework.inl.gov
       3             : //*
       4             : //* All rights reserved, see COPYRIGHT for full restrictions
       5             : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
       6             : //*
       7             : //* Licensed under LGPL 2.1, please see LICENSE for details
       8             : //* https://www.gnu.org/licenses/lgpl-2.1.html
       9             : 
      10             : #include "ConeRayStudy.h"
      11             : 
      12             : // Local includes
      13             : #include "RayTracingAngularQuadrature.h"
      14             : 
      15             : registerMooseObject("RayTracingApp", ConeRayStudy);
      16             : 
      17             : InputParameters
      18          90 : ConeRayStudy::validParams()
      19             : {
      20          90 :   auto params = RepeatableRayStudyBase::validParams();
      21             : 
      22          90 :   params.addClassDescription(
      23             :       "Ray study that spawns Rays in the direction of a cone from a given set of starting points.");
      24             : 
      25         180 :   params.addRequiredParam<std::vector<Point>>("start_points", "The point(s) of the cone(s).");
      26         180 :   params.addRequiredParam<std::vector<Point>>(
      27             :       "directions", "The direction(s) of the cone(s) (points down the center of each cone).");
      28         180 :   params.addParam<std::vector<Real>>("scaling_factors",
      29             :                                      "Scaling factors for each cone (if any). Defaults to 1.");
      30             : 
      31         180 :   params.addRequiredParam<std::vector<Real>>("half_cone_angles",
      32             :                                              "Angle of the half-cones in degrees (must be <= 90)");
      33         180 :   params.addParam<std::vector<unsigned int>>(
      34             :       "polar_quad_orders",
      35             :       "Order of the polar quadrature for each cone. Polar angle is between ray and the "
      36             :       "cone direction. Must be positive. This will default to 2 for all cones if not provided.");
      37         180 :   params.addParam<std::vector<unsigned int>>(
      38             :       "azimuthal_quad_orders",
      39             :       "Order of the azimuthal quadrature per quadrant for each cone. The azimuthal angle is "
      40             :       "measured in a plane perpendicular to the cone direction. This will default to 30 if not "
      41             :       "provided. Must be positive.");
      42             : 
      43         180 :   params.addRequiredParam<std::string>("ray_data_name",
      44             :                                        "The name of the Ray data that the angular quadrature "
      45             :                                        "weights and factors will be filled into for "
      46             :                                        "properly weighting the line source per Ray.");
      47             : 
      48             :   // Here we will not use Ray registration because we create so many Rays that
      49             :   // we don't want to keep track of them by name
      50          90 :   params.set<bool>("_use_ray_registration") = false;
      51             : 
      52          90 :   return params;
      53           0 : }
      54             : 
      55          52 : ConeRayStudy::ConeRayStudy(const InputParameters & parameters)
      56             :   : RepeatableRayStudyBase(parameters),
      57          52 :     _start_points(getParam<std::vector<Point>>("start_points")),
      58         104 :     _directions(getParam<std::vector<Point>>("directions")),
      59         108 :     _scaling_factors(isParamValid("scaling_factors")
      60          52 :                          ? getParam<std::vector<Real>>("scaling_factors")
      61             :                          : std::vector<Real>(_start_points.size(), 1)), // default to 1
      62         104 :     _half_cone_angles(getParam<std::vector<Real>>("half_cone_angles")),
      63         108 :     _polar_quad_orders(isParamValid("polar_quad_orders")
      64          52 :                            ? getParam<std::vector<unsigned int>>("polar_quad_orders")
      65          52 :                            : std::vector<unsigned int>(_start_points.size(), 2)), // default to 2
      66          56 :     _azimuthal_quad_orders(
      67          52 :         isParamValid("azimuthal_quad_orders")
      68          52 :             ? getParam<std::vector<unsigned int>>("azimuthal_quad_orders")
      69          52 :             : std::vector<unsigned int>(_start_points.size(), 30)), // default to 30
      70         156 :     _ray_data_index(registerRayData(getParam<std::string>("ray_data_name")))
      71             : {
      72          52 :   if (_directions.size() != _start_points.size())
      73           2 :     paramError("directions", "Not the same size as start_points.");
      74             : 
      75          50 :   if (_scaling_factors.size() != _start_points.size())
      76           2 :     paramError("scaling_factors", "Not the same size as start_points.");
      77             : 
      78          48 :   if (_half_cone_angles.size() != _start_points.size())
      79           2 :     paramError("half_cone_angles", "Not the same size as start_points.");
      80          90 :   for (const auto val : _half_cone_angles)
      81          46 :     if (val <= 0 || val > 90)
      82           2 :       paramError("half_cone_angles", "Must be > 0 and <= 90 degrees");
      83             : 
      84          44 :   if (_polar_quad_orders.size() != _start_points.size())
      85           2 :     paramError("polar_quad_orders", "Not the same size as start_points.");
      86             : 
      87          42 :   if (_azimuthal_quad_orders.size() != _start_points.size())
      88           2 :     paramError("azimuthal_quad_orders", "Not the same size as start_points.");
      89             : 
      90          40 :   if (_mesh.dimension() == 1)
      91           2 :     mooseError("Does not support 1D.");
      92          38 : }
      93             : 
      94             : void
      95          28 : ConeRayStudy::defineRays()
      96             : {
      97             :   // Loop through each cone
      98          56 :   for (std::size_t i = 0; i < _start_points.size(); ++i)
      99             :   {
     100             :     // Setup the angular quadrature and rotate it about the direction
     101             :     // (the direction points down the middle of the cone)
     102          28 :     RayTracingAngularQuadrature aq(_mesh.dimension(),
     103             :                                    _polar_quad_orders[i],
     104             :                                    _azimuthal_quad_orders[i],
     105          28 :                                    std::cos(_half_cone_angles[i] * M_PI / 180.),
     106          28 :                                    1);
     107          28 :     aq.rotate(_directions[i].unit());
     108             : 
     109             :     // For all angles in the angular quadrature, spawn a Ray
     110        1078 :     for (std::size_t l = 0; l < aq.numDirections(); ++l)
     111             :     {
     112             :       // Get a Ray from the study to initialize
     113        1050 :       std::shared_ptr<Ray> ray = acquireReplicatedRay();
     114             : 
     115             :       // Start from our cone point in the rotated angular quadrature direction
     116             :       // Note here that we do not need to set the starting element - all Rays
     117             :       // at this point are replicated across all processors and will be
     118             :       // properly claimed (moved to the starting proc with the correct starting elem)
     119        1050 :       ray->setStart(_start_points[i]);
     120        1050 :       ray->setStartingDirection(aq.getDirection(l));
     121             : 
     122             :       // Add the angular quadrature weight and scaling factor as data on the Ray for weighting
     123             :       //
     124             :       // In the 2D case, the 3D directions projected into the 2D plane may overlap. Therefore,
     125             :       // we could have multiple weights/sins for a single direction, which is why here we grab
     126             :       // the total weight.
     127             :       //
     128             :       // The angular quadrature weights sum to 2pi, so scale such that they scale to 1
     129        1050 :       ray->data(_ray_data_index) = _scaling_factors[i] * aq.getTotalWeight(l) / (2. * M_PI);
     130             : 
     131             :       // Done with this Ray - move it to be traced later on
     132        1050 :       _rays.emplace_back(std::move(ray));
     133             :     }
     134          28 :   }
     135          28 : }