Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Add NamedAdditionalOutput for strain weakened plastic parameters #1763

Merged

Conversation

anne-glerum
Copy link
Contributor

A visualization output that shows the cohesion and friction angle as weakened by the current value of the tracked strain.

@naliboff any comments? I hope this way it doesn't interfere with your plans to separate the viscosity calculation out into several functions too much.

Copy link
Member

@gassmoeller gassmoeller left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Only a few comments. Have you checked that the tests give reasonable output? Also please add an entry to the changelog so that people are aware of the new output.

{
std::vector<std::string> names;
names.push_back("weakened_cohesions");
names.push_back("weakened_friction_angles");
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I am not quite happy with the term weakened, what if people implement a healing algorithm that might make the values higher than the initial values? Would you be ok with current_...?

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I agree with @gassmoeller here. While the original intent was only to do strain weakening, it would be good to leave the door open for implementing strain hardening and not having to rename these terms. Another term that could be used is modified....

:
NamedAdditionalMaterialOutputs<dim>(make_plastic_additional_outputs_names()),
cohesions(n_points, -1.0),
friction_angles(n_points, -1.0)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please replace the -1 by numbers::signaling_nan<double>() in both lines. -1 was a leftover from early days in the seismic velocities, and we would like to get rid of it.

@@ -747,6 +747,46 @@ namespace aspect
return vs;
}

namespace
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Which material models will use this property? I am asking, because I would prefer this to be moved in some place other than the interface (it is only relevant for plastic material models). If there is only one (or several, but they are all derived from one base model), then please move all this code into the header and source file of that class (see e.g. the handling of MeltOutputs in melt.h and melt.cc).

const double strain_fraction = ( cut_off_strain_ii - start_strain_weakening_intervals[j] ) /
( start_strain_weakening_intervals[j] - end_strain_weakening_intervals[j] );
const double weakened_coh = cohesions[j] + ( cohesions[j] - cohesions[j] * cohesion_strain_weakening_factors[j] ) * strain_fraction;
const double weakened_phi = angles_internal_friction[j] + ( angles_internal_friction[j] - angles_internal_friction[j] * friction_strain_weakening_factors[j] ) * strain_fraction;
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Same here, I would prefer current over weakened

@@ -353,6 +365,7 @@ namespace aspect
// of compositional field viscosities is consistent with any averaging scheme.
out.viscosities[i] = average_value(composition, composition_viscosities, viscosity_averaging);


Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

remove this empty line

@gassmoeller
Copy link
Member

/run-tests

@naliboff
Copy link
Contributor

@anne-glerum - nice, thank you for adding this in! This does not interfere with any plans I had for visco_plastic. Honestly, I think it easier to just make a new material model that is both significantly restructured to accommodate all the features we have discussed (viscoelasticity, new brittle and viscous deformation mechanisms, brittle vs total strain, etc).

@gassmoeller
Copy link
Member

I agree. visco_plastic has reached a complexity that is ok for users, who have previously used diffusion_dislocation or drucker_prager, but additional complexity should be moved into a new material model (e.g. visco_elasto_plastic ?).

@anne-glerum
Copy link
Contributor Author

I've moved the code in interface to the visco_plastic header and source files. There could be another plastic model for which this output is interesting however, namely the dynamic_friction model, @alarshi?

@anne-glerum
Copy link
Contributor Author

Also, the output of the tests looks good: the cohesion and friction is weakened by the specified initial strain by the right amount in the right area. I'll run a 'real' model over time today to be sure.

@alarshi
Copy link
Contributor

alarshi commented May 24, 2017

@anne-glerum: Cool! I will try this after merge. Thank you!

Copy link
Member

@gassmoeller gassmoeller left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Two tiny things I missed yesterday, otherwise good to go.

@@ -747,6 +747,7 @@ namespace aspect
return vs;
}


Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Please remove this line

evaluate(const MaterialModel::MaterialModelInputs<dim> &in,
MaterialModel::MaterialModelOutputs<dim> &out) const
calculate_weakening(const double &strain_ii,
const unsigned int &j) const
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

hand over the doubles by-value (i.e. remove the &), for simple data types like double this is equally fast (or even faster), and the chance is lower that you accidentally change them inside of the function.

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

same for the unsigned int of course

@anne-glerum
Copy link
Contributor Author

If the dynamic_friction material model and the new, say visco_elasto_plastic, material model would also use this, should the PlasticAdditionalOutput not stay in interface.cc?

@gassmoeller
Copy link
Member

I do not feel very strongly about it. Usually we try to only put stuff in the interface class that is of importance for all (or a very significant fraction of all) material models. The fact that SeismicOutputs is in there is mostly a relic, because it was always there. If dynamic_friction and visco_elasto_plastic include the header of visco_plastic the current form should work as well, and it also represents a relationship between them. Even though they are not derived from each other, they do similar things, which is model visco_plastic behavior. I just try to avoid putting every new Output into the interface, which would grow uncontrollably over time, and it is very hard to get rid of anything in the interface (it took us years to get rid of the reference_... functions).

@naliboff
Copy link
Contributor

I do not feel strongly on this point either, but one reason to keep these functions within visco_plastic is there a multitude of permutations on how to implement strain weakening. These permutations can be based on the type of tracked strain (brittle vs total, etc) the yield mechanism (Drucker Prager, velocity-weakening, ..) etc. In many cases each material model may have a type of strain weakening similar to the current implementation, but slightly different nonetheless.

@anne-glerum
Copy link
Contributor Author

anne-glerum commented May 30, 2017

Alright, then I think it's done :) Just checked it on a rift model and calculating the current values in paraview gives the same results.

Create weakening function

Move PlasticAdditionalOuput to interface

Implement weakening function

Also use new function in viscosity calculation

Correct averaging and unweakened values

Add new output to test

Reformulate after rebase on master

Add new output to weakening tests

Add function description
@gassmoeller
Copy link
Member

Looks good.

@gassmoeller gassmoeller merged commit 97b6687 into geodynamics:master May 31, 2017
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

None yet

4 participants