three_d_mesh_dist.cc File Reference

#include "generic.h"
#include "poisson.h"
#include "meshes/simple_cubic_mesh.h"

Classes
class	ThreeDPoissonProblem< ELEMENT >
	Poisson problem. More...

Namespaces
	TanhSolnForPoisson
	Namespace for exact solution for Poisson equation with "sharp step".

Functions
void	TanhSolnForPoisson::get_exact_u (const Vector< double > &x, Vector< double > &u)
	Exact solution as a Vector. More...
void	TanhSolnForPoisson::get_exact_u (const Vector< double > &x, double &u)
	Exact solution as a scalar. More...
void	TanhSolnForPoisson::get_source (const Vector< double > &x, double &source)
	Source function to make it an exact solution. More...
void	parallel_test (const unsigned &n_refine_first, const unsigned &n_refine_once_distributed, const bool &redistribute, const bool &solve, const bool &doc)
int	main (int argc, char **argv)
	Driver code for 3D Poisson problem. More...

Function Documentation

main()

int main	(	int argc	,
		char **	argv
	)

Driver code for 3D Poisson problem.

{
#ifdef OOMPH_HAS_MPI
 // Setup MPI helpers
 MPI_Helpers::init(argc,argv);
#endif
 
 // Store command line arguments
 CommandLineArgs::setup(argc,argv);
 
 
 // Default min. number of refinements before we distribute
 unsigned n_refine_first_min=1;
 
 // Default max. number of refinements before we distribute
 unsigned n_refine_first_max=1;
 
 // Default min. number of uniform refinements after distribution
 unsigned n_uniform_min=1; 
 
 // Default max. number of uniform refinements after distribution
 unsigned n_uniform_max=2; 
 
 // Redistribute after uniform refinement?
 bool redistribute=true;
                 
 
        
 if (CommandLineArgs::Argc==1)
  {
   std::cout << "Using default settings for flags" << std::endl;
  }
 else if (CommandLineArgs::Argc==6)
  {
 
    // Min. number of refinements before we distribute
    n_refine_first_min=atoi(argv[1]);
 
    // Max. number of refinements before we distribute
    n_refine_first_max=atoi(argv[2]);
 
    // Min. number of uniform refinements after distribution
    n_uniform_min=atoi(argv[3]); 
 
    // Max. number of uniform refinements after distribution
    n_uniform_max=atoi(argv[4]); 
 
    // Redistribute after uniform refinement?
    redistribute=atoi(argv[5]); 
  }
 else
  {
   oomph_info 
    << "Wrong number of input args. Specify none or five " 
    << "(n_refine_first_min/max, n_uniform_min/max, redistributon flag)" 
    << std::endl;
   exit(1);
  }
 
 
 
 // Tell us what we're going to do
 oomph_info << "Run with  n_refine_first_min/max, n_uniform_min/max, redistribute: " 
            << n_refine_first_min << " "
            << n_refine_first_max << " "
            << n_uniform_min   << " " 
            << n_uniform_max   << " " 
            << redistribute  << " " 
            << std::endl;
 
 
 // switch off output for all but one processor
// if (my_rank>0)
//  {
//   oomph_info.stream_pt() = &oomph_nullstream;
//  }
 
 
 // Solve and self test?
 bool solve=true; // (andy) only testing mesh distribution here
 bool doc=true;
 
 // Write headers for analysis
 ofstream some_file;
 some_file.open("RESLT/three_d_initial_refinement_header.dat");
 some_file << "VARIABLES= \"n<SUB>refine first</SUB>\"," 
           << "\"n<SUB>dof orig</SUB>\","
           << "\"n<SUB>dof</SUB>\","
           << "\"t<SUB>serial adapt</SUB>\""
           << std::endl;
 some_file.close();    
 
 
 some_file.open("RESLT/three_d_initial_distr_header.dat");
 some_file << "VARIABLES= \"n_<SUB>p</SUB>\","
           << "\"n<SUB>dof</SUB>\","
           << "\"t<SUB>distr</SUB>\","
           << "\"average n<SUB>halo nodes</SUB>\","
           << "\"max. n<SUB>halo nodes</SUB>\","
           << "\"min. n<SUB>halo nodes</SUB>\","
           << "\"average efficiency\","
           << "\"max. efficiency\","
           << "\"min. efficiency\""
           << std::endl;
 some_file.close();    
 
 
 some_file.open("RESLT/three_d_distributed_refinement_header.dat");
 some_file << "VARIABLES= \"n_<SUB>p</SUB>\","
           << "\"n<SUB>dof orig</SUB>\","
           << "\"n<SUB>dof</SUB>\","
           << "\"t<SUB>distributed refinement</SUB>\","
           << "\"average n<SUB>halo nodes</SUB>\","
           << "\"max. n<SUB>halo nodes</SUB>\","
           << "\"min. n<SUB>halo nodes</SUB>\","
           << "\"average efficiency\","
           << "\"max. efficiency\","
           << "\"min. efficiency\""
           << std::endl;
 some_file.close();   
 
 
 some_file.open("RESLT/three_d_redistribution_header.dat");
 some_file << "VARIABLES= \"n_<SUB>p</SUB>\","
           << "\"n<SUB>dof</SUB>\","
           << "\"t<SUB>redistribution</SUB>\","
           << "\"average n<SUB>halo nodes</SUB>\","
           << "\"max. n<SUB>halo nodes</SUB>\","
           << "\"min. n<SUB>halo nodes</SUB>\","
           << "\"average efficiency\","
           << "\"max. efficiency\","
           << "\"min. efficiency\""
           << std::endl;
 some_file.close();   
 
 
 
 some_file.open("RESLT/three_d_distributed_unrefinement_header.dat");
 some_file << "VARIABLES= \"n_<SUB>p</SUB>\","
           << "\"n<SUB>dof orig</SUB>\","
           << "\"n<SUB>dof</SUB>\","
           << "\"t<SUB>distributed unrefinement</SUB>\","
           << "\"average n<SUB>halo nodes</SUB>\","
           << "\"max. n<SUB>halo nodes</SUB>\","
           << "\"min. n<SUB>halo nodes</SUB>\","
           << "\"average efficiency\","
           << "\"max. efficiency\","
           << "\"min. efficiency\""
           << std::endl;
 some_file.close();   
 
 
 // Test the procedures in parallel mode
 for (unsigned i=n_refine_first_min;i<=n_refine_first_max;i++)
 
  {
   for (unsigned j=n_uniform_min;j<=n_uniform_max;j++)
    {
     parallel_test(i,j,redistribute,solve,doc);
    }
  }
 
#ifdef OOMPH_HAS_MPI
 MPI_Helpers::finalize();
#endif
 
}

