LKM-code-base / RotatingMHD

GNU General Public License v3.0
1 stars 0 forks source link

[master] Preconditioner parameters are not being read from the prm file #76

Closed j507 closed 3 years ago

j507 commented 3 years ago

I merged amg_preconditioner to coriolis_acceleration which leads to

Triangulation:
 Number of initial active cells           = 11200

Spatial discretization:
 Number of velocity degrees of freedom    = 98560
 Number of pressure degrees of freedom    = 12768
 Number of temperature degrees of freedom = 49280
 Number of total degrees of freedom       = 160608

  Heat Equation: Setting up matrices... done!
  Heat Equation: Setting up vectors... done!
  Heat Equation: Assembling constant matrices... done!

  Navier Stokes: Setting up matrices... done!
  Navier Stokes: Setting up vectors... done!
  Navier Stokes: Assembling velocity mass and stiffness matrices... done!
  Navier Stokes: Assembling pressure mass and stiffness matrices... done!

  Navier Stokes: Assembling Poisson pre-step's right hand side... done!
    Right-hand side's L2-norm = 4.493372e+00
  Navier Stokes: Solving the Poisson pre-step... done!
    Number of CG iterations: 54, Final residual: 2.812297e-06.

  Heat Equation: Assembling advection matrix... done!
  Heat Equation: Assembling right hand side... done!
    Right-hand side's L2-norm = 1.283326e+02
  Heat Equation: Solving linear system... done!
    Number of GMRES iterations: 2, Final residual: 1.072787e-04.

  Navier Stokes: Assembling advection matrix... done!
  Navier Stokes: Assembling diffusion step's right hand side... done!
    Right-hand side's L2-norm = 6.734343e+02
  Navier Stokes: Solving the diffusion step... done!
    Number of GMRES iterations: 5, Final residual: 5.530679e-04.
  Navier Stokes: Assembling the projection step's right hand side... done!
    Right-hand side's L2-norm = 3.243004e+00

Thread 1 "Christensen" received signal SIGSEGV, Segmentation fault.
0x00007fffee49275c in Teuchos::RCPNodeHandle::strong_count (this=0x68)
    at /usr/local/trilinos-12.18/include/Teuchos_RCPNode.hpp:926
