Example 18 - Unsteady Navier-Stokes Equations with DiffSystem



This example shows how the transient nonlinear problem from example 13 can be solved using the new (and experimental) DiffSystem class framework

Basic include files
        #include "equation_systems.h"
        #include "error_vector.h"
        #include "getpot.h"
        #include "exodusII_io.h"
        #include "kelly_error_estimator.h"
        #include "mesh.h"
        #include "mesh_generation.h"
        #include "mesh_refinement.h"
        #include "uniform_refinement_estimator.h"
        
Some (older) compilers do not offer full stream functionality, OStringStream works around this.
        #include "o_string_stream.h"
        
The systems and solvers we may use
        #include "naviersystem.h"
        #include "diff_solver.h"
        #include "euler_solver.h"
        #include "steady_solver.h"
        
Bring in everything from the libMesh namespace
        using namespace libMesh;
        
The main program.
        int main (int argc, char** argv)
        {
Initialize libMesh.
          LibMeshInit init (argc, argv);
        
This example fails without at least double precision FP
        #ifdef LIBMESH_DEFAULT_SINGLE_PRECISION
          libmesh_example_assert(false, "--disable-singleprecision");
        #endif
        
        #ifndef LIBMESH_ENABLE_AMR
          libmesh_example_assert(false, "--enable-amr");
        #else
        
Trilinos solver NaNs by default on the zero pressure block. We'll skip this example for now.
          if (libMesh::default_solver_package() == TRILINOS_SOLVERS)
            {
              std::cout << "We skip example 18 when using the Trilinos solvers.\n"
                        << std::endl;
              return 0;
            }
        
Parse the input file
          GetPot infile("ex18.in");
        
Read in parameters from the input file
          const Real global_tolerance          = infile("global_tolerance", 0.);
          const unsigned int nelem_target      = infile("n_elements", 400);
          const bool transient                 = infile("transient", true);
          const Real deltat                    = infile("deltat", 0.005);
          unsigned int n_timesteps             = infile("n_timesteps", 20);
          const unsigned int write_interval    = infile("write_interval", 5);
          const unsigned int coarsegridsize    = infile("coarsegridsize", 1);
          const unsigned int coarserefinements = infile("coarserefinements", 0);
          const unsigned int max_adaptivesteps = infile("max_adaptivesteps", 10);
          const unsigned int dim               = infile("dimension", 2);
        
Skip higher-dimensional examples on a lower-dimensional libMesh build
          libmesh_example_assert(dim <= LIBMESH_DIM, "2D/3D support");
          
We have only defined 2 and 3 dimensional problems
          libmesh_assert (dim == 2 || dim == 3);
        
Create a mesh.
          Mesh mesh;
          
And an object to refine it
          MeshRefinement mesh_refinement(mesh);
          mesh_refinement.coarsen_by_parents() = true;
          mesh_refinement.absolute_global_tolerance() = global_tolerance;
          mesh_refinement.nelem_target() = nelem_target;
          mesh_refinement.refine_fraction() = 0.3;
          mesh_refinement.coarsen_fraction() = 0.3;
          mesh_refinement.coarsen_threshold() = 0.1;
        
Use the MeshTools::Generation mesh generator to create a uniform grid on the square [-1,1]^D. We instruct the mesh generator to build a mesh of 8x8 \p Quad9 elements in 2D, or \p Hex27 elements in 3D. Building these higher-order elements allows us to use higher-order approximation, as in example 3.
          if (dim == 2)
            MeshTools::Generation::build_square (mesh,
                                                 coarsegridsize,
                                                 coarsegridsize,
                                                 0., 1.,
                                                 0., 1.,
                                                 QUAD9);
          else if (dim == 3)
            MeshTools::Generation::build_cube (mesh,
                                               coarsegridsize,
                                               coarsegridsize,
                                               coarsegridsize,
                                               0., 1.,
                                               0., 1.,
                                               0., 1.,
                                               HEX27);
        
          mesh_refinement.uniformly_refine(coarserefinements);
        
Print information about the mesh to the screen.
          mesh.print_info();
        
Create an equation systems object.
          EquationSystems equation_systems (mesh);
        
Declare the system "Navier-Stokes" and its variables.
          NavierSystem & system = 
            equation_systems.add_system<NavierSystem> ("Navier-Stokes");
        
Solve this as a time-dependent or steady system
          if (transient)
            system.time_solver =
              AutoPtr<TimeSolver>(new EulerSolver(system));
          else
            {
              system.time_solver =
                AutoPtr<TimeSolver>(new SteadySolver(system));
              libmesh_assert(n_timesteps == 1);
            }
        
Initialize the system
          equation_systems.init ();
        
Set the time stepping options
          system.deltat = deltat;
        
And the nonlinear solver options
          DiffSolver &solver = *(system.time_solver->diff_solver().get());
          solver.quiet = infile("solver_quiet", true);
          solver.max_nonlinear_iterations =
            infile("max_nonlinear_iterations", 15);
          solver.relative_step_tolerance =
            infile("relative_step_tolerance", 1.e-3);
          solver.relative_residual_tolerance =
            infile("relative_residual_tolerance", 0.0);
          solver.absolute_residual_tolerance =
            infile("absolute_residual_tolerance", 0.0);
            
And the linear solver options
          solver.max_linear_iterations =
            infile("max_linear_iterations", 50000);
          solver.initial_linear_tolerance =
            infile("initial_linear_tolerance", 1.e-3);
        
Print information about the system to the screen.
          equation_systems.print_info();
        
Now we begin the timestep loop to compute the time-accurate solution of the equations.
          for (unsigned int t_step=0; t_step != n_timesteps; ++t_step)
            {
A pretty update message
              std::cout << "\n\nSolving time step " << t_step << ", time = "
                        << system.time << std::endl;
        
Adaptively solve the timestep
              unsigned int a_step = 0;
              for (; a_step != max_adaptivesteps; ++a_step)
                {
                  system.solve();
        
                  ErrorVector error;
        
                  AutoPtr<ErrorEstimator> error_estimator;
        
To solve to a tolerance in this problem we need a better estimator than Kelly
                  if (global_tolerance != 0.)
                    {
We can't adapt to both a tolerance and a mesh size at once
                      libmesh_assert (nelem_target == 0);
        
                      UniformRefinementEstimator *u =
                        new UniformRefinementEstimator;
        
The lid-driven cavity problem isn't in H1, so lets estimate L2 error
                      u->error_norm = L2;
        
                      error_estimator.reset(u);
                    }
                  else
                    {
If we aren't adapting to a tolerance we need a target mesh size
                      libmesh_assert (nelem_target > 0);
        
Kelly is a lousy estimator to use for a problem not in H1 - if we were doing more than a few timesteps we'd need to turn off or limit the maximum level of our adaptivity eventually
                      error_estimator.reset(new KellyErrorEstimator);
                    }
        
Calculate error based on u and v (and w?) but not p
                  std::vector<Real> weights(2,1.0);  // u, v
                  if (dim == 3)
                    weights.push_back(1.0);          // w
                  weights.push_back(0.0);            // p
Keep the same default norm type.
                  std::vector<FEMNormType>
        	    norms(1, error_estimator->error_norm.type(0));
        	  error_estimator->error_norm = SystemNorm(norms, weights);
        
                  error_estimator->estimate_error(system, error);
        
Print out status at each adaptive step.
                  Real global_error = error.l2_norm();
                  std::cout << "adaptive step " << a_step << ": ";
                  if (global_tolerance != 0.)
                    std::cout << "global_error = " << global_error
                              << " with ";
                  std::cout << mesh.n_active_elem()
                            << " active elements and "
                            << equation_systems.n_active_dofs()
                            << " active dofs." << std::endl;
                  if (global_tolerance != 0.)
                    std::cout << "worst element error = " << error.maximum()
                              << ", mean = " << error.mean() << std::endl;
        
                  if (global_tolerance != 0.)
                    {
If we've reached our desired tolerance, we don't need any more adaptive steps
                      if (global_error < global_tolerance)
                        break;
                      mesh_refinement.flag_elements_by_error_tolerance(error);
                    }
                  else
                    {
If flag_elements_by_nelem_target returns true, this should be our last adaptive step.
                      if (mesh_refinement.flag_elements_by_nelem_target(error))
                        {
                          mesh_refinement.refine_and_coarsen_elements();
                          equation_systems.reinit();
                          a_step = max_adaptivesteps;
                          break;
                        }
                    }
        
Carry out the adaptive mesh refinement/coarsening
                  mesh_refinement.refine_and_coarsen_elements();
                  equation_systems.reinit();
                }
Do one last solve if necessary
              if (a_step == max_adaptivesteps)
                {
                  system.solve();
                }
        
Advance to the next timestep in a transient problem
              system.time_solver->advance_timestep();
        
        #ifdef LIBMESH_HAVE_EXODUS_API
Write out this timestep if we're requested to
              if ((t_step+1)%write_interval == 0)
                {
                  OStringStream file_name;
        
We write the file in the ExodusII format.
                  file_name << "out_";
                  OSSRealzeroright(file_name,3,0, t_step + 1);
                  file_name << ".exd";
        
                  ExodusII_IO(mesh).write_equation_systems (file_name.str(),
                                                      equation_systems);
                }
        #endif // #ifdef LIBMESH_HAVE_EXODUS_API
            }
        #endif // #ifndef LIBMESH_ENABLE_AMR
          
All done.
          return 0;
        }



