mesh_data.c

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

// $Id: mesh_data.C 2789 2008-04-13 02:24:40Z roystgnr $

// The libMesh Finite Element Library.
// Copyright (C) 2002-2007  Benjamin S. Kirk, John W. Peterson
  
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
  
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
// Lesser General Public License for more details.
  
// You should have received a copy of the GNU Lesser General Public
// License along with this library; if not, write to the Free Software
// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA



// C++ includes
#include <sstream>

// Local includes
#include "mesh_data.h"
#include "mesh_base.h"
#include "libmesh_logging.h"
#include "elem.h"


//------------------------------------------------------
// MeshData functions
MeshData::MeshData(const MeshBase& m) :
  _mesh               (m),
  _data_descriptor    (""),
  _node_id_map_closed (false),
  _node_data_closed   (false),
  _elem_id_map_closed (false),
  _elem_data_closed   (false),
  _active             (false),
  _compatibility_mode (false),
  _unv_header         (NULL)
{
}




MeshData::~MeshData()
{
  clear();
}





void MeshData::activate (const std::string& descriptor)
{
#ifdef DEBUG
  if (_compatibility_mode)
      std::cerr << "WARNING: MeshData was in compatibility mode, now being activated."
		<< std::endl;
#endif

  _compatibility_mode = false;
  _active = true;
  _data_descriptor = descriptor;
}





void MeshData::enable_compatibility_mode (const std::string& descriptor)
{
  if (!_active)
    {
      _compatibility_mode = true;
      _active = false;
      // do as if the id maps are already closed
      _node_id_map_closed = true;
      _elem_id_map_closed = true;
      _data_descriptor = descriptor;
      // we can safely clear the id maps
      _node_id.clear();
      _id_node.clear();
      _elem_id.clear();
      _id_elem.clear();
    }
#ifdef DEBUG
  else
      std::cerr << "WARNING: MeshData was in compatibility mode, now being activated."
		<< std::endl;
#endif
}






void MeshData::clear()
{
  _data_descriptor    = "";
  _node_data.clear();
  _elem_data.clear();
  _node_data_closed   = false;
  _elem_data_closed   = false;
}





void  MeshData::slim (const bool node_id_map,
		      const bool elem_id_map)
{
  if (this->active())
    {

      if (node_id_map)
        {
	  // dumb check
	  libmesh_assert (_node_id_map_closed);

	  _node_id_map_closed = false;
	  _node_id.clear();
	  _id_node.clear();
	}

      if (elem_id_map)
        {
	  // dumb check
	  libmesh_assert (_elem_id_map_closed);

	  _elem_id_map_closed = false;
	  _elem_id.clear();
	  _id_elem.clear();
	}
    }

#ifdef DEBUG
  else if (this->compatibility_mode())
    {
      std::cerr << "WARNING: No need for MeshData::slim() in compatibility mode." << std::endl;
    }
#endif
}




void MeshData::translate (const MeshBase& out_mesh,
			  std::vector<Number>& values,
			  std::vector<std::string>& names) const
{
  libmesh_assert (_active || _compatibility_mode);

  START_LOG("translate()", "MeshData");

  const unsigned int n_comp = this->n_val_per_node();

  // transfer our nodal data to a vector
  // that may be written concurrently
  // with the \p out_mesh.
  {
    // reserve memory for the nodal data
    values.reserve(n_comp*out_mesh.n_nodes());
    
    // iterate over the mesh's nodes
    MeshBase::const_node_iterator       nodes_it  = out_mesh.nodes_begin();
    const MeshBase::const_node_iterator nodes_end = out_mesh.nodes_end();

    // Do not use the \p get_data() method, but the operator()
    // method, since this returns by default a zero value,
    // when there is no nodal data.
    for (; nodes_it != nodes_end; ++nodes_it)
      {
	const Node* node = *nodes_it;
	
	for (unsigned int c= 0; c<n_comp; c++)
	  values.push_back(this->operator()(node, c));
      }
  }


  
  // Now we have the data, nicely stored in \p values.
  // It remains to give names to the data, then write to
  // file.
  {
    names.reserve(n_comp);
    
    // this naming scheme only works up to n_comp=100
    // (at least for gmv-accepted variable names)
    libmesh_assert(n_comp < 100);

    for (unsigned int n=0; n<n_comp; n++)
      {
	std::ostringstream name_buf;
	name_buf << "bc_" << n;
	names.push_back(name_buf.str());
      }
  }

  STOP_LOG("translate()", "MeshData");
}




