Skip to content

Commit

Permalink
(A) WIP of bouncing algorithm and clean-up of monolithic algorithm
Browse files Browse the repository at this point in the history
  • Loading branch information
@j507
j507 committed Jul 11, 2024
1 parent 3f2fb3a commit 121ddde
Show file tree
Hide file tree
Showing 4 changed files with 408 additions and 291 deletions.
5 changes: 4 additions & 1 deletion include/gCP/fe_field.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -309,7 +309,10 @@ dealii::types::global_dof_index get_n_plastic_slip_dofs() const;
*/
unsigned int get_n_components() const;

std::vector<double> get_l2_norms(
double get_l2_norm(
const dealii::LinearAlgebraTrilinos::MPI::BlockVector &vector) const;

dealii::Vector<double> get_sub_l2_norms(
const dealii::LinearAlgebraTrilinos::MPI::BlockVector &vector) const;

/*!
Expand Down
53 changes: 18 additions & 35 deletions include/gCP/gradient_crystal_plasticity.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -57,7 +57,7 @@ class GradientCrystalPlasticitySolver
void set_macroscopic_strain(
const dealii::SymmetricTensor<2,dim> macroscopic_strain);

std::tuple<bool,unsigned int> solve_nonlinear_system(
void solve_nonlinear_system(
const bool flag_skip_extrapolation = false);

std::shared_ptr<const Kinematics::ElasticStrain<dim>>
Expand Down Expand Up @@ -89,28 +89,6 @@ class GradientCrystalPlasticitySolver

private:

struct SolverData
{
SolverData()
:
flag_successful_convergence(false),
flag_compute_active_set(true),
nonlinear_iteration(0),
residual_l2_norms(3,0.),
old_residual_l2_norms(3,0.)
{};

bool flag_successful_convergence;

bool flag_compute_active_set;

unsigned int nonlinear_iteration;

std::vector<double> residual_l2_norms;

std::vector<double> old_residual_l2_norms;
};

const RunTimeParameters::SolverParameters &parameters;

const dealii::DiscreteTime &discrete_time;
Expand Down Expand Up @@ -248,18 +226,20 @@ class GradientCrystalPlasticitySolver
void copy_local_to_global_quadrature_point_history(
const gCP::AssemblyData::QuadraturePointHistory::Copy &){};

void monolithic_algorithm(SolverData &solver_data);
void monolithic_algorithm();

void bouncing_algorithm(SolverData &solver_data);
void bouncing_algorithm();

void embracing_algorihtm(SolverData &solver_data);
void embracing_algorihtm();

unsigned int solve_linearized_system(
const RunTimeParameters::KrylovParameters &krylov_parameters,
const unsigned int block_id,
const double right_hand_side_l2_norm);

void update_trial_solution(const double relaxation_parameter);
void update_trial_solution(
const double relaxation_parameter,
const unsigned int block_id = 0);

void update_trial_solution(const std::vector<double>
relaxation_parameter);
Expand All @@ -269,6 +249,7 @@ class GradientCrystalPlasticitySolver
double line_search_algorithm(
dealii::Vector<double> residual_l2_norms,
const unsigned int block_id);

void store_trial_solution(
const bool flag_store_initial_trial_solution = false);

Expand All @@ -278,17 +259,19 @@ class GradientCrystalPlasticitySolver
void extrapolate_initial_trial_solution(
const bool flag_skip_extrapolation = false);

void distribute_affine_constraints_to_trial_solution();

void update_and_output_nonlinear_solver_logger(
const std::vector<double> residual_l2_norms);
const dealii::Vector<double> residual_l2_norms);

void update_and_output_nonlinear_solver_logger(
const unsigned int nonlinear_iteration,
const unsigned int n_krylov_iterations,
const unsigned int n_line_search_iterations,
const std::vector<double> newton_update_l2_norms,
const std::vector<double> residual_l2_norms,
const double order_of_convergence,
const double relaxation_parameter = 1.0);
const unsigned int nonlinear_iteration,
const unsigned int n_krylov_iterations,
const unsigned int n_line_search_iterations,
const dealii::Vector<double> newton_update_l2_norms,
const dealii::Vector<double> residual_l2_norms,
const double order_of_convergence,
const double relaxation_parameter = 1.0);

/*!
* @note Only for debugging purposes
Expand Down
78 changes: 46 additions & 32 deletions source/fe_field.cc
View file
Original file line number Diff line number Diff line change
Expand Up @@ -584,7 +584,26 @@ void FEField<dim>::reset_all_affine_constraints()


template <int dim>
std::vector<double> FEField<dim>::get_l2_norms(
double FEField<dim>::get_l2_norm(
const dealii::LinearAlgebraTrilinos::MPI::BlockVector &vector) const
{
if (vector.has_ghost_elements())
{
dealii::LinearAlgebraTrilinos::MPI::BlockVector tmp =
get_distributed_vector_instance(vector);

return (tmp.l2_norm());
}
else
{
return (vector.l2_norm());
}
}



template <int dim>
dealii::Vector<double> FEField<dim>::get_sub_l2_norms(
const dealii::LinearAlgebraTrilinos::MPI::BlockVector &vector) const
{
dealii::LinearAlgebraTrilinos::MPI::BlockVector distributed_vector_tmp;
Expand All @@ -593,9 +612,9 @@ std::vector<double> FEField<dim>::get_l2_norms(

distributed_vector_tmp = vector;

std::vector<double> l2_norms(3, 0.);
dealii::Vector<double> l2_norms(2);

l2_norms[0] = distributed_vector_tmp.l2_norm();
const double l2_norm = distributed_vector_tmp.l2_norm();

if (flag_use_single_block)
{
Expand Down Expand Up @@ -636,42 +655,37 @@ std::vector<double> FEField<dim>::get_l2_norms(
scalar_squared_entries,
MPI_COMM_WORLD);

const double control_l2_norm =
std::sqrt(vector_squared_entries + scalar_squared_entries);

auto to_string =
[](const double number)
{
std::ostringstream out;
out.precision(10);
out << std::scientific << number;
return std::move(out).str();
};

const double factor = 1e4;

Assert(
std::fabs(l2_norms[0] - control_l2_norm) <
std::numeric_limits<double>::epsilon() * factor,
dealii::ExcMessage("The norms do not match (" +
to_string(l2_norms[0]) + ", " + to_string(control_l2_norm) +
")"));
l2_norms[0] = std::sqrt(vector_squared_entries);

(void)control_l2_norm;

(void)to_string;

l2_norms[1] = std::sqrt(vector_squared_entries);

l2_norms[2] = std::sqrt(scalar_squared_entries);
l2_norms[1] = std::sqrt(scalar_squared_entries);
}
else
{
l2_norms[1] = distributed_vector_tmp.block(0).l2_norm();
l2_norms[0] = distributed_vector_tmp.block(0).l2_norm();

l2_norms[2] = distributed_vector_tmp.block(1).l2_norm();
l2_norms[1] = distributed_vector_tmp.block(1).l2_norm();
}

auto to_string =
[](const double number)
{
std::ostringstream out;
out.precision(10);
out << std::scientific << number;
return std::move(out).str();
};

const double factor = 1e4;

Assert(
std::fabs(l2_norm - l2_norms.l2_norm()) <
std::numeric_limits<double>::epsilon() * factor,
dealii::ExcMessage("The norms do not match (" +
to_string(l2_norm) + ", " + to_string(l2_norms.l2_norm()) +
")"));

(void)to_string;

return l2_norms;
}

Expand Down
Loading

0 comments on commit 121ddde

Please sign in to comment.