Line data    Source code

       1             : // The libMesh Finite Element Library.
       2             : // Copyright (C) 2002-2025 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
       3             : 
       4             : // This library is free software; you can redistribute it and/or
       5             : // modify it under the terms of the GNU Lesser General Public
       6             : // License as published by the Free Software Foundation; either
       7             : // version 2.1 of the License, or (at your option) any later version.
       8             : 
       9             : // This library is distributed in the hope that it will be useful,
      10             : // but WITHOUT ANY WARRANTY; without even the implied warranty of
      11             : // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
      12             : // Lesser General Public License for more details.
      13             : 
      14             : // You should have received a copy of the GNU Lesser General Public
      15             : // License along with this library; if not, write to the Free Software
      16             : // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
      17             : 
      18             : 
      19             : 
      20             : // C++ includes
      21             : 
      22             : // Local includes
      23             : #include "libmesh/elem.h"
      24             : #include "libmesh/mesh_refinement.h"
      25             : #include "libmesh/remote_elem.h"
      26             : 
      27             : namespace libMesh
      28             : {
      29             : 
      30             : 
      31             : //--------------------------------------------------------------------
      32             : // Elem methods
      33             : 
      34             : /**
      35             :  * The following functions only apply when
      36             :  * AMR is enabled and thus are not present
      37             :  * otherwise.
      38             :  */
      39             : #ifdef LIBMESH_ENABLE_AMR
      40             : 
      41   311395842 : void Elem::set_p_level(unsigned int p)
      42             : {
      43             :   // Maintain the parent's p level as the minimum of it's children
      44   327601301 :   if (this->parent() != nullptr)
      45             :     {
      46    64413442 :       unsigned int parent_p_level = this->parent()->p_level();
      47             : 
      48             :       // If our new p level is less than our parents, our parents drops
      49    64413442 :       if (parent_p_level > p)
      50             :         {
      51         241 :           this->parent()->set_p_level(p);
      52             : 
      53             :           // And we should keep track of the drop, in case we need to
      54             :           // do a projection later.
      55          32 :           this->parent()->set_p_refinement_flag(Elem::JUST_COARSENED);
      56             :         }
      57             :       // If we are the lowest p level and it increases, so might
      58             :       // our parent's, but we have to check every other child to see
      59    64413201 :       else if (parent_p_level == _p_level && _p_level < p)
      60             :         {
      61       24167 :           _p_level = cast_int<unsigned char>(p);
      62       24167 :           parent_p_level = cast_int<unsigned char>(p);
      63      177007 :           for (auto & c : this->parent()->child_ref_range())
      64      152840 :             parent_p_level = std::min(parent_p_level,
      65      184916 :                                       c.p_level());
      66             : 
      67             :           // When its children all have a higher p level, the parent's
      68             :           // should rise
      69       27091 :           if (parent_p_level > this->parent()->p_level())
      70             :             {
      71        1825 :               this->parent()->set_p_level(parent_p_level);
      72             : 
      73             :               // And we should keep track of the rise, in case we need to
      74             :               // do a projection later.
      75         386 :               this->parent()->set_p_refinement_flag(Elem::JUST_REFINED);
      76             :             }
      77             : 
      78       24167 :           return;
      79             :         }
      80             :     }
      81             : 
      82    46387513 :   this->hack_p_level(p);
      83             : }
      84             : 
      85             : 
      86             : 
      87     1177529 : void Elem::refine (MeshRefinement & mesh_refinement)
      88             : {
      89       72858 :   libmesh_assert_equal_to (this->refinement_flag(), Elem::REFINE);
      90       72858 :   libmesh_assert (this->active());
      91             : 
      92     1177529 :   const unsigned int nc = this->n_children();
      93             : 
      94             :   // Create my children if necessary
      95     1177529 :   if (!_children)
      96             :     {
      97     2177412 :       _children = std::make_unique<Elem *[]>(nc);
      98             : 
      99      144562 :       unsigned int parent_p_level = this->p_level();
     100     1160987 :       const unsigned int nei = this->n_extra_integers();
     101     7384256 :       for (unsigned int c = 0; c != nc; c++)
     102             :         {
     103     6603575 :           auto current_child = Elem::build(this->type(), this);
     104     6603579 :           _children[c] = current_child.get();
     105             : 
     106      380306 :           current_child->set_refinement_flag(Elem::JUST_REFINED);
     107     6223269 :           current_child->set_p_level(parent_p_level);
     108      760616 :           current_child->set_p_refinement_flag(this->p_refinement_flag());
     109             : 
     110    51446778 :           for (auto cnode : current_child->node_index_range())
     111             :             {
     112             :               Node * node =
     113    47827533 :                 mesh_refinement.add_node(*this, c, cnode,
     114     2604024 :                                          current_child->processor_id());
     115    56148367 :               node->set_n_systems (this->n_systems());
     116    45223509 :               current_child->set_node(cnode, node);
     117             :             }
     118             : 
     119     6223269 :           Elem * added_child = mesh_refinement.add_elem (std::move(current_child));
     120     7911803 :           added_child->set_n_systems(this->n_systems());
     121      380306 :           libmesh_assert_equal_to (added_child->n_extra_integers(),
     122             :                                    this->n_extra_integers());
     123     6237111 :           for (unsigned int i=0; i != nei; ++i)
     124       12734 :             added_child->set_extra_integer(i, this->get_extra_integer(i));
     125     5462653 :         }
     126             :     }
     127             :   else
     128             :     {
     129        1156 :       unsigned int parent_p_level = this->p_level();
     130       82710 :       for (unsigned int c = 0; c != nc; c++)
     131             :         {
     132        2312 :           Elem * current_child = this->child_ptr(c);
     133       66168 :           if (current_child != remote_elem)
     134             :             {
     135        2240 :               libmesh_assert(current_child->subactive());
     136        2240 :               current_child->set_refinement_flag(Elem::JUST_REFINED);
     137       62121 :               current_child->set_p_level(parent_p_level);
     138        4480 :               current_child->set_p_refinement_flag(this->p_refinement_flag());
     139             :             }
     140             :         }
     141             :     }
     142             : 
     143             :   // Un-set my refinement flag now
     144       72858 :   this->set_refinement_flag(Elem::INACTIVE);
     145             : 
     146             :   // Leave the p refinement flag set - we will need that later to get
     147             :   // projection operations correct
     148             :   // this->set_p_refinement_flag(Elem::INACTIVE);
     149             : 
     150             : #ifndef NDEBUG
     151      455476 :   for (unsigned int c = 0; c != nc; c++)
     152      382618 :     if (this->child_ptr(c) != remote_elem)
     153             :       {
     154      382546 :         libmesh_assert_equal_to (this->child_ptr(c)->parent(), this);
     155      382546 :         libmesh_assert(this->child_ptr(c)->active());
     156             :       }
     157             : #endif
     158       72858 :   libmesh_assert (this->ancestor());
     159     1177529 : }
     160             : 
     161             : 
     162             : 
     163      345905 : void Elem::coarsen()
     164             : {
     165       17078 :   libmesh_assert_equal_to (this->refinement_flag(), Elem::COARSEN_INACTIVE);
     166       17078 :   libmesh_assert (!this->active());
     167             : 
     168             :   // We no longer delete children until MeshRefinement::contract()
     169             : 
     170       17078 :   unsigned int parent_p_level = 0;
     171             : 
     172      345905 :   const unsigned int n_n = this->n_nodes();
     173             : 
     174             :   // re-compute hanging node nodal locations
     175     1732257 :   for (unsigned int c = 0, nc = this->n_children(); c != nc; ++c)
     176             :     {
     177       68504 :       Elem * mychild = this->child_ptr(c);
     178     1386352 :       if (mychild == remote_elem)
     179       29746 :         continue;
     180     7287133 :       for (auto cnode : mychild->node_index_range())
     181             :         {
     182      381288 :           Point new_pos;
     183      381288 :           bool calculated_new_pos = false;
     184             : 
     185    39104974 :           for (unsigned int n=0; n<n_n; n++)
     186             :             {
     187             :               // The value from the embedding matrix
     188    33174375 :               const Real em_val = this->embedding_matrix(c,cnode,n);
     189             : 
     190             :               // The node location is somewhere between existing vertices
     191    33174375 :               if ((em_val != 0.) && (em_val != 1.))
     192             :                 {
     193    11455713 :                   new_pos.add_scaled (this->point(n), em_val);
     194      809208 :                   calculated_new_pos = true;
     195             :                 }
     196             :             }
     197             : 
     198     5930599 :           if (calculated_new_pos)
     199             :             {
     200             :               //Move the existing node back into it's original location
     201    15974880 :               for (unsigned int i=0; i<LIBMESH_DIM; i++)
     202             :                 {
     203     1446864 :                   Point & child_node = mychild->point(cnode);
     204    11981160 :                   child_node(i)=new_pos(i);
     205             :                 }
     206             :             }
     207             :         }
     208             :     }
     209             : 
     210     1732257 :   for (auto & mychild : this->child_ref_range())
     211             :     {
     212     1386352 :       if (&mychild == remote_elem)
     213       29746 :         continue;
     214       68432 :       libmesh_assert_equal_to (mychild.refinement_flag(), Elem::COARSEN);
     215       68432 :       mychild.set_refinement_flag(Elem::INACTIVE);
     216     1424966 :       if (mychild.p_level() > parent_p_level)
     217           0 :         parent_p_level = mychild.p_level();
     218             :     }
     219             : 
     220       17078 :   this->set_refinement_flag(Elem::JUST_COARSENED);
     221      345905 :   this->set_p_level(parent_p_level);
     222             : 
     223       17078 :   libmesh_assert (this->active());
     224      345905 : }
     225             : 
     226             : 
     227             : 
     228     5625044 : void Elem::contract()
     229             : {
     230             :   // Subactive elements get deleted entirely, not contracted
     231      330550 :   libmesh_assert (this->active());
     232             : 
     233             :   // Active contracted elements no longer can have children
     234      330550 :   _children.reset(nullptr);
     235             : 
     236     5625044 :   if (this->refinement_flag() == Elem::JUST_COARSENED)
     237       16498 :     this->set_refinement_flag(Elem::DO_NOTHING);
     238     5625044 : }
     239             : 
     240             : #endif // #ifdef LIBMESH_ENABLE_AMR
     241             : 
     242             : 
     243             : } // namespace libMesh