Line data    Source code

       1             : /********************************************************************/
       2             : /*                  SOFTWARE COPYRIGHT NOTIFICATION                 */
       3             : /*                             Cardinal                             */
       4             : /*                                                                  */
       5             : /*                  (c) 2021 UChicago Argonne, LLC                  */
       6             : /*                        ALL RIGHTS RESERVED                       */
       7             : /*                                                                  */
       8             : /*                 Prepared by UChicago Argonne, LLC                */
       9             : /*               Under Contract No. DE-AC02-06CH11357               */
      10             : /*                With the U. S. Department of Energy               */
      11             : /*                                                                  */
      12             : /*             Prepared by Battelle Energy Alliance, LLC            */
      13             : /*               Under Contract No. DE-AC07-05ID14517               */
      14             : /*                With the U. S. Department of Energy               */
      15             : /*                                                                  */
      16             : /*                 See LICENSE for full restrictions                */
      17             : /********************************************************************/
      18             : 
      19             : #ifdef ENABLE_NEK_COUPLING
      20             : 
      21             : #include "NekPointValue.h"
      22             : #include "pointInterpolation.hpp"
      23             : 
      24             : registerMooseObject("CardinalApp", NekPointValue);
      25             : 
      26             : InputParameters
      27         495 : NekPointValue::validParams()
      28             : {
      29         495 :   InputParameters params = NekFieldPostprocessor::validParams();
      30         990 :   params.addRequiredParam<Point>("point", "The physical point where the field will be evaluated");
      31         495 :   params.addClassDescription("Uses NekRS's pointInterpolation to query the NekRS solution at a "
      32             :                              "point (does not need to be a grid point).");
      33         495 :   return params;
      34           0 : }
      35             : 
      36         160 : NekPointValue::NekPointValue(const InputParameters & parameters)
      37         319 :   : NekFieldPostprocessor(parameters), _point(getParam<Point>("point")), _value(0)
      38             : {
      39         159 : }
      40             : 
      41             : void
      42         151 : NekPointValue::execute()
      43             : {
      44         151 :   std::vector<dfloat> x = {_point(0) / nekrs::referenceLength()};
      45         151 :   std::vector<dfloat> y = {_point(1) / nekrs::referenceLength()};
      46         151 :   std::vector<dfloat> z = {_point(2) / nekrs::referenceLength()};
      47         151 :   int n = x.size();
      48             : 
      49         151 :   auto nrs = nekrs::nrsPtr();
      50             : 
      51             :   // set the points to be interpolated; we include this every time we call the
      52             :   // interpolator, in case the mesh is moving. TODO: auto-detect for efficiency
      53         151 :   auto interp = pointInterpolation_t(nrs->meshV, platform->comm.mpiComm());
      54         151 :   interp.setPoints(x, y, z);
      55             :   const auto verbosity = pointInterpolation_t::VerbosityLevel::Basic;
      56         151 :   interp.find(verbosity);
      57             : 
      58             :   // if this slot is not used for coupling, we are responsible for copying it to
      59             :   // device before calling the interp function. Otherwise, Cardinal handles copying
      60             :   // to device automatically.
      61         151 :   if (_field == field::usrwrk00 && !_nek_problem->isUsrWrkSlotReservedForCoupling(0))
      62           1 :     _nek_problem->copyIndividualScratchSlot(0);
      63         151 :   if (_field == field::usrwrk01 && !_nek_problem->isUsrWrkSlotReservedForCoupling(1))
      64           1 :     _nek_problem->copyIndividualScratchSlot(1);
      65         151 :   if (_field == field::usrwrk02 && !_nek_problem->isUsrWrkSlotReservedForCoupling(2))
      66           1 :     _nek_problem->copyIndividualScratchSlot(2);
      67             : 
      68             :   // interpolate the field onto those points
      69         151 :   occa::memory o_interpolated;
      70             :   int n_values = n;
      71         151 :   switch (_field)
      72             :   {
      73          64 :     case field::velocity_component:
      74             :     case field::velocity_x:
      75             :     case field::velocity_y:
      76             :     case field::velocity_z:
      77             :     case field::velocity:
      78             :     case field::velocity_x_squared:
      79             :     case field::velocity_y_squared:
      80             :     case field::velocity_z_squared:
      81          64 :       n_values = n * 3; // nVFields
      82          64 :       o_interpolated = platform->device.malloc<dfloat>(n_values);
      83          64 :       interp.eval(n_values, nrs->fieldOffset, nrs->scalar->o_U, n, o_interpolated);
      84             :       break;
      85           8 :     case field::pressure:
      86           8 :       o_interpolated = platform->device.malloc<dfloat>(n);
      87           8 :       interp.eval(1, nrs->fieldOffset, nrs->fluid->o_P, n, o_interpolated);
      88             :       break;
      89          44 :     case field::temperature:
      90             :     case field::scalar01:
      91             :     case field::scalar02:
      92             :     case field::scalar03:
      93          44 :       n_values = n * nrs->Nscalar;
      94          44 :       o_interpolated = platform->device.malloc<dfloat>(n_values);
      95          44 :       interp.eval(n_values, nekrs::scalarFieldOffset(), nrs->scalar->o_S, n, o_interpolated);
      96             :       break;
      97             :     case field::unity:
      98             :       break;
      99          27 :     case field::usrwrk00:
     100             :     case field::usrwrk01:
     101             :     case field::usrwrk02:
     102          27 :       n_values = n * _nek_problem->nUsrWrkSlots();
     103          27 :       o_interpolated = platform->device.malloc<dfloat>(n_values);
     104          27 :       interp.eval(n_values, nekrs::fieldOffset(), nrs->bc.o_usrwrk, n, o_interpolated);
     105             :       break;
     106           0 :     default:
     107           0 :       mooseError("Unhandled NekFieldEnum in NekPointValue!");
     108             :   }
     109             : 
     110             :   // the interpolation happens on device, so we need to copy it back to the host
     111         151 :   std::vector<dfloat> interpolated(n_values);
     112         151 :   o_interpolated.copyTo(interpolated.data(), n_values);
     113             : 
     114             :   // because of NekRS's way of storing solution, we need extra steps to actually
     115             :   // return what the user wants
     116         151 :   switch (_field)
     117             :   {
     118             :     case field::velocity_component:
     119           8 :       _value = interpolated[0] * _velocity_direction(0) + interpolated[1] * _velocity_direction(1) +
     120           8 :                interpolated[2] * _velocity_direction(2);
     121           8 :       break;
     122             :     case field::velocity_x:
     123           8 :       _value = interpolated[0];
     124           8 :       break;
     125             :     case field::velocity_y:
     126           8 :       _value = interpolated[1];
     127           8 :       break;
     128             :     case field::velocity_z:
     129           8 :       _value = interpolated[2];
     130           8 :       break;
     131             :     case field::velocity:
     132           8 :       _value = std::sqrt(interpolated[0] * interpolated[0] + interpolated[1] * interpolated[1] +
     133           8 :                          interpolated[2] * interpolated[2]);
     134           8 :       break;
     135             :     case field::velocity_x_squared:
     136           8 :       _value = interpolated[0] * interpolated[0];
     137           8 :       break;
     138             :     case field::velocity_y_squared:
     139           8 :       _value = interpolated[1] * interpolated[1];
     140           8 :       break;
     141             :     case field::velocity_z_squared:
     142           8 :       _value = interpolated[2] * interpolated[2];
     143           8 :       break;
     144             :     case field::pressure:
     145           8 :       _value = interpolated[0];
     146           8 :       break;
     147             :     case field::temperature:
     148           8 :       _value = interpolated[0];
     149           8 :       break;
     150          12 :     case field::scalar01:
     151          12 :       _value = interpolated[nekrs::scalarSlot(1)];
     152          12 :       break;
     153          12 :     case field::scalar02:
     154          12 :       _value = interpolated[nekrs::scalarSlot(2)];
     155          12 :       break;
     156          12 :     case field::scalar03:
     157          12 :       _value = interpolated[nekrs::scalarSlot(3)];
     158          12 :       break;
     159           8 :     case field::unity:
     160           8 :       _value = 1;
     161           8 :       break;
     162             :     case field::usrwrk00:
     163           9 :       _value = interpolated[0];
     164           9 :       break;
     165             :     case field::usrwrk01:
     166           9 :       _value = interpolated[1];
     167           9 :       break;
     168             :     case field::usrwrk02:
     169           9 :       _value = interpolated[2];
     170           9 :       break;
     171           0 :     default:
     172           0 :       mooseError("Unhandled NekFieldEnum in NekPointValue!");
     173             :   }
     174             : 
     175         151 :   _value = _value * nekrs::nondimensionalDivisor(_field) + nekrs::nondimensionalAdditive(_field);
     176         148 : }
     177             : 
     178             : Real
     179         148 : NekPointValue::getValue() const
     180             : {
     181         148 :   return _value;
     182             : }
     183             : 
     184             : #endif