The program without comments:

 
  
  #include "equation_systems.h"
  #include "error_vector.h"
  #include "getpot.h"
  #include "exodusII_io.h"
  #include "kelly_error_estimator.h"
  #include "mesh.h"
  #include "mesh_generation.h"
  #include "mesh_refinement.h"
  #include "uniform_refinement_estimator.h"
  
  #include "o_string_stream.h"
  
  #include "naviersystem.h"
  #include "diff_solver.h"
  #include "euler_solver.h"
  #include "steady_solver.h"
  
  using namespace libMesh;
  
  int main (int argc, char** argv)
  {
    LibMeshInit init (argc, argv);
  
  #ifdef LIBMESH_DEFAULT_SINGLE_PRECISION
    libmesh_example_assert(false, "--disable-singleprecision");
  #endif
  
  #ifndef LIBMESH_ENABLE_AMR
    libmesh_example_assert(false, "--enable-amr");
  #else
  
    if (libMesh::default_solver_package() == TRILINOS_SOLVERS)
      {
        std::cout << "We skip example 18 when using the Trilinos solvers.\n"
                  << std::endl;
        return 0;
      }
  
    GetPot infile("ex18.in");
  
    const Real global_tolerance          = infile("global_tolerance", 0.);
    const unsigned int nelem_target      = infile("n_elements", 400);
    const bool transient                 = infile("transient", true);
    const Real deltat                    = infile("deltat", 0.005);
    unsigned int n_timesteps             = infile("n_timesteps", 20);
    const unsigned int write_interval    = infile("write_interval", 5);
    const unsigned int coarsegridsize    = infile("coarsegridsize", 1);
    const unsigned int coarserefinements = infile("coarserefinements", 0);
    const unsigned int max_adaptivesteps = infile("max_adaptivesteps", 10);
    const unsigned int dim               = infile("dimension", 2);
  
    libmesh_example_assert(dim <= LIBMESH_DIM, "2D/3D support");
    
    libmesh_assert (dim == 2 || dim == 3);
  
    Mesh mesh;
    
    MeshRefinement mesh_refinement(mesh);
    mesh_refinement.coarsen_by_parents() = true;
    mesh_refinement.absolute_global_tolerance() = global_tolerance;
    mesh_refinement.nelem_target() = nelem_target;
    mesh_refinement.refine_fraction() = 0.3;
    mesh_refinement.coarsen_fraction() = 0.3;
    mesh_refinement.coarsen_threshold() = 0.1;
  
    if (dim == 2)
      MeshTools::Generation::build_square (mesh,
                                           coarsegridsize,
                                           coarsegridsize,
                                           0., 1.,
                                           0., 1.,
                                           QUAD9);
    else if (dim == 3)
      MeshTools::Generation::build_cube (mesh,
                                         coarsegridsize,
                                         coarsegridsize,
                                         coarsegridsize,
                                         0., 1.,
                                         0., 1.,
                                         0., 1.,
                                         HEX27);
  
    mesh_refinement.uniformly_refine(coarserefinements);
  
    mesh.print_info();
  
    EquationSystems equation_systems (mesh);
  
    NavierSystem & system = 
      equation_systems.add_system<NavierSystem> ("Navier-Stokes");
  
    if (transient)
      system.time_solver =
        AutoPtr<TimeSolver>(new EulerSolver(system));
    else
      {
        system.time_solver =
          AutoPtr<TimeSolver>(new SteadySolver(system));
        libmesh_assert(n_timesteps == 1);
      }
  
    equation_systems.init ();
  
    system.deltat = deltat;
  
    DiffSolver &solver = *(system.time_solver->diff_solver().get());
    solver.quiet = infile("solver_quiet", true);
    solver.max_nonlinear_iterations =
      infile("max_nonlinear_iterations", 15);
    solver.relative_step_tolerance =
      infile("relative_step_tolerance", 1.e-3);
    solver.relative_residual_tolerance =
      infile("relative_residual_tolerance", 0.0);
    solver.absolute_residual_tolerance =
      infile("absolute_residual_tolerance", 0.0);
      
    solver.max_linear_iterations =
      infile("max_linear_iterations", 50000);
    solver.initial_linear_tolerance =
      infile("initial_linear_tolerance", 1.e-3);
  
    equation_systems.print_info();
  
    for (unsigned int t_step=0; t_step != n_timesteps; ++t_step)
      {
        std::cout << "\n\nSolving time step " << t_step << ", time = "
                  << system.time << std::endl;
  
        unsigned int a_step = 0;
        for (; a_step != max_adaptivesteps; ++a_step)
          {
            system.solve();
  
            ErrorVector error;
  
            AutoPtr<ErrorEstimator> error_estimator;
  
            if (global_tolerance != 0.)
              {
                libmesh_assert (nelem_target == 0);
  
                UniformRefinementEstimator *u =
                  new UniformRefinementEstimator;
  
                u->error_norm = L2;
  
                error_estimator.reset(u);
              }
            else
              {
                libmesh_assert (nelem_target > 0);
  
                error_estimator.reset(new KellyErrorEstimator);
              }
  
  	  std::vector<Real> weights(2,1.0);  // u, v
            if (dim == 3)
              weights.push_back(1.0);          // w
            weights.push_back(0.0);            // p
  	  std::vector<FEMNormType>
  	    norms(1, error_estimator->error_norm.type(0));
  	  error_estimator->error_norm = SystemNorm(norms, weights);
  
            error_estimator->estimate_error(system, error);
  
            Real global_error = error.l2_norm();
            std::cout << "adaptive step " << a_step << ": ";
            if (global_tolerance != 0.)
              std::cout << "global_error = " << global_error
                        << " with ";
            std::cout << mesh.n_active_elem()
                      << " active elements and "
                      << equation_systems.n_active_dofs()
                      << " active dofs." << std::endl;
            if (global_tolerance != 0.)
              std::cout << "worst element error = " << error.maximum()
                        << ", mean = " << error.mean() << std::endl;
  
            if (global_tolerance != 0.)
              {
                if (global_error < global_tolerance)
                  break;
                mesh_refinement.flag_elements_by_error_tolerance(error);
              }
            else
              {
                if (mesh_refinement.flag_elements_by_nelem_target(error))
                  {
                    mesh_refinement.refine_and_coarsen_elements();
                    equation_systems.reinit();
                    a_step = max_adaptivesteps;
                    break;
                  }
              }
  
            mesh_refinement.refine_and_coarsen_elements();
            equation_systems.reinit();
          }
        if (a_step == max_adaptivesteps)
          {
            system.solve();
          }
  
        system.time_solver->advance_timestep();
  
  #ifdef LIBMESH_HAVE_EXODUS_API
        if ((t_step+1)%write_interval == 0)
          {
            OStringStream file_name;
  
            file_name << "out_";
            OSSRealzeroright(file_name,3,0, t_step + 1);
            file_name << ".exd";
  
            ExodusII_IO(mesh).write_equation_systems (file_name.str(),
                                                equation_systems);
          }
  #endif // #ifdef LIBMESH_HAVE_EXODUS_API
      }
  #endif // #ifndef LIBMESH_ENABLE_AMR
    
    return 0;
  }



