meshdiff.C File Reference

Go to the source code of this file.

Functions
int	main (int argc, char **argv)

Variables
unsigned char	dim = 2

Function Documentation

main()

int main	(	int	argc,
		char **	argv
	)

Definition at line 37 of file meshdiff.C.

References libMesh::ExactSolution::attach_reference_solution(), libMesh::System::calculate_norm(), libMesh::ExactSolution::compute_error(), dim, libMesh::System::get_dof_map(), libMesh::H1, libMesh::ExactSolution::h1_error(), libMesh::H2, libMesh::ExactSolution::h2_error(), libMesh::System::has_variable(), libMesh::index_range(), libMesh::TriangleWrapper::init(), libMesh::MeshFunction::init(), libMesh::invalid_uint, libMesh::L2, libMesh::ExactSolution::l2_error(), libMesh::System::n_vars(), libMesh::out, libMesh::System::project_solution(), libMesh::System::solution, libMesh::System::variable_name(), libMesh::System::variable_number(), and libMesh::NameBasedIO::write_equation_systems().

{
  LibMeshInit init(argc, argv);

  Mesh coarse_mesh(init.comm(), dim), fine_mesh(init.comm(), dim);
  EquationSystems coarse_es(coarse_mesh), fine_es(fine_mesh);

  libMesh::out << "Usage: " << argv[0]
               << " coarsemesh coarsesolution finemesh finesolution [outputdiff]" << std::endl;

  if (argc < 5)
    libmesh_error();

  coarse_mesh.read(argv[1]);
  libMesh::out << "Loaded coarse mesh " << argv[1] << std::endl;
  coarse_es.read(argv[2]);
  libMesh::out << "Loaded coarse solution " << argv[2] << std::endl;
  fine_mesh.read(argv[3]);
  libMesh::out << "Loaded fine mesh " << argv[3] << std::endl;
  fine_es.read(argv[4]);
  libMesh::out << "Loaded fine solution " << argv[4] << std::endl;

  ExactSolution exact_sol(coarse_es);
  exact_sol.attach_reference_solution(&fine_es);

  std::vector<std::string> sysnames;
  sysnames.reserve(coarse_es.n_systems());

  for (unsigned int i = 0; i != coarse_es.n_systems(); ++i)
    if (fine_es.has_system(coarse_es.get_system(i).name()))
      sysnames.push_back(coarse_es.get_system(i).name());
    else
      libMesh::out << "Coarse solution system "
                   << coarse_es.get_system(i).name()
                   << " not found in fine solution!" << std::endl;

  for (unsigned int i = 0; i != fine_es.n_systems(); ++i)
    if (!coarse_es.has_system(fine_es.get_system(i).name()))
      libMesh::out << "Fine solution system "
                   << fine_es.get_system(i).name()
                   << " not found in coarse solution!" << std::endl;

  if (!coarse_es.n_systems() && !fine_es.n_systems())
    libMesh::out << "No systems found in fine or coarse solution!"
                 << std::endl;

  for (auto i : index_range(sysnames))
    {
      const std::string sysname = sysnames[i];
      const System & coarse_sys = coarse_es.get_system(sysname);
      const System & fine_sys = fine_es.get_system(sysname);

      for (unsigned int j = 0; j != coarse_sys.n_vars(); ++j)
        {
          const std::string varname = coarse_sys.variable_name(j);

          if (!fine_sys.has_variable(varname))
            {
              libMesh::out << "Fine solution system " << sysname
                           << " variable " << varname
                           << " not found in coarse solution!" << std::endl;
              continue;
            }
          const unsigned int j2 = fine_sys.variable_number(varname);

          exact_sol.compute_error(sysname, varname);

          libMesh::out << "Errors in system " << sysname << ", variable " << varname << ":" << std::endl;
          libMesh::out << "L2 error: " << exact_sol.l2_error(sysname, varname)
                       << ", fine norm: " << coarse_sys.calculate_norm(*coarse_sys.solution, j, L2)
                       << ", coarse norm: " << fine_sys.calculate_norm(*fine_sys.solution, j2, L2) << std::endl;
          libMesh::out << "H1 error: " << exact_sol.h1_error(sysname, varname)
                       << ", fine norm: " << coarse_sys.calculate_norm(*coarse_sys.solution, j, H1)
                       << ", coarse norm: " << fine_sys.calculate_norm(*fine_sys.solution, j2, H1) << std::endl;
          libMesh::out << "H2 error: " << exact_sol.h2_error(sysname, varname)
                       << ", fine norm: " << coarse_sys.calculate_norm(*coarse_sys.solution, j, H2)
                       << ", coarse norm: " << fine_sys.calculate_norm(*fine_sys.solution, j2, H2) << std::endl;
        }

      for (unsigned int j = 0; j != fine_sys.n_vars(); ++j)
        {
          const std::string varname = fine_sys.variable_name(j);

          if (!coarse_sys.has_variable(varname))
            {
              libMesh::out << "Coarse solution system " << sysname
                           << " variable " << varname
                           << " not found in fine solution!" << std::endl;
              continue;
            }
        }

      if (!coarse_sys.n_vars() && !fine_sys.n_vars())
        libMesh::out << "No variables found in fine or coarse solution system "
                     << sysname << '!' << std::endl;
    }

  if (argc > 5)
    {
      libMesh::out << "Writing diff solution " << argv[5] << std::endl;

      for (auto i : index_range(sysnames))
        {
          const std::string sysname = sysnames[i];
          const System & coarse_sys = coarse_es.get_system(sysname);
          const System & fine_sys = fine_es.get_system(sysname);

          std::unique_ptr<NumericVector<Number>> fine_solution = fine_sys.solution->clone();
          fine_sys.solution->zero();

          std::vector<unsigned int>
            var_remapping(fine_sys.n_vars(), libMesh::invalid_uint);

          for (unsigned int j = 0; j != fine_sys.n_vars(); ++j)
            {
              const std::string varname = fine_sys.variable_name(j);

              if (coarse_sys.has_variable(varname))
                var_remapping[j] = coarse_sys.variable_number(varname);
            }

          MeshFunction coarse_solution
            (coarse_es, *coarse_sys.solution,
             coarse_sys.get_dof_map(), std::move(var_remapping));
          coarse_solution.init();
          const DenseVector<Number> oom_value(fine_sys.n_vars(), 0);
          coarse_solution.enable_out_of_mesh_mode(oom_value);

          fine_sys.project_solution(&coarse_solution);
          *fine_sys.solution -= *fine_solution;
        }

      NameBasedIO(fine_mesh).write_equation_systems (argv[5], fine_es);
    }

  return 0;
}

Variable Documentation

dim

unsigned char dim = 2

Definition at line 35 of file meshdiff.C.

Referenced by main().