926     if (node_) {

so I tried setting the AdditionalData of the PreconditionILU to its default values and noticed it was not being updated at all. Hopefully it is just that the values set in the constructor are just bad for this specific problem.

j507 commented 3 years ago

The preconditioner produced by the build_preconditioner method in Christensen leads to the segmentation fault. I replaced the preconditioner

          {
            build_preconditioner(correction_step_preconditioner,
                                 projection_mass_matrix,
                                 solver_parameters.preconditioner_parameters_ptr,
                                 (pressure->fe_degree > 1? true: false));
          }

          #ifdef USE_PETSC_LA
            LinearAlgebra::SolverCG solver(solver_control,
                                           mpi_communicator);
          #else
            LinearAlgebra::SolverCG solver(solver_control);
          #endif

          LinearAlgebra::MPI::PreconditionILU     preconditioner;
          preconditioner.initialize(projection_mass_matrix);

          try
          {
            solver.solve(projection_mass_matrix,
                         distributed_pressure,
                         correction_step_rhs,
                         preconditioner);
          }

in the projection step and correction step and it works. Curiously this does not happen in the other applications.

sebglane commented 3 years ago

Did you try to debug this using the GNU Debugger?

j507 commented 3 years ago

The terminal output from the first post is from the GNU Debugger (gdb ./Christensen)

j507 commented 3 years ago

To troubleshoot the error I created the coriolis_plus_master branch. It is basically that origin/coriolis_acceleration (22de0f0) with origin/master (710c0b7) merged into it.

Summarizing:

I added a PreconditionerILU to solve_projection_step and correction_step

  const typename RunTimeParameters::LinearSolverParameters &solver_parameters
    = parameters.projection_step_solver_parameters;
  if (reinit_prec)
  {
    build_preconditioner(projection_step_preconditioner,
                         phi_laplace_matrix,
                         solver_parameters.preconditioner_parameters_ptr,
                         (phi->fe_degree > 1? true: false));
  }

  LinearAlgebra::MPI::PreconditionILU     preconditioner;
  preconditioner.initialize(phi_laplace_matrix);

  SolverControl solver_control(
    solver_parameters.n_maximum_iterations,
    std::max(solver_parameters.relative_tolerance * projection_step_rhs.l2_norm(),
             solver_parameters.absolute_tolerance));

  #ifdef USE_PETSC_LA
    LinearAlgebra::SolverCG solver(solver_control,
                                   mpi_communicator);
  #else
    LinearAlgebra::SolverCG solver(solver_control);
  #endif

  try
  {
    solver.solve(phi_laplace_matrix,
                 distributed_phi,
                 projection_step_rhs,
                 /*preconditioner*/*projection_step_preconditioner);
  }

if ones uses preconditioner on both solve calls, the problem runs. At least until the 600th time step, where the AztecOO error occurs.

j507 commented 3 years ago

Update: I merged fix_pressure_constraint into my local coriolis_plus_master branch and the segmentation error still occurs. Here is a longer output from the GNU Debugger:

Thread 1 "Christensen" received signal SIGSEGV, Segmentation fault.
0x00007fffee49275c in Teuchos::RCPNodeHandle::strong_count (this=0x68)
    at /usr/local/trilinos-12.18/include/Teuchos_RCPNode.hpp:926
926     if (node_) {
(gdb) up
#1  0x00007fffee80bfca in Teuchos::RCP<Epetra_Operator>::strong_count (this=0x60)
    at /usr/local/trilinos-12.18/include/Teuchos_RCP.hpp:456
456   return node_.strong_count();
(gdb) up
#2  0x00007fffee80428d in dealii::TrilinosWrappers::SolverBase::set_preconditioner<dealii::TrilinosWrappers::PreconditionBase> (this=0x7fffffffb570, solver=..., preconditioner=...)
    at /home/j507/Source/dealii-9.2.0/source/lac/trilinos_solver.cc:570
570 /home/j507/Source/dealii-9.2.0/source/lac/trilinos_solver.cc: No such file or directory.
(gdb) up
#3  0x00007fffee809e67 in dealii::TrilinosWrappers::SolverBase::do_solve<dealii::TrilinosWrappers::PreconditionBase> (this=0x7fffffffb570, preconditioner=...)
    at /home/j507/Source/dealii-9.2.0/source/lac/trilinos_solver.cc:446
446 in /home/j507/Source/dealii-9.2.0/source/lac/trilinos_solver.cc
(gdb) up
#4  0x00007fffee80274b in dealii::TrilinosWrappers::SolverBase::solve (this=0x7fffffffb570, A=..., 
    x=..., b=..., preconditioner=...) at /home/j507/Source/dealii-9.2.0/source/lac/trilinos_solver.cc:88
88  in /home/j507/Source/dealii-9.2.0/source/lac/trilinos_solver.cc
(gdb) up
#5  0x0000555555c102ea in RMHD::NavierStokesProjection<2>::solve_projection_step (this=0x7fffffffc5e0, 
    reinit_prec=false)
    at /home/j507/Documents/Projekte/Backup/RotatingMHD/source/navier_stokes_projection/projection_step_methods.cc:60
60      solver.solve(phi_laplace_matrix,
(gdb) up
#6  0x0000555555ba88bd in RMHD::NavierStokesProjection<2>::projection_step (this=0x7fffffffc5e0, 
    reinit_prec=false)
    at /home/j507/Documents/Projekte/Backup/RotatingMHD/source/navier_stokes_projection/solve.cc:62
62    solve_projection_step(reinit_prec);

Could you check if you can also replicate the error?

sebglane commented 3 years ago

Here is what I did:

git checkout coriolis_plus_master
git pull
git merge origin/fix_pressure_constrain
git mergetool

I had to merge some files manually using meld. I hope I didn't break anything here.

Next:

make
bash Christensen.sh 2

This is the output:

