Skip to content

Commit

Permalink
Correct bug
Browse files Browse the repository at this point in the history
  • Loading branch information
hepap committed Apr 5, 2024
1 parent 955ac1e commit 1051d53
Show file tree
Hide file tree
Showing 4 changed files with 66 additions and 75 deletions.
121 changes: 57 additions & 64 deletions source/solvers/cahn_hilliard.cc
Original file line number Diff line number Diff line change
Expand Up @@ -1327,78 +1327,71 @@ void
CahnHilliard<dim>::output_newton_update_norms(
const unsigned int display_precision)
{
auto mpi_communicator = triangulation->get_communicator();

FEValuesExtractors::Scalar phase_order(0);
FEValuesExtractors::Scalar chemical_potential(1);

ComponentMask phase_order_mask = fe->component_mask(phase_order);
ComponentMask chemical_potential_mask = fe->component_mask(chemical_potential);

const std::vector<IndexSet> index_set_phase_order =
DoFTools::locally_owned_dofs_per_component(dof_handler, phase_order_mask);
const std::vector<IndexSet> index_set_chemical_potential =
DoFTools::locally_owned_dofs_per_component(dof_handler, chemical_potential_mask);

double local_sum = 0.0;
double local_max = std::numeric_limits<double>::lowest();


for (auto j = index_set_phase_order[0].begin();
j != index_set_phase_order[0].end();
j++)
{
double dof_newton_update = newton_update[*j];
auto mpi_communicator = triangulation->get_communicator();

local_sum += dof_newton_update * dof_newton_update;
FEValuesExtractors::Scalar phase_order(0);
FEValuesExtractors::Scalar chemical_potential(1);

local_max = std::max(local_max, dof_newton_update);
}
ComponentMask phase_order_mask = fe->component_mask(phase_order);
ComponentMask chemical_potential_mask =
fe->component_mask(chemical_potential);

double global_phase_order_l2_norm =
std::sqrt(Utilities::MPI::sum(local_sum, mpi_communicator));
double global_phase_order_linfty_norm =
Utilities::MPI::max(local_max, mpi_communicator);
const std::vector<IndexSet> index_set_phase_order =
DoFTools::locally_owned_dofs_per_component(dof_handler, phase_order_mask);
const std::vector<IndexSet> index_set_chemical_potential =
DoFTools::locally_owned_dofs_per_component(dof_handler,
chemical_potential_mask);

local_sum = 0.0;
local_max = std::numeric_limits<double>::lowest();
double local_sum = 0.0;
double local_max = std::numeric_limits<double>::lowest();

for (auto j = index_set_chemical_potential[1].begin();
j != index_set_chemical_potential[1].end();
j++)
{
double dof_newton_update = newton_update[*j];

local_sum += dof_newton_update * dof_newton_update;
for (auto j = index_set_phase_order[0].begin();
j != index_set_phase_order[0].end();
j++)
{
double dof_newton_update = newton_update[*j];

local_sum += dof_newton_update * dof_newton_update;

local_max = std::max(local_max, dof_newton_update);
}
local_max = std::max(local_max, std::abs(dof_newton_update));
}


double global_phase_order_l2_norm =
std::sqrt(Utilities::MPI::sum(local_sum, mpi_communicator));
double global_phase_order_linfty_norm =
Utilities::MPI::max(local_max, mpi_communicator);

local_sum = 0.0;
local_max = std::numeric_limits<double>::lowest();

for (auto j = index_set_chemical_potential[1].begin();
j != index_set_chemical_potential[1].end();
j++)
{
double dof_newton_update = newton_update[*j];

local_sum += dof_newton_update * dof_newton_update;

local_max = std::max(local_max, std::abs(dof_newton_update));
}

double global_chemical_potential_l2_norm =
std::sqrt(Utilities::MPI::sum(local_sum, mpi_communicator));
double global_chemical_potential_linfty_norm =
Utilities::MPI::max(local_max, mpi_communicator);

this->pcout << std::setprecision(display_precision)
<< "\n\t||du||_L2 = " << std::setw(6)
<< newton_update.l2_norm() << std::setw(6)
<< "\t||du||_Linfty = "
<< std::setprecision(display_precision)
<< newton_update.linfty_norm() << std::endl;

this->pcout << std::setprecision(display_precision)
<< "\n\t||dphi||_L2 = " << std::setw(6)
<< global_phase_order_l2_norm << std::setw(6)
<< "\t||dphi||_Linfty = "
<< std::setprecision(display_precision)
<< global_phase_order_linfty_norm << std::endl;
this->pcout << std::setprecision(display_precision)
<< "\t||deta||_L2 = " << std::setw(6) << global_chemical_potential_l2_norm
<< std::setw(6) << "\t||deta||_Linfty = "
<< std::setprecision(display_precision)
<< global_chemical_potential_linfty_norm << std::endl;

double global_chemical_potential_l2_norm =
std::sqrt(Utilities::MPI::sum(local_sum, mpi_communicator));
double global_chemical_potential_linfty_norm =
Utilities::MPI::max(local_max, mpi_communicator);

this->pcout << std::setprecision(display_precision)
<< "\n\t||dphi||_L2 = " << std::setw(6)
<< global_phase_order_l2_norm << std::setw(6)
<< "\t||dphi||_Linfty = " << std::setprecision(display_precision)
<< global_phase_order_linfty_norm << std::endl;
this->pcout << std::setprecision(display_precision)
<< "\t||deta||_L2 = " << std::setw(6)
<< global_chemical_potential_l2_norm << std::setw(6)
<< "\t||deta||_Linfty = " << std::setprecision(display_precision)
<< global_chemical_potential_linfty_norm << std::endl;
}