The console output of the program:

Compiling C++ (in optimized mode) ex18.C...
Compiling C++ (in optimized mode) naviersystem.C...
Linking ex18-opt...
***************************************************************
* Running Example  mpirun -np 2 ./ex18-opt -pc_type bjacobi -sub_pc_type ilu -sub_pc_factor_levels 4 -sub_pc_factor_zeropivot 0 -ksp_right_pc -log_summary
***************************************************************
 
 Mesh Information:
  mesh_dimension()=2
  spatial_dimension()=3
  n_nodes()=1681
    n_local_nodes()=863
  n_elem()=400
    n_local_elem()=200
    n_active_elem()=400
  n_subdomains()=1
  n_processors()=2
  processor_id()=0

 EquationSystems
  n_systems()=1
   System "Navier-Stokes"
    Type "Implicit"
    Variables="u" "v" "p" 
    Finite Element Types="LAGRANGE" "LAGRANGE" "LAGRANGE" 
    Approximation Orders="SECOND" "SECOND" "FIRST" 
    n_dofs()=3803
    n_local_dofs()=1958
    n_constrained_dofs()=0
    n_vectors()=2



Solving time step 0, time = 0


Solving time step 1, time = 0.005


Solving time step 2, time = 0.01


Solving time step 3, time = 0.015


