Program Listing for File elasticity.h¶
↰ Return to documentation for file (include/elasticity.h)
// ---------------------------------------------------------------------
//
// Copyright (C) 2024 by Luca Heltai
//
// This file is part of the reduced_lagrange_multipliers application, based on
// the deal.II library.
//
// The reduced_lagrange_multipliers application is free software; you can use
// it, redistribute it, and/or modify it under the terms of the Apache-2.0
// License WITH LLVM-exception as published by the Free Software Foundation;
// either version 3.0 of the License, or (at your option) any later version. The
// full text of the license can be found in the file LICENSE.md at the top level
// of the reduced_lagrange_multipliers distribution.
//
// ---------------------------------------------------------------------
/* ---------------------------------------------------------------------
*/
#ifndef dealii_distributed_lagrange_multiplier_elasticity_h
#define dealii_distributed_lagrange_multiplier_elasticity_h
#include <deal.II/base/function.h>
#include <deal.II/base/parsed_convergence_table.h>
#include <deal.II/base/quadrature_lib.h>
#include <deal.II/base/timer.h>
#include <deal.II/lac/block_linear_operator.h>
#include <deal.II/lac/generic_linear_algebra.h>
#include <deal.II/lac/linear_operator.h>
#include <deal.II/lac/linear_operator_tools.h>
#define FORCE_USE_OF_TRILINOS
namespace LA
{
#if defined(DEAL_II_WITH_PETSC) && !defined(DEAL_II_PETSC_WITH_COMPLEX) && \
!(defined(DEAL_II_WITH_TRILINOS) && defined(FORCE_USE_OF_TRILINOS))
using namespace dealii::LinearAlgebraPETSc;
# define USE_PETSC_LA
#elif defined(DEAL_II_WITH_TRILINOS)
using namespace dealii::LinearAlgebraTrilinos;
#else
# error DEAL_II_WITH_PETSC or DEAL_II_WITH_TRILINOS required
#endif
} // namespace LA
#include <deal.II/base/conditional_ostream.h>
#include <deal.II/base/index_set.h>
#include <deal.II/base/parameter_acceptor.h>
#include <deal.II/base/parsed_function.h>
#include <deal.II/base/utilities.h>
#include <deal.II/base/work_stream.h>
#include <deal.II/distributed/grid_refinement.h>
#include <deal.II/distributed/solution_transfer.h>
#include <deal.II/distributed/tria.h>
#include <deal.II/dofs/dof_handler.h>
#include <deal.II/dofs/dof_renumbering.h>
#include <deal.II/dofs/dof_tools.h>
#include <deal.II/fe/fe_nothing.h>
#include <deal.II/fe/fe_q.h>
#include <deal.II/fe/fe_system.h>
#include <deal.II/fe/fe_values.h>
#include <deal.II/fe/mapping_fe_field.h>
#include <deal.II/fe/mapping_q.h>
#include <deal.II/grid/grid_generator.h>
#include <deal.II/grid/grid_in.h>
#include <deal.II/grid/grid_out.h>
#include <deal.II/grid/grid_refinement.h>
#include <deal.II/grid/grid_tools.h>
#include <deal.II/grid/grid_tools_cache.h>
#include <deal.II/grid/manifold_lib.h>
#include <deal.II/grid/tria.h>
#include <deal.II/grid/tria_accessor.h>
#include <deal.II/grid/tria_iterator.h>
#include <deal.II/lac/affine_constraints.h>
#include <deal.II/lac/dynamic_sparsity_pattern.h>
#include <deal.II/lac/full_matrix.h>
#include <deal.II/lac/petsc_precondition.h>
#include <deal.II/lac/petsc_solver.h>
#include <deal.II/lac/petsc_sparse_matrix.h>
#include <deal.II/lac/petsc_vector.h>
#include <deal.II/lac/solver_cg.h>
#include <deal.II/lac/solver_gmres.h>
#include <deal.II/lac/solver_minres.h>
#include <deal.II/lac/sparsity_tools.h>
#include <deal.II/lac/trilinos_precondition.h>
#include <deal.II/lac/trilinos_solver.h>
#include <deal.II/lac/trilinos_sparse_matrix.h>
#include <deal.II/lac/trilinos_vector.h>
#include <deal.II/physics/elasticity/standard_tensors.h>
// #include <deal.II/trilinos/parameter_acceptor.h>
#include <deal.II/lac/vector.h>
#include <deal.II/meshworker/dof_info.h>
#include <deal.II/meshworker/integration_info.h>
#include <deal.II/meshworker/loop.h>
#include <deal.II/meshworker/scratch_data.h>
#include <deal.II/meshworker/simple.h>
#include <deal.II/numerics/data_out.h>
#include <deal.II/numerics/data_out_faces.h>
#include <deal.II/numerics/error_estimator.h>
#include <deal.II/numerics/vector_tools.h>
// #include <deal.II/numerics/matrix_tools.h>
#include <deal.II/opencascade/manifold_lib.h>
#include <deal.II/opencascade/utilities.h>
#include "inclusions.h"
#ifdef DEAL_II_WITH_OPENCASCADE
# include <TopoDS.hxx>
#endif
#include <deal.II/base/hdf5.h>
#include <cmath>
#include <filesystem>
#include <fstream>
#include <iomanip>
#include <iostream>
#include <map>
#include <memory>
#include <set>
#include <string>
#include <system_error>
#include "elasticity_problem_parameters.h"
#include "material_properties.h"
template <int dim>
class RigidBodyMotion : public Function<dim>
{
public:
RigidBodyMotion(const unsigned int type_);
virtual double
value(const Point<dim> &p, const unsigned int component) const override;
private:
const unsigned int type;
};
template <int dim, int spacedim = dim>
class ElasticityProblem : public Subscriptor
{
public:
ElasticityProblem(const ElasticityProblemParameters<dim, spacedim> &par);
void
make_grid();
void
setup_fe();
void
setup_dofs();
void
assemble_elasticity_system();
void
assemble_forcing_terms();
void
compute_system_rhs();
void
assemble_coupling();
void
run();
void
run_static();
void
run_quasistatic();
void
run_newmark();
IndexSet
assemble_coupling_sparsity(DynamicSparsityPattern &dsp);
void
solve();
void
solve_static();
void
solve_quasistatic();
void
solve_newmark();
void
refine_and_transfer();
void
refine_and_transfer_around_inclusions();
void
execute_actual_refine_and_transfer();
std::string
output_solution() const;
void
output_results() const;
void
print_parameters() const;
void
compute_internal_and_boundary_stress(bool) const;
void
output_lambda() const;
std::string
output_stresses() const;
// void
// compute_face_stress();
const ElasticityProblemParameters<dim, spacedim> ∥
MPI_Comm mpi_communicator;
ConditionalOStream pcout;
mutable TimerOutput computing_timer;
parallel::distributed::Triangulation<spacedim> tria;
std::unique_ptr<FiniteElement<spacedim>> fe;
Inclusions<spacedim> inclusions;
std::unique_ptr<Quadrature<spacedim>> quadrature;
std::unique_ptr<Quadrature<spacedim - 1>> face_quadrature_formula;
DoFHandler<spacedim> dh;
std::vector<IndexSet> owned_dofs;
std::vector<IndexSet> relevant_dofs;
AffineConstraints<double> constraints;
AffineConstraints<double>
inclusion_constraints;
AffineConstraints<double> mean_value_constraints;
LA::MPI::SparseMatrix stiffness_matrix;
LA::MPI::SparseMatrix
newmark_matrix;
LA::MPI::SparseMatrix mass_matrix;
LA::MPI::SparseMatrix coupling_matrix;
LA::MPI::SparseMatrix damping_matrix;
LA::MPI::PreconditionAMG prec_A;
LA::MPI::PreconditionAMG prec_newmark;
LA::MPI::PreconditionAMG prec_C;
LA::MPI::PreconditionAMG prec_M;
LA::MPI::SparseMatrix inclusion_matrix;
LA::MPI::BlockVector solution;
LA::MPI::BlockVector velocity;
LA::MPI::BlockVector acceleration;
LA::MPI::BlockVector predictor;
LA::MPI::BlockVector corrector;
LA::MPI::BlockVector locally_relevant_solution;
LA::MPI::BlockVector force_rhs;
LA::MPI::BlockVector bc_rhs;
LA::MPI::BlockVector neumann_bc_rhs;
LA::MPI::BlockVector system_rhs;
std::vector<std::vector<BoundingBox<spacedim>>> global_bounding_boxes;
unsigned int cycle = 0;
unsigned int time_step = 0;
FEValuesExtractors::Vector displacement;
std::map<types::boundary_id, Tensor<1, spacedim>> forces;
std::map<types::boundary_id, Tensor<1, spacedim>>
average_displacements;
std::map<types::boundary_id, Tensor<1, spacedim>>
average_normals;
std::map<types::boundary_id, double>
areas;
TrilinosWrappers::MPI::Vector sigma_n;
double current_time = 0.0;
class Postprocessor;
};
#endif