References oomph::CommandLineArgs::Argc, i, j, oomph::oomph_info, parallel_test(), Flag_definition::setup(), and solve.

parallel_test()

void parallel_test	(	const unsigned &	n_refine_first,
		const unsigned &	n_refine_once_distributed,
		const bool &	redistribute,
		const bool &	solve,
		const bool &	doc
	)

///////////////////////////////////////////////////////////////////// ///////////////////////////////////////////////////////////////////// ///////////////////////////////////////////////////////////////////// Test mesh distribution in parallel, starting with n_uniform_first uniform refinments before starting the distribution. Then refine and unrefine uniformly as specified. Re-distribute after every uniform refinement if requested by the boolean flag.

{
 
 oomph_info 
  << std::endl  << std::endl
  << "Running with " << n_refine_first << " initial serial refinements and \n"
  << "with " << n_refine_once_distributed 
  << " subsequent distributed refinements"
  << std::endl << std::endl;
 
 double av_efficiency; 
 double max_efficiency;
 double min_efficiency;
 double av_number_halo_nodes;
 unsigned max_number_halo_nodes;
 unsigned min_number_halo_nodes;
 
 char filename[200];
 std::ofstream some_file;
 
 DocInfo doc_info;
 doc_info.set_directory("RESLT");
 
 DocInfo quiet_doc_info;
 quiet_doc_info.set_directory("RESLT");
 
 // Create the problem
 //-------------------
 ThreeDPoissonProblem<RefineableQPoissonElement<3,3> >* problem_pt=
  new ThreeDPoissonProblem<RefineableQPoissonElement<3,3> >
  (&TanhSolnForPoisson::get_source);
 
 // Storage for number of processors and current processor
 int n_proc=problem_pt->communicator_pt()->nproc();
 int my_rank=problem_pt->communicator_pt()->my_rank(); 
 
 // Initial number of dofs in base mesh
 unsigned n_dof_orig=problem_pt->ndof();
 unsigned n_dof=problem_pt->ndof();
 
 
 // Start with a few rounds of uniform refinement
 //----------------------------------------------
 clock_t t_start = clock();
 for (unsigned i=0;i<n_refine_first;i++)
  {
   // Refine uniformly
   problem_pt->refine_uniformly();
  } 
 clock_t t_end = clock();
 
 
 // Doc stats
 if (my_rank==0)
  {
   double t_tot=double(t_end-t_start)/CLOCKS_PER_SEC;
   oomph_info << "Time for initial adaptation [ndof=" << problem_pt->ndof() 
              << "]: " << t_tot << std::endl;
 
   n_dof=problem_pt->ndof();
   sprintf(
    filename,
    "%s/three_d_initial_refinement_nrefinement_%i_ndofbefore_%i_ndofafter_%i.dat",
    doc_info.directory().c_str(),
    n_refine_first,
    n_dof_orig,
    n_dof);
   some_file.open(filename);
   some_file << n_refine_first << " " 
             << n_dof_orig << " " 
             << n_dof << " " 
             << t_tot << " "
             << std::endl;
   some_file.close();    
 }
 
 
 
 // Now distribute the problem
 //---------------------------
 t_start = clock();
 if (doc)
  {
   std::ifstream input_file;
   std::ofstream output_file;
 
   // Get partition from file
   unsigned n_partition=problem_pt->mesh_pt()->nelement();
   Vector<unsigned> element_partition(n_partition,0);
   sprintf(filename,"three_d_mesh_dist_partition.dat");
   input_file.open(filename);
   std::string input_string;
   for (unsigned e=0;e<n_partition;e++)
    {
     getline(input_file,input_string,'\n');
     element_partition[e]=atoi(input_string.c_str());
    }
 
//   Vector<unsigned> out_element_partition;
   bool report_stats=false;
   problem_pt->distribute(element_partition,doc_info,report_stats);
//                           out_element_partition);
 
//    sprintf(filename,"out_three_d_mesh_dist_partition.dat");
//    output_file.open(filename);
//    for (unsigned e=0;e<n_partition;e++)
//     {
//      output_file << out_element_partition[e] << std::endl;
//     }
 
//   problem_pt->distribute (doc_info,report_stats);
  }
 else
  {
   problem_pt->distribute();
  }
 t_end = clock();
 
 
 // Assess the quality of the distribution
 //---------------------------------------
 problem_pt->mesh_pt()->get_halo_node_stats(av_number_halo_nodes,
                                            max_number_halo_nodes,
                                            min_number_halo_nodes);
 problem_pt->mesh_pt()->get_efficiency_of_mesh_distribution
  (av_efficiency,max_efficiency,min_efficiency);
 
 
 // Doc stats
 if (my_rank==0)
  {
   double t_tot=double(t_end-t_start)/CLOCKS_PER_SEC;
   oomph_info << "Time for problem distribution [ndof=" << problem_pt->ndof() 
              << "]: " << t_tot << std::endl;
 
   sprintf(
    filename,
    "%s/three_d_initial_distr_np_%i_nrefinement_%i_ndof_%i.dat",
    doc_info.directory().c_str(),
    n_proc,
    n_refine_first,
    n_dof);
   some_file.open(filename);
   some_file << n_proc << " " 
             << n_dof << " " 
             << t_tot << " "
             << av_number_halo_nodes << " " 
             << max_number_halo_nodes << " " 
             << min_number_halo_nodes  << " " 
             << av_efficiency << " " 
             << max_efficiency << " " 
             << min_efficiency << " " 
             << std::endl;
   some_file.close();    
  }
 
 
 // Check things
 //-------------
 problem_pt->check_halo_schemes(quiet_doc_info);
 
 
 // Solve?
 //-------
 if (solve)
  {
 
#ifdef OOMPH_HAS_TRILINOS
 
  // Create a Trilinos solver
  TrilinosAztecOOSolver* linear_solver_pt = new TrilinosAztecOOSolver;
 
  // Create the Trilinos ML preconditioner
  TrilinosMLPreconditioner* preconditioner_pt = new TrilinosMLPreconditioner;
 
  // Set the preconditioner pointer
  linear_solver_pt->preconditioner_pt() = preconditioner_pt;
 
  // Choose solver type (GMRES)
  linear_solver_pt->solver_type()=TrilinosAztecOOSolver::GMRES;
 
  // Set linear solver
  problem_pt->linear_solver_pt() = linear_solver_pt;
  linear_solver_pt->disable_doc_time();
    
#endif
 
   // Solve
   problem_pt->newton_solve();
   
   if (doc)
    {
     // Doc the solution
     problem_pt->doc_solution(doc_info);
    }
  }
 
 
   // distributed... so now check halo schemes, 
   // as there seems to be a problem after the uniform refinement...
   // andy trying to find bug
//   DocInfo tmp_doc_info;
 
   // Check files into DIST directory
//   tmp_doc_info.set_directory("DIST");
//   problem_pt->check_halo_schemes(tmp_doc_info); // andy checking for bug 
//   problem_pt->mesh_pt()->doc_mesh_distribution(tmp_doc_info);
 
 
 // Now refine the distributed problem uniformly
 //---------------------------------------------
 for (unsigned i=0;i<n_refine_once_distributed;i++)
  {
   // Number of dofs before uniform refinement
   n_dof_orig=problem_pt->ndof();
 
   // Refine uniformly
   //-----------------
   t_start = clock();   
   problem_pt->refine_uniformly();
   t_end = clock();
 
 
   // Assess efficiency
   //------------------
   problem_pt->mesh_pt()->get_halo_node_stats(av_number_halo_nodes,
                                              max_number_halo_nodes,
                                              min_number_halo_nodes);
   problem_pt->mesh_pt()->get_efficiency_of_mesh_distribution
    (av_efficiency,max_efficiency,min_efficiency);
 
 
   if (my_rank==0)
    {
     double t_tot=double(t_end-t_start)/CLOCKS_PER_SEC;
     oomph_info << "Time for distributed adaptation [ndof=" 
                << problem_pt->ndof() 
                << "]: " << t_tot << std::endl;
 
     // Number of dofs achieved after uniform refinement
     n_dof=problem_pt->ndof();
 
     // Snippet indicating if  redistribution of haloes is used or not
     string snippet="";
     if (!redistribute) snippet="out";
 
     sprintf(
      filename,
      "%s/three_d_distributed_refinement_np_%i_ninitialrefinement_%i_ntotalrefinement_%i_ndof_%i_with%sredistribution.dat",
      doc_info.directory().c_str(),
      n_proc,
      n_refine_first,
      n_refine_first+i+1,
      n_dof,
      snippet.c_str());
     some_file.open(filename);
     some_file << n_proc << " " 
               << n_dof_orig << " " 
               << n_dof << " " 
               << t_tot << " "
               << av_number_halo_nodes << " " 
               << max_number_halo_nodes << " " 
               << min_number_halo_nodes  << " " 
               << av_efficiency << " " 
               << max_efficiency << " " 
               << min_efficiency << " " 
               << std::endl;
     some_file.close();   
 
    }
 
 
 
   // Redistribute (if required)
   //---------------------------
   if (redistribute)
    {
      t_start = clock();   
      bool report_stats=false;
      doc_info.disable_doc();
      problem_pt->prune_halo_elements_and_nodes(doc_info, 
                                                report_stats);
      t_end = clock();
      double t_tot=double(t_end-t_start)/CLOCKS_PER_SEC;
     
 
      if (my_rank==0)
       {
 
        oomph_info << "Av., min. max. efficiency before redistribution " 
                   << av_efficiency << " " 
                   << min_efficiency << " " 
                   << max_efficiency << " " 
                   << std::endl;
        sprintf(
         filename,
         "%s/three_d_redistribution_np_%i_ninitialrefinement_%i_ntotalrefinement_%i_ndof_%i.dat",
         doc_info.directory().c_str(),
         n_proc,
         n_refine_first,
         n_refine_first+i+1,
         n_dof);
        some_file.open(filename);
        some_file << n_proc << " " 
                  << n_dof << " " 
                  << t_tot << " "
                  << av_number_halo_nodes << " " 
                  << max_number_halo_nodes << " " 
                  << min_number_halo_nodes  << " " 
                  << av_efficiency << " " 
                  << max_efficiency << " " 
                  << min_efficiency << " " 
                  << std::endl;
       }
 
 
      // Re-assess efficiency
      //----------------------
      problem_pt->mesh_pt()->get_halo_node_stats(av_number_halo_nodes,
                                                 max_number_halo_nodes,
                                                 min_number_halo_nodes);
      problem_pt->mesh_pt()->
       get_efficiency_of_mesh_distribution(av_efficiency,
                                           max_efficiency,
                                           min_efficiency);
      
      
      if (my_rank==0)
       {
        oomph_info << "Time for redistribution  [ndof=" 
                   << problem_pt->ndof() 
                   << "]: " << t_tot << std::endl;
 
 
        some_file << n_proc << " " 
                  << n_dof << " " 
                  << t_tot << " "
                  << av_number_halo_nodes << " " 
                  << max_number_halo_nodes << " " 
                  << min_number_halo_nodes  << " " 
                  << av_efficiency << " " 
                  << max_efficiency << " " 
                  << min_efficiency << " " 
                  << std::endl;
        some_file.close();   
       }
            
      problem_pt->mesh_pt()->
       get_efficiency_of_mesh_distribution(av_efficiency,
                                           max_efficiency,
                                           min_efficiency);
      oomph_info << "Av., min., max. efficiency after redistribution " 
                 << av_efficiency << " " 
                 << min_efficiency << " " 
                 << max_efficiency << " " 
                 << std::endl;
    }
   
     
   // Check things
   //-------------
   problem_pt->check_halo_schemes(quiet_doc_info);
 
 
   // Solve again
   //------------
   if (solve)
    {
     // Solve again
     problem_pt->newton_solve();
     
     if (doc)
      {
       // Doc the solution
       problem_pt->doc_solution(doc_info);
      }
    }
 
 }
 
 
 // Unrefine the distributed problem uniformly (only if not redistributed)
 //-----------------------------------------------------------------------
 if (!redistribute)
  {
   // Number of dofs before unrefinement
   n_dof_orig=problem_pt->ndof();
   
   
   t_start = clock();
   problem_pt->unrefine_uniformly();
   t_end = clock();
   
   if (my_rank==0)
    {
     double t_tot=double(t_end-t_start)/CLOCKS_PER_SEC;
     oomph_info << "Time for distributed unrefinement [ndof=" 
                << problem_pt->ndof() 
                << "]: " << t_tot << std::endl;
     
     // Number of dofs after unrefinement
     n_dof=problem_pt->ndof();
     
     sprintf(
      filename,
      "%s/three_d_distributed_unrefinement_np_%i_ninitialrefinement_%i_ntotalrefinement_%i_ndof_%i.dat",
      doc_info.directory().c_str(),
      n_proc,
      n_refine_first,
      n_refine_first+n_refine_once_distributed,
      n_dof);
     some_file.open(filename);
     some_file << n_proc << " " 
               << n_dof_orig << " " 
               << n_dof << " " 
               << t_tot << " "
               << av_number_halo_nodes << " " 
               << max_number_halo_nodes << " " 
               << min_number_halo_nodes  << " " 
               << av_efficiency << " " 
               << max_efficiency << " " 
               << min_efficiency << " " 
               << std::endl;
     some_file.close();   
     
    }
   
      
   // Check things
   problem_pt->check_halo_schemes(quiet_doc_info);
   
   // Solve again
   if (solve)
    {
     // Solve again
     problem_pt->newton_solve();
     
     
     if (doc)
      {
       // Doc the solution
       problem_pt->doc_solution(doc_info);
      }
    }
  }
 
}

References oomph::Problem::communicator_pt(), oomph::DocInfo::directory(), oomph::DocInfo::disable_doc(), oomph::LinearSolver::disable_doc_time(), ThreeDPoissonProblem< ELEMENT >::doc_solution(), e(), MergeRestartFiles::filename, TanhSolnForPoisson::get_source(), i, oomph::Problem::linear_solver_pt(), ThreeDPoissonProblem< ELEMENT >::mesh_pt(), oomph::OomphCommunicator::my_rank(), oomph::Problem::ndof(), oomph::Problem::newton_solve(), oomph::OomphCommunicator::nproc(), oomph::oomph_info, oomph::IterativeLinearSolver::preconditioner_pt(), oomph::Problem::refine_uniformly(), oomph::DocInfo::set_directory(), solve, oomph::TrilinosAztecOOSolver::solver_type(), oomph::Global_string_for_annotation::string(), and oomph::Problem::unrefine_uniformly().

Referenced by main().