Solving time step 4, time = 0.02


Solving time step 5, time = 0.025


Solving time step 6, time = 0.03


Solving time step 7, time = 0.035


Solving time step 8, time = 0.04


Solving time step 9, time = 0.045


Solving time step 10, time = 0.05


Solving time step 11, time = 0.055


Solving time step 12, time = 0.06


Solving time step 13, time = 0.065


Solving time step 14, time = 0.07
************************************************************************************************************************
***             WIDEN YOUR WINDOW TO 120 CHARACTERS.  Use 'enscript -r -fCourier9' to print this document            ***
************************************************************************************************************************

---------------------------------------------- PETSc Performance Summary: ----------------------------------------------

./ex18-opt on a gcc-4.5-l named daedalus with 2 processors, by roystgnr Tue Feb 22 12:21:48 2011
Using Petsc Release Version 3.1.0, Patch 5, Mon Sep 27 11:51:54 CDT 2010

                         Max       Max/Min        Avg      Total 
Time (sec):           1.322e+01      1.00007   1.322e+01
Objects:              4.840e+02      1.00000   4.840e+02
Flops:                8.464e+09      1.14869   7.916e+09  1.583e+10
Flops/sec:            6.403e+08      1.14877   5.988e+08  1.198e+09
MPI Messages:         2.178e+03      1.00000   2.178e+03  4.356e+03
MPI Message Lengths:  5.076e+06      1.01961   2.308e+03  1.005e+07
MPI Reductions:       4.729e+03      1.00000

