Line data    Source code

       1             : 
       2             : /**********************************************************************/
       3             : /*                    DO NOT MODIFY THIS HEADER                       */
       4             : /*             Swift, a Fourier spectral solver for MOOSE             */
       5             : /*                                                                    */
       6             : /*            Copyright 2024 Battelle Energy Alliance, LLC            */
       7             : /*                        ALL RIGHTS RESERVED                         */
       8             : /**********************************************************************/
       9             : 
      10             : #include "LBMDirichletWallBC.h"
      11             : 
      12             : using namespace torch::indexing;
      13             : 
      14             : registerMooseObject("SwiftApp", LBMDirichletWallBC);
      15             : 
      16             : InputParameters
      17           0 : LBMDirichletWallBC::validParams()
      18             : {
      19           0 :   InputParameters params = LBMBoundaryCondition::validParams();
      20           0 :   params.addRequiredParam<TensorInputBufferName>("f_old", "Old state distribution function");
      21           0 :   params.addClassDescription("LBMDirichletWallBC object");
      22           0 :   params.addRequiredParam<TensorInputBufferName>("velocity", "Fluid velocity");
      23           0 :   params.addParam<std::string>("value",
      24             :                                "0.0"
      25             :                                "Value at the boundary");
      26           0 :   return params;
      27           0 : }
      28             : 
      29           0 : LBMDirichletWallBC::LBMDirichletWallBC(const InputParameters & parameters)
      30             :   : LBMBoundaryCondition(parameters),
      31           0 :     _f_old(_lb_problem.getBufferOld(getParam<TensorInputBufferName>("f_old"), 1)),
      32           0 :     _velocity(getInputBuffer("velocity")),
      33           0 :     _value(_lb_problem.getConstant<Real>(getParam<std::string>("value")))
      34             : {
      35           0 :   computeBoundaryNormals();
      36           0 : }
      37             : 
      38             : void
      39           0 : LBMDirichletWallBC::computeBoundaryNormals()
      40             : {
      41           0 :   if (_lb_problem.isBinaryMedia())
      42             :   {
      43             :     const torch::Tensor & binary_media = _lb_problem.getBinaryMedia();
      44           0 :     _binary_mesh = binary_media.clone();
      45             : 
      46             :     /* boolean travelling tensor will indicate which directions at every boundary lattice need to be
      47             :      * computed. This is done by rolling binary media around in each direction and finding where the
      48             :      * zeros are. */
      49             : 
      50             :     torch::Tensor boolean_travelling_tensor =
      51           0 :         torch::ones(_shape_q, MooseTensor::intTensorOptions());
      52             : 
      53           0 :     for (int64_t ic = 1; ic < _stencil._q; ic++)
      54             :     {
      55           0 :       int64_t ex = _stencil._ex[ic].item<int64_t>();
      56           0 :       int64_t ey = _stencil._ey[ic].item<int64_t>();
      57           0 :       int64_t ez = _stencil._ez[ic].item<int64_t>();
      58           0 :       torch::Tensor shifted_mesh = torch::roll(binary_media, {ex, ey, ez}, {0, 1, 2});
      59           0 :       torch::Tensor adjacent_to_boundary = (shifted_mesh == 0) & (binary_media == 1);
      60           0 :       auto boolean_travelling_tensor_slice = torch::ones_like(adjacent_to_boundary);
      61           0 :       boolean_travelling_tensor_slice.masked_fill_(adjacent_to_boundary, 0);
      62           0 :       boolean_travelling_tensor.index_put_({Slice(), Slice(), Slice(), ic},
      63             :                                            boolean_travelling_tensor_slice);
      64           0 :       _binary_mesh.masked_fill_(adjacent_to_boundary, 2);
      65             :     }
      66             : 
      67             :     /* normal vectors are computed by summing up microscopic velocity vectors e_i at boundaries
      68             :      * where these vectors are pointing away from the boundary */
      69             : 
      70             :     // boolean_travelling_tensor = 1 - boolean_travelling_tensor;
      71             :     // boolean_travelling_tensor.unsqueeze_(-1);
      72             :     // torch::Tensor e_xyz = torch::stack({_ex, _ey, _ez}).permute({1, 2, 3, 4, 0});
      73             : 
      74             :     // torch::Tensor normals = torch::einsum("abcqm,ijkqn->abcn", {boolean_travelling_tensor,
      75             :     // e_xyz})
      76             :     //                             .to(MooseTensor::floatTensorOptions());
      77             : 
      78             :     // _boundary_normals = normals / torch::norm(normals, 2, -1).unsqueeze(-1);
      79             :     // _boundary_normals = torch::where(
      80             :     //     torch::isnan(_boundary_normals), torch::zeros_like(_boundary_normals),
      81             :     //     _boundary_normals);
      82             : 
      83             :     // /* tangent vectors simply t_i = e_i - (e_i dot n) * n where n is unit normal vector to the
      84             :     //  * boundary */
      85             : 
      86             :     // _e_xyz = e_xyz.to(MooseTensor::floatTensorOptions());
      87             :     // _boundary_tangent_vectors =
      88             :     //     _e_xyz - torch::einsum("ijkqm,abcdm->abcq", {_e_xyz, _boundary_normals.unsqueeze(-2)})
      89             :     //                      .unsqueeze(-1) *
      90             :     //                  _boundary_normals.unsqueeze(-2);
      91             :   }
      92             :   else
      93           0 :     mooseError("Binary media must be avialble to use LBMDirichletWallBC boundary condition.");
      94           0 : }
      95             : 
      96             : void
      97           0 : LBMDirichletWallBC::wallBoundary()
      98             : {
      99             : 
     100           0 :   if (_lb_problem.getTotalSteps() == 0)
     101             :   {
     102           0 :     _boundary_mask = (_binary_mesh.unsqueeze(-1).expand_as(_u) == 2);
     103           0 :     _boundary_mask = _boundary_mask.to(torch::kBool);
     104             :   }
     105             : 
     106           0 :   torch::Tensor f_bounce_back = torch::ones_like(_u) * _w * _value;
     107           0 :   _u.index_put_({_boundary_mask}, f_bounce_back.index({_boundary_mask}));
     108           0 : }