cell_hex20.c

一个用来实现偏微分方程中网格的计算库

  libmesh_assert (i < this->n_sides());

  if (proxy)
    {
      AutoPtr<Elem> ap(new Side<Quad8,Hex20>(this,i));
      return ap;
    }

  else
    {
      AutoPtr<Elem> face(new Quad8);

      // Think of a unit cube: (-1,1) x (-1,1)x (1,1)
      switch (i)
	{
	case 0:  // the face at z=0
	  {
	    face->set_node(0) = this->get_node(0);
	    face->set_node(1) = this->get_node(3);
	    face->set_node(2) = this->get_node(2);
	    face->set_node(3) = this->get_node(1);
	    face->set_node(4) = this->get_node(11);
	    face->set_node(5) = this->get_node(10);
	    face->set_node(6) = this->get_node(9);
	    face->set_node(7) = this->get_node(8);

	    return face;
	  }
	case 1:  // the face at y = 0
	  {
	    face->set_node(0) = this->get_node(0);
	    face->set_node(1) = this->get_node(1);
	    face->set_node(2) = this->get_node(5);
	    face->set_node(3) = this->get_node(4);
	    face->set_node(4) = this->get_node(8);
	    face->set_node(5) = this->get_node(13);
	    face->set_node(6) = this->get_node(16);
	    face->set_node(7) = this->get_node(12);

	    return face;
	  }
	case 2:  // the face at x=1
	  {
	    face->set_node(0) = this->get_node(1);
	    face->set_node(1) = this->get_node(2);
	    face->set_node(2) = this->get_node(6);
	    face->set_node(3) = this->get_node(5);
	    face->set_node(4) = this->get_node(9);
	    face->set_node(5) = this->get_node(14);
	    face->set_node(6) = this->get_node(17);
	    face->set_node(7) = this->get_node(13);

	    return face;
	  }
	case 3: // the face at y=1
	  {
	    face->set_node(0) = this->get_node(2);
	    face->set_node(1) = this->get_node(3);
	    face->set_node(2) = this->get_node(7);
	    face->set_node(3) = this->get_node(6);
	    face->set_node(4) = this->get_node(10);
	    face->set_node(5) = this->get_node(15);
	    face->set_node(6) = this->get_node(18);
	    face->set_node(7) = this->get_node(14);

	    return face;
	  }
	case 4: // the face at x=0
	  {
	    face->set_node(0) = this->get_node(3);
	    face->set_node(1) = this->get_node(0);
	    face->set_node(2) = this->get_node(4);
	    face->set_node(3) = this->get_node(7);
	    face->set_node(4) = this->get_node(11);
	    face->set_node(5) = this->get_node(12);
	    face->set_node(6) = this->get_node(19);
	    face->set_node(7) = this->get_node(15);

	    return face;
	  }
	case 5: // the face at z=1
	  {
	    face->set_node(0) = this->get_node(4);
	    face->set_node(1) = this->get_node(5);
	    face->set_node(2) = this->get_node(6);
	    face->set_node(3) = this->get_node(7);
	    face->set_node(4) = this->get_node(16);
	    face->set_node(5) = this->get_node(17);
	    face->set_node(6) = this->get_node(18);
	    face->set_node(7) = this->get_node(19);

	    return face;
	  }
	default:
	  {
	    libmesh_error();
	    return face;
	  }
	}
    }
  
  // We'll never get here.
  libmesh_error();
  AutoPtr<Elem> ap(NULL);  return ap;
}



AutoPtr<Elem> Hex20::build_edge (const unsigned int i) const
{
  libmesh_assert (i < this->n_edges());

  AutoPtr<Elem> ap(new SideEdge<Edge3,Hex20>(this,i));
  return ap;
}



void Hex20::connectivity(const unsigned int sc,
			 const IOPackage iop,
			 std::vector<unsigned int>& conn) const
{
  libmesh_assert (_nodes != NULL);
  libmesh_assert (sc < this->n_sub_elem());
  libmesh_assert (iop != INVALID_IO_PACKAGE);

  
  switch (iop)
    {
    case TECPLOT:
      {
	switch (sc)
	  {
	  case 0:
	    conn.resize(8);
	    conn[0] = this->node(0)+1;
	    conn[1] = this->node(1)+1;
	    conn[2] = this->node(2)+1;
	    conn[3] = this->node(3)+1;
	    conn[4] = this->node(4)+1;
	    conn[5] = this->node(5)+1;
	    conn[6] = this->node(6)+1;
	    conn[7] = this->node(7)+1;

	    return;