template std::pair<Tensor<1, 2>, Tensor<1, 2>>
Expand Down
4 changes: 2 additions & 2 deletions source/solvers/navier_stokes_base.cc
Original file line number Diff line number Diff line change
Expand Up @@ -2516,7 +2516,7 @@ NavierStokesBase<dim, VectorType, DofsType>::output_newton_update_norms(

local_sum += dof_newton_update * dof_newton_update;

local_max = std::max(local_max, dof_newton_update);
local_max = std::max(local_max, std::abs(dof_newton_update));
}
}

Expand All @@ -2536,7 +2536,7 @@ NavierStokesBase<dim, VectorType, DofsType>::output_newton_update_norms(

local_sum += dof_newton_update * dof_newton_update;

local_max = std::max(local_max, dof_newton_update);
local_max = std::max(local_max, std::abs(dof_newton_update));
}

double global_pressure_l2_norm =
Expand Down
8 changes: 4 additions & 4 deletions tests/core/non_linear_test_system_01.h
Original file line number Diff line number Diff line change
Expand Up @@ -130,10 +130,10 @@ class NonLinearProblemTestClass : public PhysicsSolver<Vector<double>>
output_newton_update_norms(const unsigned int display_precision) override
{
deallog << std::setprecision(display_precision)
<< "\t||dx||_L2 = " << std::setw(6) << newton_update.l2_norm()
<< std::setw(6)
<< "\t||dx||_Linfty = " << std::setprecision(display_precision)
<< newton_update.linfty_norm() << std::endl;
<< "\t||dx||_L2 = " << std::setw(6) << newton_update.l2_norm()
<< std::setw(6)
<< "\t||dx||_Linfty = " << std::setprecision(display_precision)
<< newton_update.linfty_norm() << std::endl;
};
virtual Vector<double> &
get_present_solution() override
Expand Down
8 changes: 3 additions & 5 deletions tests/core/vector_problem.cc
Original file line number Diff line number Diff line change
Expand Up @@ -105,7 +105,8 @@ test()
{
correction_norm += dummy_solution[*j] * dummy_solution[*j];

max_correction = std::max(max_correction, dummy_solution[*j]);
max_correction =
std::max(max_correction, std::abs(dummy_solution[*j]));
}
}

Expand All @@ -120,10 +121,7 @@ test()
{
correction_norm += dummy_solution[*j] * dummy_solution[*j];

if (dummy_solution[*j] > max_correction)
{
max_correction = dummy_solution[*j];
}
max_correction = std::max(max_correction, std::abs(dummy_solution[*j]));
}

deallog << "||p||_L2 : " << std::sqrt(correction_norm) << std::endl;
Expand Down

0 comments on commit 1051d53

Please sign in to comment.