[ 58%] Built target rotatingMHD
[ 66%] Built target Christensen
[ 66%] Built target Couette
[ 69%] Built target ThermalTGV
[ 73%] Built target GuermondNeumannBC
[ 77%] Built target MIT
[ 81%] Built target Guermond
[ 84%] Built target TGV
[ 88%] Built target step-35
[ 92%] Built target timestepping_coefficients
[ 96%] Built target DFG
[100%] Built target test_preconditioners_dealii
rm: cannot remove '*.pvtu': No such file or directory
rm: cannot remove '*.vtu': No such file or directory
+-----------------------------------------+--------------------+
| Problem parameters                                        |
+-----------------------------------------+--------------------+
| Problem type                             | rotating_boussinesq  |
| Spatial dimension                        | 3                    |
| Mapping                                  | MappingQ<3>(2)       |
| Mapping - Apply to interior cells        | true                 |
| Finite Element - Velocity                | FE_Q<3>(2)^3         |
| Finite Element - Pressure                | FE_Q<3>(1)           |
| Finite Element - Temperature             | FE_Q<3>(2)           |
| Verbose                                  | false                |
+-----------------------------------------+--------------------+
| Output control parameters                                 |
+-----------------------------------------+--------------------+
| Graphical output frequency               | 1                    |
| Terminal output frequency                | 1                    |
| Graphical output directory               | ChristensenResults/  |
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
| Dimensionless numbers                                     |
+-----------------------------------------+--------------------+
| Prandtl number                           | 1                    |
| Rayleigh number                          | 100000               |
| Ekman number                             | 0.001                |
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
| Refinement control parameters                             |
+-----------------------------------------+--------------------+
| Adaptive mesh refinement                 | False                |
| Number of initial adapt. refinements     | 0                    |
| Number of initial global refinements     | 5                    |
| Number of initial boundary refinements   | 2                    |
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
| Timestepping parameters                                      |
+-----------------------------------------+--------------------+
| IMEX scheme                             | BDF2               |
| Maximum number of time steps            | 10                 |
| Adaptive timestepping                   | false              |
| Initial time step                       | 0.0001             |
| Start time                              | 0                  |
| Final time                              | 1.2                |
| Verbose                                 | true               |
+-----------------------------------------+--------------------+
| Navier-Stokes discretization parameters                   |
+-----------------------------------------+--------------------+
| Incremental pressure-correction scheme   | rotational           |
| Convective term weak form                | skew-symmetric       |
| Convective temporal form                 | semi-implicit        |
| Preconditioner update frequency          | 10                   |
+-----------------------------------------+--------------------+
| Linear solver parameters                                  |
+-----------------------------------------+--------------------+
| Maximum number of iterations             | 1000                 |
| Relative tolerance                       | 1e-06                |
| Absolute tolerance                       | 1e-09                |
+-----------------------------------------+--------------------+
| Precondition ILU Parameters                               |
+-----------------------------------------+--------------------+
| Fill-in level                            | 1                    |
| Overlap                                  | 1                    |
| Relative tolerance                       | 1                    |
| Absolute tolerance                       | 0                    |
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
| Linear solver parameters                                  |
+-----------------------------------------+--------------------+
| Maximum number of iterations             | 1000                 |
| Relative tolerance                       | 1e-06                |
| Absolute tolerance                       | 1e-09                |
+-----------------------------------------+--------------------+
| Precondition ILU Parameters                               |
+-----------------------------------------+--------------------+
| Fill-in level                            | 1                    |
| Overlap                                  | 1                    |
| Relative tolerance                       | 1                    |
| Absolute tolerance                       | 0                    |
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
| Linear solver parameters                                  |
+-----------------------------------------+--------------------+
| Maximum number of iterations             | 1000                 |
| Relative tolerance                       | 1e-06                |
| Absolute tolerance                       | 1e-09                |
+-----------------------------------------+--------------------+
| Precondition ILU Parameters                               |
+-----------------------------------------+--------------------+
| Fill-in level                            | 1                    |
| Overlap                                  | 1                    |
| Relative tolerance                       | 1                    |
| Absolute tolerance                       | 0                    |
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
| Linear solver parameters                                  |
+-----------------------------------------+--------------------+
| Maximum number of iterations             | 1000                 |
| Relative tolerance                       | 1e-06                |
| Absolute tolerance                       | 1e-09                |
+-----------------------------------------+--------------------+
| Precondition ILU Parameters                               |
+-----------------------------------------+--------------------+
| Fill-in level                            | 1                    |
| Overlap                                  | 1                    |
| Relative tolerance                       | 1                    |
| Absolute tolerance                       | 0                    |
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
| Heat equation solver parameters                           |
+-----------------------------------------+--------------------+
| Convective term weak form                | skew-symmetric       |
| Convective temporal form                 | semi-implicit        |
| Preconditioner update frequency          | 10                   |
+-----------------------------------------+--------------------+
| Linear solver parameters                                  |
+-----------------------------------------+--------------------+
| Maximum number of iterations             | 1000                 |
| Relative tolerance                       | 1e-06                |
| Absolute tolerance                       | 1e-09                |
+-----------------------------------------+--------------------+
| Precondition ILU Parameters                               |
+-----------------------------------------+--------------------+
| Fill-in level                            | 1                    |
| Overlap                                  | 1                    |
| Relative tolerance                       | 1                    |
| Absolute tolerance                       | 0                    |
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+
+-----------------------------------------+--------------------+

