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 "VectorEMRobinBC.h"
      11             : #include "ElectromagneticEnums.h"
      12             : #include "ElectromagneticConstants.h"
      13             : #include "Function.h"
      14             : #include <complex>
      15             : 
      16             : registerMooseObject("ElectromagneticsApp", VectorEMRobinBC);
      17             : 
      18             : InputParameters
      19         258 : VectorEMRobinBC::validParams()
      20             : {
      21         258 :   InputParameters params = VectorIntegratedBC::validParams();
      22         258 :   params.addClassDescription("First order Robin-style Absorbing/Port BC for vector variables.");
      23         516 :   params.addParam<FunctionName>("beta", 1.0, "Scalar wave number.");
      24         516 :   MooseEnum component("real imaginary");
      25         516 :   params.addParam<MooseEnum>(
      26             :       "component", component, "Variable field component (real or imaginary).");
      27         516 :   params.addRequiredCoupledVar("coupled_field", "Coupled field variable.");
      28         516 :   params.addParam<FunctionName>("real_incoming", 0.0, "Real incoming field vector.");
      29         516 :   params.addParam<FunctionName>("imag_incoming", 0.0, "Imaginary incoming field vector.");
      30         516 :   MooseEnum mode("absorbing port", "port");
      31         516 :   params.addParam<MooseEnum>("mode", mode, "Mode of operation for VectorEMRobinBC.");
      32         258 :   return params;
      33         258 : }
      34             : 
      35         136 : VectorEMRobinBC::VectorEMRobinBC(const InputParameters & parameters)
      36             :   : VectorIntegratedBC(parameters),
      37             : 
      38         136 :     _beta(getFunction("beta")),
      39             : 
      40         272 :     _component(getParam<MooseEnum>("component")),
      41             : 
      42         136 :     _coupled_val(coupledVectorValue("coupled_field")),
      43         136 :     _coupled_var_num(coupled("coupled_field")),
      44             : 
      45         136 :     _inc_real(getFunction("real_incoming")),
      46         136 :     _inc_imag(getFunction("imag_incoming")),
      47             : 
      48         272 :     _mode(getParam<MooseEnum>("mode")),
      49             : 
      50         136 :     _real_incoming_was_set(parameters.isParamSetByUser("real_incoming")),
      51         272 :     _imag_incoming_was_set(parameters.isParamSetByUser("imag_incoming"))
      52             : {
      53         136 :   if (_mode == EM::ABSORBING && (_real_incoming_was_set || _imag_incoming_was_set))
      54             :   {
      55           2 :     mooseError(
      56             :         "In ",
      57           2 :         _name,
      58             :         ", mode was set to Absorbing, while an incoming profile function (used for Port BCs) was "
      59             :         "defined. Either remove the profile function parameters, or set your BC to Port mode!");
      60             :   }
      61         134 : }
      62             : 
      63             : Real
      64     8650256 : VectorEMRobinBC::computeQpResidual()
      65             : {
      66             :   // Initialize E components
      67     8650256 :   std::complex<double> field_0(0, 0);
      68     8650256 :   std::complex<double> field_1(0, 0);
      69     8650256 :   std::complex<double> field_2(0, 0);
      70             : 
      71             :   // Create E and ncrossE for residual based on component parameter
      72     8650256 :   if (_component == EM::REAL)
      73             :   {
      74     4325128 :     field_0.real(_u[_qp](0));
      75     4325128 :     field_0.imag(_coupled_val[_qp](0));
      76             : 
      77     4325128 :     field_1.real(_u[_qp](1));
      78     4325128 :     field_1.imag(_coupled_val[_qp](1));
      79             : 
      80     4325128 :     field_2.real(_u[_qp](2));
      81     4325128 :     field_2.imag(_coupled_val[_qp](2));
      82             :   }
      83             :   else
      84             :   {
      85     4325128 :     field_0.real(_coupled_val[_qp](0));
      86     4325128 :     field_0.imag(_u[_qp](0));
      87             : 
      88     4325128 :     field_1.real(_coupled_val[_qp](1));
      89     4325128 :     field_1.imag(_u[_qp](1));
      90             : 
      91     4325128 :     field_2.real(_coupled_val[_qp](2));
      92     4325128 :     field_2.imag(_u[_qp](2));
      93             :   }
      94             :   VectorValue<std::complex<double>> field(field_0, field_1, field_2);
      95             : 
      96             :   // Initialize proper zero defaults for the case when real and imag incoming functions are not set
      97             :   VectorValue<std::complex<double>> field_inc(std::complex<double>(0.0, 0.0),
      98             :                                               std::complex<double>(0.0, 0.0),
      99             :                                               std::complex<double>(0.0, 0.0));
     100             :   VectorValue<std::complex<double>> curl_inc(std::complex<double>(0.0, 0.0),
     101             :                                              std::complex<double>(0.0, 0.0),
     102             :                                              std::complex<double>(0.0, 0.0));
     103             : 
     104             :   /**
     105             :    * Stops us calling un-implemented methods within the scalar ConstantFunction
     106             :    * default for _inc_real and _inc_imag function variables, since we don't have
     107             :    * an automatic vector function constant default in MOOSE FEProblemBase right now.
     108             :    */
     109     8650256 :   if (_real_incoming_was_set && _imag_incoming_was_set)
     110             :   {
     111             :     // Creating vector, curl for field_inc before residual and Jacobian contributions
     112     8637440 :     RealVectorValue inc_real_value = _inc_real.vectorValue(_t, _q_point[_qp]);
     113     8637440 :     RealVectorValue inc_imag_value = _inc_imag.vectorValue(_t, _q_point[_qp]);
     114     8637440 :     RealVectorValue inc_real_curl = _inc_real.curl(_t, _q_point[_qp]);
     115     8637440 :     RealVectorValue inc_imag_curl = _inc_imag.curl(_t, _q_point[_qp]);
     116             : 
     117     8637440 :     std::complex<double> field_inc_0(inc_real_value(0), inc_imag_value(0));
     118     8637440 :     std::complex<double> field_inc_1(inc_real_value(1), inc_imag_value(1));
     119     8637440 :     std::complex<double> field_inc_2(inc_real_value(2), inc_imag_value(2));
     120     8637440 :     field_inc = VectorValue<std::complex<double>>(field_inc_0, field_inc_1, field_inc_2);
     121             : 
     122     8637440 :     std::complex<double> curl_inc_0(inc_real_curl(0), inc_imag_curl(0));
     123     8637440 :     std::complex<double> curl_inc_1(inc_real_curl(1), inc_imag_curl(1));
     124     8637440 :     std::complex<double> curl_inc_2(inc_real_curl(2), inc_imag_curl(2));
     125     8637440 :     curl_inc = VectorValue<std::complex<double>>(curl_inc_0, curl_inc_1, curl_inc_2);
     126             :   }
     127             : 
     128             :   // Do some error checking
     129             :   mooseAssert(_beta.value(_t, _q_point[_qp]) > 0,
     130             :               "Wave number expected to be positive, calculated to be "
     131             :                   << _beta.value(_t, _q_point[_qp]));
     132             : 
     133             :   std::complex<double> u_inc_dot_test = 0.0;
     134     8650256 :   switch (_mode)
     135             :   {
     136     8637440 :     case EM::PORT:
     137             :       // Calculate incoming wave contribution to BC residual
     138     8637440 :       u_inc_dot_test = _test[_i][_qp].cross(_normals[_qp]) * curl_inc +
     139     8637440 :                        EM::j * _beta.value(_t, _q_point[_qp]) *
     140     8637440 :                            (_test[_i][_qp].cross(_normals[_qp]) * _normals[_qp].cross(field_inc));
     141     8637440 :       break;
     142             :     case EM::ABSORBING:
     143             :       break;
     144             :   }
     145             : 
     146             :   // Calculate solution field contribution to BC residual (first order version)
     147     8650256 :   std::complex<double> p_dot_test = EM::j * _beta.value(_t, _q_point[_qp]) *
     148     8650256 :                                     _test[_i][_qp].cross(_normals[_qp]) *
     149     8650256 :                                     _normals[_qp].cross(field);
     150             : 
     151             :   std::complex<double> diff = u_inc_dot_test - p_dot_test;
     152             :   Real res = 0.0;
     153     8650256 :   switch (_component)
     154             :   {
     155             :     case EM::REAL:
     156             :       res = diff.real();
     157     4325128 :       break;
     158             :     case EM::IMAGINARY:
     159             :       res = diff.imag();
     160     4325128 :       break;
     161             :   }
     162     8650256 :   return res;
     163             : }
     164             : 
     165             : Real
     166       70424 : VectorEMRobinBC::computeQpJacobian()
     167             : {
     168       70424 :   return 0.0;
     169             : }
     170             : 
     171             : Real
     172       70424 : VectorEMRobinBC::computeQpOffDiagJacobian(unsigned int jvar)
     173             : {
     174       70424 :   Real off_diag_jac = _beta.value(_t, _q_point[_qp]) * _test[_i][_qp].cross(_normals[_qp]) *
     175       70424 :                       _normals[_qp].cross(_phi[_j][_qp]);
     176             : 
     177       70424 :   if (_component == EM::REAL && jvar == _coupled_var_num)
     178             :     return off_diag_jac;
     179       35212 :   else if (_component == EM::IMAGINARY && jvar == _coupled_var_num)
     180       35212 :     return -off_diag_jac;
     181             :   else
     182             :     return 0.0;
     183             : }