Flop counting convention: 1 flop = 1 real number operation of type (multiply/divide/add/subtract)
                            e.g., VecAXPY() for real vectors of length N --> 2N flops
                            and VecAXPY() for complex vectors of length N --> 8N flops

Summary of Stages:   ----- Time ------  ----- Flops -----  --- Messages ---  -- Message Lengths --  -- Reductions --
                        Avg     %Total     Avg     %Total   counts   %Total     Avg         %Total   counts   %Total 
 0:      Main Stage: 1.3218e+01 100.0%  1.5832e+10 100.0%  4.356e+03 100.0%  2.308e+03      100.0%  4.632e+03  97.9% 

------------------------------------------------------------------------------------------------------------------------
See the 'Profiling' chapter of the users' manual for details on interpreting output.
Phase summary info:
   Count: number of times phase was executed
   Time and Flops: Max - maximum over all processors
                   Ratio - ratio of maximum to minimum over all processors
   Mess: number of messages sent
   Avg. len: average message length
   Reduct: number of global reductions
   Global: entire computation
   Stage: stages of a computation. Set stages with PetscLogStagePush() and PetscLogStagePop().
      %T - percent time in this phase         %F - percent flops in this phase
      %M - percent messages in this phase     %L - percent message lengths in this phase
      %R - percent reductions in this phase
   Total Mflop/s: 10e-6 * (sum of flops over all processors)/(max time over all processors)
------------------------------------------------------------------------------------------------------------------------
Event                Count      Time (sec)     Flops                             --- Global ---  --- Stage ---   Total
                   Max Ratio  Max     Ratio   Max  Ratio  Mess   Avg len Reduct  %T %F %M %L %R  %T %F %M %L %R Mflop/s
------------------------------------------------------------------------------------------------------------------------