+----------+----------+----------+----------+----------+----------+
|    C1    |    C2    |    C3    |    C4    |    C5    |    C6    |
+----------+----------+----------+----------+----------+----------+
|  2.0e+03 |  1.0e+00 |  1.0e+05 |  1.0e+00 |  0.0e+00 |  1.0e+03 |
+----------+----------+----------+----------+----------+----------+

Triangulation:
 Number of initial active cells           = 11200

Spatial discretization:
 Number of velocity degrees of freedom    = 98560
 Number of pressure degrees of freedom    = 12768
 Number of temperature degrees of freedom = 49280
 Number of total degrees of freedom       = 160608

+--------------------+-----------------------------------------+
| Boundary conditions of the Velocity entity                   |
+--------------------+-----------------------------------------+
| Dirichlet boundary conditions                                |
| Boundary id       |    Type name                             |
| 0                 | ZeroFunction<2, double>                  |
| 1                 | ZeroFunction<2, double>                  |
+--------------------+-----------------------------------------+
---------------------------------------------------------
TimerOutput objects finalize timed values printed to the
screen by communicating over MPI in their destructors.
Since an exception is currently uncaught, this
synchronization (and subsequent output) will be skipped
to avoid a possible deadlock.
---------------------------------------------------------

----------------------------------------------------
Exception on processing: 

----------------------------------------------------
Exception on processing: 

--------------------------------------------------------
An error occurred in line <507> of file </home/sg/Documents/Projects/RotatingMHD/source/entities_structs.cc> in function
    void RMHD::Entities::ScalarEntity<dim>::apply_boundary_conditions() [with int dim = 2]
The violated condition was: 
    boundary_conditions.is_regularity_guaranteed()
Additional information: 
    No boundary conditions were set for the "Pressure" entity
--------------------------------------------------------

Aborting!
----------------------------------------------------
--------------------------------------------------------
An error occurred in line <507> of file </home/sg/Documents/Projects/RotatingMHD/source/entities_structs.cc> in function
    void RMHD::Entities::ScalarEntity<dim>::apply_boundary_conditions() [with int dim = 2]
The violated condition was: 
    boundary_conditions.is_regularity_guaranteed()
Additional information: 
    No boundary conditions were set for the "Pressure" entity
--------------------------------------------------------
sebglane commented 3 years ago

I guess I found the bug regarding the parsing of the preconditioner parameters. Please see my comment in #79.

j507 commented 3 years ago

The problem was caused ny my manual call of navier_stokes.setup in order to output the pressure field computed by the Poisson pre-step. So when the first call of navier_stokes.solve() occurs

template <int dim>
void NavierStokesProjection<dim>::solve()
{
  if (velocity->solution.size() != diffusion_step_rhs.size())
  {
    setup();

    diffusion_step(true);

    projection_step(true);

    pressure_correction(true);

    flag_matrices_were_updated = false;
  }
  else
  {
    diffusion_step(time_stepping.get_step_number() %
                   parameters.preconditioner_update_frequency == 0);

    projection_step(false);

    pressure_correction(false);
  }

  phi->update_solution_vectors();
}

the else scope is being executed, i.e. the preconditioner was not being initialized. Sorry for the trouble.