<commit_before>// The libMesh Finite Element Library.
// Copyright (C) 2002-2016 <NAME>, <NAME>, <NAME>

// 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



// <h1>Eigenproblems Example 3 - Can you "hear the shape" of a drum?</h1>
// \author <NAME>
// \date 2012
//
// The sound that a drum makes is determined by it's resonant frequencies,
// which are given by the eigenvalues of the Laplacian. "Can One Hear the
// Shape of a Drum?" was the title of an article by <NAME>
// in the American Mathematical Monthly in 1966, where he raised the question:
// If we know all the eigenvalues of a drum, can we uniquely determine it's shape?
// This question was resolved in 1992, when Gordon, Webb, and Wolpert constructed
// a pair of regions in 2D that have different shapes but identical eigenvalues.
// So the answer to Kac's question is no: the spectrum of the Laplacian does
// not uniquely determine the shape of the domain.

// In this example, we compute the first few eigenvalues of the two domains proposed
// by Gordon, Webb and Wolpert. This amounts to solving a generalized eigenvalue
// problem in each case. The computed eigenvalues are stored in drum1_evals.txt and
// drum2_evals.txt. We can compare these to the (highly accurate) values reported in:
// <NAME>, "Eigenmodes of Isospectral Drums", SIAM Review, Vol. 39, No. 1, pp. 1-17, 1997.
//
// The first five eigenvalues from Driscoll are listed below (the author states that
// "all digits shown are believed to be correct"):
//  2.53794399980
//  3.65550971352
//  5.17555935622
//  6.53755744376
//  7.24807786256


// C++ include files
#include <fstream>

// libMesh include files.
#include "libmesh/libmesh.h"
#include "libmesh/mesh.h"
#include "libmesh/mesh_generation.h"
#include "libmesh/exodusII_io.h"
#include "libmesh/condensed_eigen_system.h"
#include "libmesh/equation_systems.h"
#include "libmesh/fe.h"
#include "libmesh/quadrature_gauss.h"
#include "libmesh/dense_matrix.h"
#include "libmesh/sparse_matrix.h"
#include "libmesh/numeric_vector.h"
#include "libmesh/dof_map.h"
#include "libmesh/fe_interface.h"
#include "libmesh/getpot.h"
#include "libmesh/elem.h"
#include "libmesh/zero_function.h"
#include "libmesh/dirichlet_boundaries.h"

#define BOUNDARY_ID 100

// Bring in everything from the libMesh namespace
using namespace libMesh;


// Function prototype.  This is the function that will assemble
// the eigen system. Here, we will simply assemble a mass matrix.
void assemble_matrices(EquationSystems & es,
                       const std::string & system_name);

// We store the Dirichlet dofs in a set in order to impose the boundary conditions
void get_dirichlet_dofs(EquationSystems & es,
                        const std::string & system_name,
                        std::set<unsigned int> & global_dirichlet_dofs_set);


