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 "CrossIC.h"
      11             : 
      12             : // Portions of this code Copyright 2007-2009 Roy Stogner
      13             : //
      14             : // Permission is hereby granted, free of charge, to any person obtaining
      15             : // a copy of this software and associated documentation files (the
      16             : //"Software"), to deal in the Software without restriction, including
      17             : // without limitation the rights to use, copy, modify, merge, publish,
      18             : // distribute, sublicense, and/or sell copies of the Software, and to
      19             : // permit persons to whom the Software is furnished to do so, subject to
      20             : // the following conditions:
      21             : //
      22             : // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
      23             : // EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
      24             : // MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
      25             : // IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
      26             : // CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
      27             : // TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
      28             : // SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
      29             : 
      30             : registerMooseObject("PhaseFieldApp", CrossIC);
      31             : 
      32             : InputParameters
      33         507 : CrossIC::validParams()
      34             : {
      35         507 :   InputParameters params = C1ICBase::validParams();
      36         507 :   params.addClassDescription("Cross-shaped initial condition");
      37        1014 :   params.addParam<Real>("x1", 0.0, "The x coordinate of the lower left-hand corner of the box");
      38        1014 :   params.addParam<Real>("y1", 0.0, "The y coordinate of the lower left-hand corner of the box");
      39        1014 :   params.addParam<Real>("x2", 1.0, "The x coordinate of the upper right-hand corner of the box");
      40        1014 :   params.addParam<Real>("y2", 1.0, "The y coordinate of the upper right-hand corner of the box");
      41         507 :   return params;
      42           0 : }
      43             : 
      44         266 : CrossIC::CrossIC(const InputParameters & parameters)
      45             :   : C1ICBase(parameters),
      46         266 :     _x1(parameters.get<Real>("x1")),
      47         266 :     _y1(parameters.get<Real>("y1")),
      48         266 :     _x2(parameters.get<Real>("x2")),
      49         532 :     _y2(parameters.get<Real>("y2"))
      50             : {
      51         266 : }
      52             : 
      53             : Real
      54      937642 : CrossIC::value(const Point & p)
      55             : {
      56      937642 :   const Real x = (p(0) - _x1) / (_x2 - _x1);
      57      937642 :   const Real y = (p(1) - _y1) / (_y2 - _y1);
      58             : 
      59      937642 :   const Real cmax = _average + _amplitude;
      60      937642 :   const Real cmin = _average - _amplitude;
      61             : 
      62      937642 :   if (x > (0.5 + _length / 2.0 + _buffer + _interface))
      63             :   {
      64             :     return cmin;
      65             :   }
      66      850362 :   else if (x > 0.5 + _length / 2.0)
      67             :   {
      68       52612 :     if (y > 0.5 + _width / 2.0 + _buffer + _interface)
      69             :     {
      70             :       return cmin;
      71             :     }
      72       32120 :     else if (y > 0.5 + _width / 2.0)
      73             :     {
      74        2992 :       Real xd = x - .5 - _length / 2.0;
      75        2992 :       Real yd = y - .5 - _width / 2.0;
      76        2992 :       Real r = std::sqrt(xd * xd + yd * yd);
      77        2992 :       return interfaceValue(r);
      78             :     }
      79       29128 :     else if (y > 0.5 - _width / 2.0)
      80             :     {
      81        5644 :       return interfaceValue(x - .5 - _length / 2.0);
      82             :     }
      83       23484 :     else if (y > 0.5 - _width / 2. - _buffer - _interface)
      84             :     {
      85        2824 :       Real xd = x - .5 - _length / 2.0;
      86        2824 :       Real yd = y - .5 + _width / 2.0;
      87        2824 :       Real r = std::sqrt(xd * xd + yd * yd);
      88        2824 :       return interfaceValue(r);
      89             :     }
      90             :     else
      91             :     {
      92             :       return cmin;
      93             :     }
      94             :   }
      95      797750 :   else if (x > 0.5 + _width / 2. + 2. * _buffer + _interface)
      96             :   {
      97      158970 :     if (y > 0.5 + _width / 2 + _buffer + _interface)
      98             :       return cmin;
      99       97164 :     else if (y > 0.5 + _width / 2.0)
     100        9240 :       return interfaceValue(y - .5 - _width / 2.0);
     101       87924 :     else if (y > 0.5 - _width / 2.0)
     102             :       return cmax;
     103       71046 :     else if (y > 0.5 - _width / 2. - _buffer - _interface)
     104        8532 :       return interfaceValue(.5 - _width / 2. - y);
     105             :     else
     106             :       return cmin;
     107             :   }
     108      638780 :   else if (x > 0.5 + _width / 2.0)
     109             :   {
     110       76828 :     if (y > 0.5 + _length / 2. + _buffer + _interface)
     111             :     {
     112             :       return cmin;
     113             :     }
     114       68772 :     else if (y > 0.5 + _length / 2.0)
     115             :     {
     116        4840 :       Real xd = x - (.5 + _width / 2.0);
     117        4840 :       Real yd = y - (.5 + _length / 2.0);
     118        4840 :       Real r = std::sqrt(xd * xd + yd * yd);
     119        4840 :       return interfaceValue(r);
     120             :     }
     121       63932 :     else if (y > 0.5 + _width / 2. + 2. * _buffer + _interface)
     122             :     {
     123       14404 :       return interfaceValue(x - .5 - _width / 2.0);
     124             :     }
     125       49528 :     else if (y > 0.5 + _width / 2.0)
     126             :     {
     127        7272 :       Real xd = x - (.5 + _width / 2. + 2. * _buffer + _interface);
     128        7272 :       Real yd = y - (.5 + _width / 2. + 2. * _buffer + _interface);
     129        7272 :       Real r = std::sqrt(xd * xd + yd * yd);
     130        7272 :       return interfaceValue(2. * _buffer + _interface - r);
     131             :     }
     132       42256 :     else if (y > 0.5 - _width / 2.0)
     133             :     {
     134             :       return cmax;
     135             :     }
     136       34572 :     else if (y > 0.5 - _width / 2. - 2. * _buffer - _interface)
     137             :     {
     138        6504 :       Real xd = x - (.5 + _width / 2. + 2. * _buffer + _interface);
     139        6504 :       Real yd = y - (.5 - _width / 2. - 2. * _buffer - _interface);
     140        6504 :       Real r = std::sqrt(xd * xd + yd * yd);
     141        6504 :       return interfaceValue(2. * _buffer + _interface - r);
     142             :     }
     143       28068 :     else if (y > 0.5 - _length / 2.0)
     144             :     {
     145       14860 :       return interfaceValue(x - .5 - _width / 2.0);
     146             :     }
     147       13208 :     else if (y > 0.5 - _length / 2. - _buffer - _interface)
     148             :     {
     149        5152 :       Real xd = x - (.5 + _width / 2.0);
     150        5152 :       Real yd = y - (.5 - _length / 2.0);
     151        5152 :       Real r = std::sqrt(xd * xd + yd * yd);
     152        5152 :       return interfaceValue(r);
     153             :     }
     154             :     else
     155             :     {
     156             :       return cmin;
     157             :     }
     158             :   }
     159      561952 :   else if (x > 0.5 - _width / 2.0)
     160             :   {
     161       96352 :     if (y > 0.5 + _length / 2 + _buffer + _interface)
     162             :       return cmin;
     163       88332 :     else if (y > 0.5 + _length / 2.0)
     164        5644 :       return interfaceValue(y - .5 - _length / 2.0);
     165       82688 :     else if (y > 0.5 - _length / 2.0)
     166             :       return cmax;
     167       16367 :     else if (y > 0.5 - _length / 2. - _buffer - _interface)
     168        7348 :       return interfaceValue(.5 - _length / 2. - y);
     169             :     else
     170             :       return cmin;
     171             :   }
     172      465600 :   else if (x > 0.5 - _width / 2. - 2. * _buffer - _interface)
     173             :   {
     174       91688 :     if (y > 0.5 + _length / 2. + _buffer + _interface)
     175             :     {
     176             :       return cmin;
     177             :     }
     178       84432 :     else if (y > 0.5 + _length / 2.0)
     179             :     {
     180        4520 :       Real xd = x - (.5 - _width / 2.0);
     181        4520 :       Real yd = y - (.5 + _length / 2.0);
     182        4520 :       Real r = std::sqrt(xd * xd + yd * yd);
     183        4520 :       return interfaceValue(r);
     184             :     }
     185       79912 :     else if (y > 0.5 + _width / 2. + 2. * _buffer + _interface)
     186             :     {
     187       13700 :       return interfaceValue(.5 - _width / 2. - x);
     188             :     }
     189       66212 :     else if (y > 0.5 + _width / 2.0)
     190             :     {
     191        6696 :       Real xd = x - (.5 - _width / 2. - 2. * _buffer - _interface);
     192        6696 :       Real yd = y - (.5 + _width / 2. + 2. * _buffer + _interface);
     193        6696 :       Real r = std::sqrt(xd * xd + yd * yd);
     194        6696 :       return interfaceValue(2. * _buffer + _interface - r);
     195             :     }
     196       59516 :     else if (y > 0.5 - _width / 2.0)
     197             :     {
     198             :       return cmax;
     199             :     }
     200       50154 :     else if (y > 0.5 - _width / 2. - 2. * _buffer - _interface)
     201             :     {
     202       10116 :       Real xd = x - (.5 - _width / 2. - 2. * _buffer - _interface);
     203       10116 :       Real yd = y - (.5 - _width / 2. - 2. * _buffer - _interface);
     204       10116 :       Real r = std::sqrt(xd * xd + yd * yd);
     205       10116 :       return interfaceValue(2. * _buffer + _interface - r);
     206             :     }
     207       40038 :     else if (y > 0.5 - _length / 2.0)
     208             :     {
     209       23288 :       return interfaceValue(.5 - _width / 2. - x);
     210             :     }
     211       16750 :     else if (y > 0.5 - _length / 2. - _buffer - _interface)
     212             :     {
     213        7496 :       Real xd = x - (.5 - _width / 2.0);
     214        7496 :       Real yd = y - (.5 - _length / 2.0);
     215        7496 :       Real r = std::sqrt(xd * xd + yd * yd);
     216        7496 :       return interfaceValue(r);
     217             :     }
     218             :     else
     219             :     {
     220             :       return cmin;
     221             :     }
     222             :   }
     223      373912 :   else if (x > 0.5 - _length / 2.0)
     224             :   {
     225      206810 :     if (y > 0.5 + _width / 2 + _buffer + _interface)
     226             :       return cmin;
     227      144448 :     else if (y > 0.5 + _width / 2.0)
     228        9320 :       return interfaceValue(y - .5 - _width / 2.0);
     229      135128 :     else if (y > 0.5 - _width / 2.0)
     230             :       return cmax;
     231      113640 :     else if (y > 0.5 - _width / 2. - _buffer - _interface)
     232       14900 :       return interfaceValue(.5 - _width / 2. - y);
     233             :     else
     234             :       return cmin;
     235             :   }
     236      167102 :   else if (x > (.5 - _length / 2. - _buffer - _interface))
     237             :   {
     238       69832 :     if (y > 0.5 + _width / 2. + _buffer + _interface)
     239             :     {
     240             :       return cmin;
     241             :     }
     242       47816 :     else if (y > 0.5 + _width / 2.0)
     243             :     {
     244        3232 :       Real xd = x - (.5 - _length / 2.0);
     245        3232 :       Real yd = y - .5 - _width / 2.0;
     246        3232 :       Real r = std::sqrt(xd * xd + yd * yd);
     247        3232 :       return interfaceValue(r);
     248             :     }
     249       44584 :     else if (y > 0.5 - _width / 2.0)
     250             :     {
     251        7348 :       return interfaceValue(.5 - _length / 2. - x);
     252             :     }
     253       37236 :     else if (y > 0.5 - _width / 2. - _buffer - _interface)
     254             :     {
     255        4768 :       Real xd = x - (.5 - _length / 2.0);
     256        4768 :       Real yd = y - .5 + _width / 2.0;
     257        4768 :       Real r = std::sqrt(xd * xd + yd * yd);
     258        4768 :       return interfaceValue(r);
     259             :     }
     260             :     else
     261             :       return cmin;
     262             :   }
     263             :   else
     264             :     return cmin;
     265             : }
     266             : 
     267             : RealGradient
     268      110352 : CrossIC::gradient(const Point & p)
     269             : {
     270      110352 :   Point pxminus = p, pxplus = p, pyminus = p, pyplus = p;
     271             : 
     272      110352 :   pxminus(0) -= TOLERANCE;
     273      110352 :   pyminus(1) -= TOLERANCE;
     274      110352 :   pxplus(0) += TOLERANCE;
     275      110352 :   pyplus(1) += TOLERANCE;
     276             : 
     277      110352 :   Number uxminus = value(pxminus), uxplus = value(pxplus), uyminus = value(pyminus),
     278      110352 :          uyplus = value(pyplus);
     279             : 
     280      110352 :   return Gradient((uxplus - uxminus) / 2.0 / TOLERANCE, (uyplus - uyminus) / 2.0 / TOLERANCE);
     281             : }