--- Event Stage 0: Main Stage

KSPGMRESOrthog      1454 1.0 1.9099e-01 3.1 1.64e+08 1.1 0.0e+00 0.0e+00 1.5e+03  1  2  0  0 31   1  2  0  0 31  1672
KSPSetup              52 1.0 4.7922e-05 1.2 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00  0  0  0  0  0   0  0  0  0  0     0
KSPSolve              26 1.0 1.2732e+01 1.0 8.46e+09 1.1 3.0e+03 1.3e+03 3.1e+03 96100 69 40 65  96100 69 40 66  1243
PCSetUp               52 1.0 1.0410e+01 1.1 5.54e+09 1.2 0.0e+00 0.0e+00 7.8e+01 74 65  0  0  2  74 65  0  0  2   984
PCSetUpOnBlocks       26 1.0 1.0410e+01 1.1 5.54e+09 1.2 0.0e+00 0.0e+00 7.8e+01 74 65  0  0  2  74 65  0  0  2   984
PCApply             1513 1.0 2.0386e+00 1.1 2.52e+09 1.1 0.0e+00 0.0e+00 0.0e+00 15 30  0  0  0  15 30  0  0  0  2351
VecMDot             1454 1.0 1.6936e-01 4.3 8.22e+07 1.1 0.0e+00 0.0e+00 1.5e+03  1  1  0  0 31   1  1  0  0 31   943
VecNorm             1654 1.0 3.5646e-02 1.2 6.48e+06 1.1 0.0e+00 0.0e+00 1.7e+03  0  0  0  0 35   0  0  0  0 36   353
VecScale            1513 1.0 2.4424e-03 1.1 2.96e+06 1.1 0.0e+00 0.0e+00 0.0e+00  0  0  0  0  0   0  0  0  0  0  2356
VecCopy              193 1.0 4.4966e-04 1.3 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00  0  0  0  0  0   0  0  0  0  0     0
VecSet              1684 1.0 2.1021e-03 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00  0  0  0  0  0   0  0  0  0  0     0
VecAXPY              144 1.0 4.9233e-04 1.3 5.64e+05 1.1 0.0e+00 0.0e+00 0.0e+00  0  0  0  0  0   0  0  0  0  0  2225
VecMAXPY            1513 1.0 2.2277e-02 1.1 8.79e+07 1.1 0.0e+00 0.0e+00 0.0e+00  0  1  0  0  0   0  1  0  0  0  7667
VecAssemblyBegin     283 1.0 8.8170e-03 3.1 0.00e+00 0.0 1.0e+02 1.1e+03 8.5e+02  0  0  2  1 18   0  0  2  1 18     0
VecAssemblyEnd       283 1.0 1.8954e-04 1.1 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00  0  0  0  0  0   0  0  0  0  0     0
VecScatterBegin     1729 1.0 1.2525e-02 1.9 0.00e+00 0.0 3.2e+03 1.4e+03 0.0e+00  0  0 74 44  0   0  0 74 44  0     0
VecScatterEnd       1729 1.0 1.3028e+00292.7 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00  5  0  0  0  0   5  0  0  0  0     0
VecNormalize        1513 1.0 3.4015e-02 1.2 8.89e+06 1.1 0.0e+00 0.0e+00 1.5e+03  0  0  0  0 32   0  0  0  0 33   507
MatMult             1513 1.0 1.5012e+00 8.1 2.30e+08 1.1 3.0e+03 1.3e+03 0.0e+00  6  3 69 40  0   6  3 69 40  0   295
MatSolve            1513 1.0 2.0200e+00 1.1 2.52e+09 1.1 0.0e+00 0.0e+00 0.0e+00 15 30  0  0  0  15 30  0  0  0  2373
MatLUFactorNum        26 1.0 4.1903e+00 1.2 5.54e+09 1.2 0.0e+00 0.0e+00 0.0e+00 29 65  0  0  0  29 65  0  0  0  2445
MatILUFactorSym       26 1.0 6.2176e+00 1.1 0.00e+00 0.0 0.0e+00 0.0e+00 2.6e+01 44  0  0  0  1  44  0  0  0  1     0
MatAssemblyBegin      52 1.0 2.0273e-03 2.9 0.00e+00 0.0 7.8e+01 2.2e+04 1.0e+02  0  0  2 17  2   0  0  2 17  2     0
MatAssemblyEnd        52 1.0 9.5301e-03 1.2 0.00e+00 0.0 4.0e+00 3.3e+02 5.8e+01  0  0  0  0  1   0  0  0  0  1     0
MatGetRowIJ           26 1.0 7.6294e-06 1.3 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00  0  0  0  0  0   0  0  0  0  0     0
MatGetOrdering        26 1.0 3.6812e-04 1.0 0.00e+00 0.0 0.0e+00 0.0e+00 5.2e+01  0  0  0  0  1   0  0  0  0  1     0
MatZeroEntries        28 1.0 2.0542e-03 1.1 0.00e+00 0.0 0.0e+00 0.0e+00 0.0e+00  0  0  0  0  0   0  0  0  0  0     0
------------------------------------------------------------------------------------------------------------------------

