visco_plastic.h

Go to the documentation of this file.

// ------------------------------------------------------
//
// SPDX-License-Identifier: GPL-2.0-or-later
//
// SPDX-FileCopyrightText: Copyright (C) 2020-2026 by the ASPECT authors.
//
// This file is part of ASPECT.
//
// Detailed license information governing the source code and
// contributions can be found in the folder LICENSES and in
// LICENSE and CONTRIBUTING.md at the top level directory.
//
// ------------------------------------------------------

#ifndef _aspect_material_model_rheology_visco_plastic_h
#define _aspect_material_model_rheology_visco_plastic_h

#include <aspect/global.h>
#include <aspect/material_model/interface.h>
#include <aspect/material_model/utilities.h>
#include <aspect/material_model/rheology/strain_dependent.h>
#include <aspect/material_model/rheology/friction_models.h>
#include <aspect/material_model/rheology/diffusion_creep.h>
#include <aspect/material_model/rheology/dislocation_creep.h>
#include <aspect/material_model/rheology/frank_kamenetskii.h>
#include <aspect/material_model/rheology/peierls_creep.h>
#include <aspect/material_model/rheology/constant_viscosity_prefactors.h>
#include <aspect/material_model/rheology/compositional_viscosity_prefactors.h>
#include <aspect/material_model/rheology/drucker_prager.h>
#include <aspect/material_model/rheology/elasticity.h>
#include <aspect/material_model/rheology/grain_boundary_sliding.h>
#include <aspect/simulator_access.h>

#include<deal.II/fe/component_mask.h>

namespace aspect
{
  namespace MaterialModel
  {
    template <int dim>
    class PlasticAdditionalOutputs : public NamedAdditionalMaterialOutputs<dim>
    {
      public:
        PlasticAdditionalOutputs(const unsigned int n_points);

        std::vector<double> get_nth_output(const unsigned int idx) const override;

        std::vector<double> cohesions;

        std::vector<double> friction_angles;

        std::vector<double> yield_stresses;

        std::vector<double> yielding;

    };

    template <int dim>
    class ViscosityAdditionalOutputs : public NamedAdditionalMaterialOutputs<dim>
    {
      public:
        enum class Property
        {
          diffusion_viscosity,
          dislocation_viscosity
        };

        ViscosityAdditionalOutputs(const unsigned int n_points,
                                   const std::vector<Property> &active_properties);

        std::vector<double>
        get_nth_output(const unsigned int idx) const override;

      private:
        std::vector<Property> active_properties;

      public:
        std::vector<double> diffusion_viscosities;

        std::vector<double> dislocation_viscosities;

    };

    struct IsostrainViscosities
    {
      std::vector<double> composition_viscosities;

      std::vector<bool> composition_yielding;

      std::vector<Rheology::DruckerPragerParameters> drucker_prager_parameters;

      std::vector<double> dilation_lhs_terms;

      std::vector<double> dilation_rhs_terms;

      std::vector<double> diffusion_viscosities;

      std::vector<double> dislocation_viscosities;
    };

    namespace Rheology
    {

      template <int dim>
      class ViscoPlastic : public ::aspect::SimulatorAccess<dim>
      {
        public:
          ViscoPlastic();

          IsostrainViscosities
          calculate_isostrain_viscosities ( const MaterialModel::MaterialModelInputs<dim> &in,
                                            const unsigned int i,
                                            const std::vector<double> &volume_fractions,
                                            const std::vector<double> &phase_function_values = {},
                                            const std::vector<unsigned int> &n_phase_transitions_per_composition =
                                              std::vector<unsigned int>()) const;

          void compute_viscosity_derivatives(const unsigned int point_index,
                                             const std::vector<double> &volume_fractions,
                                             const IsostrainViscosities &isostrain_values,
                                             const MaterialModel::MaterialModelInputs<dim> &in,
                                             MaterialModel::MaterialModelOutputs<dim> &out,
                                             const std::vector<double> &phase_function_values = {},
                                             const std::vector<unsigned int> &n_phase_transitions_per_composition =
                                               std::vector<unsigned int>()) const;

          ComponentMask get_volumetric_composition_mask() const;

          static
          void
          declare_parameters (ParameterHandler &prm);

          void
          parse_parameters (ParameterHandler &prm,
                            const std::unique_ptr<std::vector<unsigned int>> &expected_n_phases_per_composition = nullptr);

          void
          create_plastic_outputs (MaterialModel::MaterialModelOutputs<dim> &out) const;

          void fill_plastic_outputs(const unsigned int point_index,
                                    const std::vector<double> &volume_fractions,
                                    const bool plastic_yielding,
                                    const MaterialModel::MaterialModelInputs<dim> &in,
                                    MaterialModel::MaterialModelOutputs<dim> &out,
                                    const IsostrainViscosities &isostrain_viscosities) const;
          void
          create_viscosity_outputs(MaterialModel::MaterialModelOutputs<dim> &out) const;

          void fill_viscosity_outputs(const unsigned int point_index,
                                      const std::vector<double> &volume_fractions,
                                      MaterialModel::MaterialModelOutputs<dim> &out,
                                      const IsostrainViscosities &isostrain_viscosities) const;

          double min_strain_rate;

          MaterialUtilities::CompositionalAveragingOperation viscosity_averaging;

          Rheology::StrainDependent<dim> strain_rheology;

          Rheology::FrictionModels<dim> friction_models;

          Rheology::Elasticity<dim> elastic_rheology;


        private:

          double ref_strain_rate;

          std::vector<double> minimum_viscosity;
          std::vector<double> maximum_viscosity;

          enum ViscosityScheme
          {
            diffusion,
            dislocation,
            frank_kamenetskii,
            composite,
            minimum_diffusion_dislocation
          };
          ViscosityScheme viscous_flow_law;

          enum YieldScheme
          {
            stress_limiter,
            drucker_prager
          };
          YieldScheme yield_mechanism;

          bool allow_negative_pressures_in_plasticity;

          bool use_adiabatic_pressure_in_creep;

          bool use_adiabatic_pressure_in_plasticity;

          std::vector<double> exponents_stress_limiter;

          double adiabatic_temperature_gradient_for_viscosity;

          Rheology::DiffusionCreep<dim> diffusion_creep;
          Rheology::DislocationCreep<dim> dislocation_creep;
          std::unique_ptr<Rheology::FrankKamenetskii<dim>> frank_kamenetskii_rheology;

          bool use_grain_boundary_sliding;

          std::unique_ptr<Rheology::GrainBoundarySliding<dim>> grain_boundary_sliding_rheology;

          bool use_peierls_creep;

          std::unique_ptr<Rheology::PeierlsCreep<dim>> peierls_creep;

          Rheology::ConstantViscosityPrefactors<dim> constant_viscosity_prefactors;

          Rheology::CompositionalViscosityPrefactors<dim> compositional_viscosity_prefactors;

          /*
           * Object for computing plastic stresses, viscosities, and additional outputs
           */
          Rheology::DruckerPrager<dim> drucker_prager_plasticity;

      };
    }
  }
}
#endif