void MeshData::close_foreign_id_maps ()
{
  if (_active)
    {
      libmesh_assert (!_elem_id.empty());
      libmesh_assert (!_id_elem.empty());
      libmesh_assert (!_node_id.empty());
      libmesh_assert (!_id_node.empty());

      _elem_id_map_closed = true;
      _node_id_map_closed = true;
    }
}





void MeshData::read (const std::string& name)
{
  START_LOG("read()", "MeshData");

  libmesh_assert (_active || _compatibility_mode);

  // the id maps have to be closed before reading
  // (note that in compatibility mode these are also true)
  libmesh_assert (_elem_id_map_closed && _node_id_map_closed);

#ifdef DEBUG
  if (this->compatibility_mode())
      std::cerr << "WARNING: MeshData in compatibility mode, node and element ids" << std::endl
		<< "         stored in file may be totally different from libMesh ids!" << std::endl;
#endif

  // Read the file based on extension
  // For now, only processor 0 should read the
  // mesh data.
  if (libMesh::processor_id() == 0)
    {
      if (name.rfind(".xta") < name.size())
	this->read_xdr (name, READ);
      
      else if (name.rfind(".xtr")  < name.size())
	this->read_xdr (name, DECODE);
      
      else if (name.rfind(".unv") < name.size())
	this->read_unv (name);
      
      else if ((name.rfind(".node") < name.size()) ||
	       (name.rfind(".ele") < name.size()))
	this->read_tetgen (name);
      
      else
	{
	  std::cerr << " ERROR: Unrecognized file extension: " << name
		    << "\n   I understand the following:\n\n"
		    << "     *.xta  -- Internal ASCII data format\n"
		    << "     *.xtr  -- Internal binary data format\n"
		    << "     *.unv  -- I-deas format\n"
		    << std::endl;
	  libmesh_error();
	  
	}    
    }
  STOP_LOG("read()", "MeshData");
}






void MeshData::write (const std::string& name)
{
  START_LOG("write()", "MeshData");

  libmesh_assert (_active || _compatibility_mode);

  // the id maps have to be closed before writing
  // (note that in compatibility mode these are also true)
  libmesh_assert (_elem_id_map_closed && _node_id_map_closed);
  
#ifdef DEBUG
  if (this->compatibility_mode())
      std::cerr << "WARNING: MeshData in compatibility mode.  Node and element ids" << std::endl
		<< "         written to file may differ from libMesh numbering" << std::endl
		<< "         next time this file is read!" << std::endl;
#endif

  // Read the file based on extension
  {
    if (name.rfind(".xta") < name.size())
      write_xdr (name, WRITE);

    else if (name.rfind(".xtr")  < name.size())
      write_xdr (name, ENCODE);

    else if (name.rfind(".unv") < name.size())
      write_unv (name);

    else
      {
	std::cerr << " ERROR: Unrecognized file extension: " << name
		  << "\n   I understand the following:\n\n"
		  << "     *.xta  -- Internal ASCII data format\n"
		  << "     *.xtr  -- Internal binary data format\n"
		  << "     *.unv  -- I-deas format\n"
		  << std::endl;
	libmesh_error();

      }    
  }
  STOP_LOG("write()", "MeshData");
}




std::string MeshData::get_info() const
{
  std::ostringstream out;

  if (this->active() || this->compatibility_mode())
    {
      out << " MeshData Information:\n";
      if (this->active())
	  out << "  object activated.\n";
      if (this->compatibility_mode())
	  out << "  object in compatibility mode.\n";
      if (this->_data_descriptor != "")
	  out << "  descriptor=" << this->_data_descriptor << '\n';
      if (this->elem_initialized())
	  out << "  Element associated data initialized.\n"
	      << "   n_val_per_elem()=" << this->n_val_per_elem() << '\n'
	      << "   n_elem_data()=" << this->n_elem_data() << '\n';
      if (this->node_initialized())
	  out << "  Node associated data initialized.\n"
	      << "   n_val_per_node()=" << this->n_val_per_node() << '\n'
	      << "   n_node_data()=" << this->n_node_data() << '\n';
    }
  else
      out << " MeshData neither active nor in compatibility mode.\n";

  return out.str();
}




void MeshData::print_info(std::ostream& os) const
{
  os << this->get_info()
     << std::endl;
}


std::ostream& operator << (std::ostream& os, const MeshData& m)
{
  m.print_info(os);
  return os;
}




const Node* MeshData::foreign_id_to_node (const unsigned int fid) const
{
  if (_active)
    {