Memory usage is given in bytes:

Object Type          Creations   Destructions     Memory  Descendants' Mem.
Reports information only for process 0.

--- Event Stage 0: Main Stage

       Krylov Solver     3              3        19872     0
      Preconditioner     3              3         2048     0
                 Vec    88             88      1211040     0
         Vec Scatter   125            125       108500     0
           Index Set   220            220       739584     0
   IS L to G Mapping    16             16       145344     0
              Matrix    29             29    260814268     0
========================================================================================================================
Average time to get PetscTime(): 1.19209e-07
Average time for MPI_Barrier(): 1.62125e-06
Average time for zero size MPI_Send(): 5.96046e-06
#PETSc Option Table entries:
-ksp_right_pc
-log_summary
-pc_type bjacobi
-sub_pc_factor_levels 4
-sub_pc_factor_zeropivot 0
-sub_pc_type ilu
#End of PETSc Option Table entries
Compiled without FORTRAN kernels
Compiled with full precision matrices (default)
sizeof(short) 2 sizeof(int) 4 sizeof(long) 8 sizeof(void*) 8 sizeof(PetscScalar) 8
Configure run at: Fri Oct 15 13:01:23 2010
Configure options: --with-debugging=false --COPTFLAGS=-O3 --CXXOPTFLAGS=-O3 --FOPTFLAGS=-O3 --with-clanguage=C++ --with-shared=1 --with-mpi-dir=/org/centers/pecos/LIBRARIES/MPICH2/mpich2-1.2.1-gcc-4.5-lucid --with-mumps=true --download-mumps=ifneeded --with-parmetis=true --download-parmetis=ifneeded --with-superlu=true --download-superlu=ifneeded --with-superludir=true --download-superlu_dist=ifneeded --with-blacs=true --download-blacs=ifneeded --with-scalapack=true --download-scalapack=ifneeded --with-hypre=true --download-hypre=ifneeded --with-blas-lib="[/org/centers/pecos/LIBRARIES/MKL/mkl-10.0.3.020-gcc-4.5-lucid/lib/em64t/libmkl_intel_lp64.so,/org/centers/pecos/LIBRARIES/MKL/mkl-10.0.3.020-gcc-4.5-lucid/lib/em64t/libmkl_sequential.so,/org/centers/pecos/LIBRARIES/MKL/mkl-10.0.3.020-gcc-4.5-lucid/lib/em64t/libmkl_core.so]" --with-lapack-lib=/org/centers/pecos/LIBRARIES/MKL/mkl-10.0.3.020-gcc-4.5-lucid/lib/em64t/libmkl_solver_lp64_sequential.a
-----------------------------------------
Libraries compiled on Fri Oct 15 13:01:23 CDT 2010 on atreides 
Machine characteristics: Linux atreides 2.6.32-25-generic #44-Ubuntu SMP Fri Sep 17 20:05:27 UTC 2010 x86_64 GNU/Linux 
Using PETSc directory: /org/centers/pecos/LIBRARIES/PETSC3/petsc-3.1-p5
Using PETSc arch: gcc-4.5-lucid-mpich2-1.2.1-cxx-opt
-----------------------------------------
Using C compiler: /org/centers/pecos/LIBRARIES/MPICH2/mpich2-1.2.1-gcc-4.5-lucid/bin/mpicxx -Wall -Wwrite-strings -Wno-strict-aliasing -O3   -fPIC   
Using Fortran compiler: /org/centers/pecos/LIBRARIES/MPICH2/mpich2-1.2.1-gcc-4.5-lucid/bin/mpif90 -fPIC -Wall -Wno-unused-variable -O3    
-----------------------------------------
Using include paths: -I/org/centers/pecos/LIBRARIES/PETSC3/petsc-3.1-p5/gcc-4.5-lucid-mpich2-1.2.1-cxx-opt/include -I/org/centers/pecos/LIBRARIES/PETSC3/petsc-3.1-p5/include -I/org/centers/pecos/LIBRARIES/PETSC3/petsc-3.1-p5/gcc-4.5-lucid-mpich2-1.2.1-cxx-opt/include -I/org/centers/pecos/LIBRARIES/MPICH2/mpich2-1.2.1-gcc-4.5-lucid/include  
------------------------------------------
Using C linker: /org/centers/pecos/LIBRARIES/MPICH2/mpich2-1.2.1-gcc-4.5-lucid/bin/mpicxx -Wall -Wwrite-strings -Wno-strict-aliasing -O3 
Using Fortran linker: /org/centers/pecos/LIBRARIES/MPICH2/mpich2-1.2.1-gcc-4.5-lucid/bin/mpif90 -fPIC -Wall -Wno-unused-variable -O3  
Using libraries: -Wl,-rpath,/org/centers/pecos/LIBRARIES/PETSC3/petsc-3.1-p5/gcc-4.5-lucid-mpich2-1.2.1-cxx-opt/lib -L/org/centers/pecos/LIBRARIES/PETSC3/petsc-3.1-p5/gcc-4.5-lucid-mpich2-1.2.1-cxx-opt/lib -lpetsc       -lX11 -Wl,-rpath,/org/centers/pecos/LIBRARIES/PETSC3/petsc-3.1-p5/gcc-4.5-lucid-mpich2-1.2.1-cxx-opt/lib -L/org/centers/pecos/LIBRARIES/PETSC3/petsc-3.1-p5/gcc-4.5-lucid-mpich2-1.2.1-cxx-opt/lib -lHYPRE -lsuperlu_dist_2.4 -lcmumps -ldmumps -lsmumps -lzmumps -lmumps_common -lpord -lparmetis -lmetis -lscalapack -lblacs -lsuperlu_4.0 -Wl,-rpath,/org/centers/pecos/LIBRARIES/MKL/mkl-10.0.3.020-gcc-4.5-lucid/lib/em64t -L/org/centers/pecos/LIBRARIES/MKL/mkl-10.0.3.020-gcc-4.5-lucid/lib/em64t -lmkl_solver_lp64_sequential -lmkl_intel_lp64 -lmkl_sequential -lmkl_core -lm -Wl,-rpath,/org/centers/pecos/LIBRARIES/MPICH2/mpich2-1.2.1-gcc-4.5-lucid/lib -L/org/centers/pecos/LIBRARIES/MPICH2/mpich2-1.2.1-gcc-4.5-lucid/lib -Wl,-rpath,/org/centers/pecos/LIBRARIES/GCC/gcc-4.5.1-lucid/lib/gcc/x86_64-unknown-linux-gnu/4.5.1 -L/org/centers/pecos/LIBRARIES/GCC/gcc-4.5.1-lucid/lib/gcc/x86_64-unknown-linux-gnu/4.5.1 -Wl,-rpath,/org/centers/pecos/LIBRARIES/GCC/gcc-4.5.1-lucid/lib64 -L/org/centers/pecos/LIBRARIES/GCC/gcc-4.5.1-lucid/lib64 -Wl,-rpath,/org/centers/pecos/LIBRARIES/GCC/gcc-4.5.1-lucid/lib -L/org/centers/pecos/LIBRARIES/GCC/gcc-4.5.1-lucid/lib -ldl -lmpich -lopa -lpthread -lrt -lgcc_s -lmpichf90 -lgfortran -lm -lm -lmpichcxx -lstdc++ -ldl -lmpich -lopa -lpthread -lrt -lgcc_s -ldl  
------------------------------------------
 
***************************************************************
* Done Running Example  ./ex18-opt
***************************************************************

Site Created By: libMesh Developers
Last modified: December 02 2011 21:19:02 UTC

Hosted By:
SourceForge.net Logo