int main (int argc, char ** argv)
{
  // Initialize libMesh and the dependent libraries.
  LibMeshInit init (argc, argv);

  // This example uses an ExodusII input file
#ifndef LIBMESH_HAVE_EXODUS_API
  libmesh_example_requires(false, "--enable-exodus");
#endif

  // This example is designed for the SLEPc eigen solver interface.
#ifndef LIBMESH_HAVE_SLEPC
  if (init.comm().rank() == 0)
    libMesh::err << "ERROR: This example requires libMesh to be\n"
                 << "compiled with SLEPc eigen solvers support!"
                 << std::endl;

  return 0;
#else

#ifdef LIBMESH_DEFAULT_SINGLE_PRECISION
  // SLEPc currently gives us a nasty crash with Real==float
  libmesh_example_requires(false, "--disable-singleprecision");
#endif

#ifdef LIBMESH_USE_COMPLEX_NUMBERS
  // SLEPc currently gives us an "inner product not well defined" with
  // Number==complex
  libmesh_example_requires(false, "--disable-complex");
#endif

  // Tell the user what we are doing.
  {
    libMesh::out << "Running " << argv[0];

    for (int i=1; i<argc; i++)
      libMesh::out << " " << argv[i];

    libMesh::out << std::endl << std::endl;
  }

  // Skip this 2D example if libMesh was compiled as 1D-only.
  libmesh_example_requires(2 <= LIBMESH_DIM, "2D support");

  // Use GetPot to parse the command line arguments
  GetPot command_line (argc, argv);

  // Read the mesh name from the command line
  std::string mesh_name = "";
  if (command_line.search(1, "-mesh_name"))
    mesh_name = command_line.next(mesh_name);

  // Also, read in the index of the eigenvector that we should plot
  // (zero-based indexing, as usual!)
  unsigned int plotting_index = 0;
  if (command_line.search(1, "-plotting_index"))
    plotting_index = command_line.next(plotting_index);

  // Finally, read in the number of eigenpairs we want to compute!
  unsigned int n_evals = 0;
  if (command_line.search(1, "-n_evals"))
    n_evals = command_line.next(n_evals);

  // Append the .e to mesh_name
  std::ostringstream mesh_name_exodus;
  mesh_name_exodus << mesh_name << "_mesh.e";

  // Create a mesh, with dimension to be overridden by the file, on
  // the default MPI communicator.
  Mesh mesh(init.comm());

  mesh.read(mesh_name_exodus.str());

  {
    // Add boundary IDs to this mesh so that we can use DirichletBoundary
  // Each processor should know about each boundary condition it can
  // see, so we loop over all elements, not just local elements.
    MeshBase::const_element_iterator       el     = mesh.elements_begin();
    const MeshBase::const_element_iterator end_el = mesh.elements_end();
    for ( ; el != end_el; ++el)
      {
        const Elem * elem = *el;

        for (unsigned int side=0; side<elem->n_sides(); side++)
          {
            if(elem->neighbor(side) == NULL)
              {
                mesh.get_boundary_info().add_side(elem, side, BOUNDARY_ID);
              }
          }
      }
  }

  // Print information about the mesh to the screen.
  mesh.print_info();

  // Create an equation systems object.
  EquationSystems equation_systems (mesh);

  // Create a CondensedEigenSystem named "Eigensystem" and (for convenience)
  // use a reference to the system we create.
  CondensedEigenSystem & eigen_system =
    equation_systems.add_system<CondensedEigenSystem> ("Eigensystem");

  // Declare the system variables.
  // Adds the variable "p" to "Eigensystem".   "p"<commit_msg>Fix indentation.<commit_after>// The libMesh Finite Element Library.
// Copyright (C) 2002-2016 <NAME>, <NAME>, <NAME>

// 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



// <h1>Eigenproblems Example 3 - Can you "hear the shape" of a drum?</h1>
// \author <NAME>
// \date 2012
//
// The sound that a drum makes is determined by it's resonant frequencies,
// which are given by the eigenvalues of the Laplacian. "Can One Hear the
// Shape of a Drum?" was the title of an article by <NAME>
// in the American Mathematical Monthly in 1966, where he raised the question:
// If we know all the eigenvalues of a drum, can we uniquely determine it's shape?
// This question was resolved in 1992, when Gordon, Webb, and Wolpert constructed
// a pair of regions in 2D that have different shapes but identical eigenvalues.
// So the answer to Kac's question is no: the spectrum of the Laplacian does
// not uniquely determine the shape of the domain.

// In this example, we compute the first few eigenvalues of the two domains proposed
// by Gordon, Webb and Wolpert. This amounts to solving a generalized eigenvalue
// problem in each case. The computed eigenvalues are stored in drum1_evals.txt and
// drum2_evals.txt. We can compare these to the (highly accurate) values reported in:
// <NAME>, "Eigenmodes of Isospectral Drums", SIAM Review, Vol. 39, No. 1, pp. 1-17, 1997.
//
// The first five eigenvalues from Driscoll are listed below (the author states that
// "all digits shown are believed to be correct"):
//  2.53794399980
//  3.65550971352
//  5.17555935622
//  6.53755744376
//  7.24807786256


// C++ include files
#include <fstream>

// libMesh include files.
#include "libmesh/libmesh.h"
#include "libmesh/mesh.h"
#include "libmesh/mesh_generation.h"
#include "libmesh/exodusII_io.h"
#include "libmesh/condensed_eigen_system.h"
#include "libmesh/equation_systems.h"
#include "libmesh/fe.h"
#include "libmesh/quadrature_gauss.h"
#include "libmesh/dense_matrix.h"
#include "libmesh/sparse_matrix.h"
#include "libmesh/numeric_vector.h"
#include "libmesh/dof_map.h"
#include "libmesh/fe_interface.h"
#include "libmesh/getpot.h"
#include "libmesh/elem.h"
#include "libmesh/zero_function.h"
#include "libmesh/dirichlet_boundaries.h"

#define BOUNDARY_ID 100

// Bring in everything from the libMesh namespace
using namespace libMesh;


// Function prototype.  This is the function that will assemble
// the eigen system. Here, we will simply assemble a mass matrix.
void assemble_matrices(EquationSystems & es,
                       const std::string & system_name);

// We store the Dirichlet dofs in a set in order to impose the boundary conditions
void get_dirichlet_dofs(EquationSystems & es,
                        const std::string & system_name,
                        std::set<unsigned int> & global_dirichlet_dofs_set);


int main (int argc, char ** argv)
{
  // Initialize libMesh and the dependent libraries.
  LibMeshInit init (argc, argv);

  // This example uses an ExodusII input file
#ifndef LIBMESH_HAVE_EXODUS_API
  libmesh_example_requires(false, "--enable-exodus");
#endif

  // This example is designed for the SLEPc eigen solver interface.
#ifndef LIBMESH_HAVE_SLEPC
  if (init.comm().rank() == 0)
    libMesh::err << "ERROR: This example requires libMesh to be\n"
                 << "compiled with SLEPc eigen solvers support!"
                 << std::endl;

  return 0;
#else

#ifdef LIBMESH_DEFAULT_SINGLE_PRECISION
  // SLEPc currently gives us a nasty crash with Real==float
  libmesh_example_requires(false, "--disable-singleprecision");
#endif

#ifdef LIBMESH_USE_COMPLEX_NUMBERS
  // SLEPc currently gives us an "inner product not well defined" with
  // Number==complex
  libmesh_example_requires(false, "--disable-complex");
#endif

  // Tell the user what we are doing.
  {
    libMesh::out << "Running " << argv[0];

    for (int i=1; i<argc; i++)
      libMesh::out << " " << argv[i];

    libMesh::out << std::endl << std::endl;
  }

  // Skip this 2D example if libMesh was compiled as 1D-only.
  libmesh_example_requires(2 <= LIBMESH_DIM, "2D support");

  // Use GetPot to parse the command line arguments
  GetPot command_line (argc, argv);

  // Read the mesh name from the command line
  std::string mesh_name = "";
  if (command_line.search(1, "-mesh_name"))
    mesh_name = command_line.next(mesh_name);

  // Also, read in the index of the eigenvector that we should plot
  // (zero-based indexing, as usual!)
  unsigned int plotting_index = 0;
  if (command_line.search(1, "-plotting_index"))
    plotting_index = command_line.next(plotting_index);

  // Finally, read in the number of eigenpairs we want to compute!
  unsigned int n_evals = 0;
  if (command_line.search(1, "-n_evals"))
    n_evals = command_line.next(n_evals);

  // Append the .e to mesh_name
  std::ostringstream mesh_name_exodus;
  mesh_name_exodus << mesh_name << "_mesh.e";

  // Create a mesh, with dimension to be overridden by the file, on
  // the default MPI communicator.
  Mesh mesh(init.comm());

  mesh.read(mesh_name_exodus.str());

  {
    // Add boundary IDs to this mesh so that we can use DirichletBoundary
    // Each processor should know about each boundary condition it can
    // see, so we loop over all elements, not just local elements.
    MeshBase::const_element_iterator       el     = mesh.elements_begin();
    const MeshBase::const_element_iterator end_el = mesh.elements_end();
    for ( ; el != end_el; ++el)
      {
        const Elem * elem = *el;

        for (unsigned int side=0; side<elem->n_sides(); side++)
          {
            if(elem->neighbor(side) == NULL)
              {
                mesh.get_boundary_info().add_side(elem, side, BOUNDARY_ID);
              }
          }
      }
  }

  // Print information about the mesh to the screen.
  mesh.print_info();

  // Create an equation systems object.
  EquationSystems equation_systems (mesh);

  // Create a CondensedEigenSystem named "Eigensystem" and (for convenience)
  // use a reference to the system we create.
  CondensedEigenSystem & eigen_system =
    equation_systems.add_system<CondensedEigenSystem> ("Eigensystem");

  // Declare the system variables.
  // Adds the variable "p" to "Eigensystem".   "p"
{'fraction_non_alphanumeric': 0.09963152951376346, 'fraction_numerical': 0.018640271656672206, 'mean_word_length': 3.3734597156398105, 'pattern_counts': {'":': 0, '<': 75, '<?xml version=': 0, '>': 29, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'dir': 'git-commits-cleaned', 'id': '1887990', 'n_tokens_mistral': 4473, 'n_tokens_neox': 4094, 'n_words': 1566}