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.

Sign up for GitHub

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

Fix missing in.current_cell parameter entry. #6240

Open
wants to merge 4 commits into
base: main
Choose a base branch
from
Open

Fix missing in.current_cell parameter entry. #6240

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
9bb9def
Fix missing in.current_cell parameter entry.
MFraters Feb 20, 2025
0b9af89
add test case
danieldouglas92 Feb 20, 2025
9e6e137
Add changelog entry for cpo in.current_cell bug.
MFraters Feb 25, 2025
daa045b
fix test results.
MFraters Feb 25, 2025
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 5 additions & 0 deletions doc/modules/changes/20250225_mfraters
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,5 @@
Fixed: Adds current cell to material model intpust
Copy link
Member

Choose a reason for hiding this comment

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

Suggested change
Fixed: Adds current cell to material model intpust
Fixed: Add the current cell to the material model inputs

in the cpo particle property to fix the Olivine: Karato
2008 mineral type. Also added a test to test this feature.
Copy link
Member

Choose a reason for hiding this comment

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

Suggested change
2008 mineral type. Also added a test to test this feature.
2008 mineral type. Also added a test for this feature.

<br>
(Menno Fraters and Daniel Douglas, 2025/02/25)
2 changes: 2 additions & 0 deletions include/aspect/particle/property/crystal_preferred_orientation.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -217,6 +217,7 @@ namespace aspect
const SymmetricTensor<2,3> &strain_rate_3d,
const Tensor<2,3> &velocity_gradient_tensor,
const Point<dim> &position,
const typename DoFHandler<dim>::active_cell_iterator &cell,
const double temperature,
const double pressure,
const Tensor<1,dim> &velocity,
Expand Down Expand Up @@ -287,6 +288,7 @@ namespace aspect
DeformationType
determine_deformation_type(const DeformationTypeSelector deformation_type_selector,
const Point<dim> &position,
const typename DoFHandler<dim>::active_cell_iterator &cell,
const double temperature,
const double pressure,
const Tensor<1,dim> &velocity,
Expand Down
6 changes: 6 additions & 0 deletions source/particle/property/crystal_preferred_orientation.cc
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -298,6 +298,7 @@ namespace aspect
}

ArrayView<double> data = particle.get_properties();
const typename DoFHandler<dim>::active_cell_iterator cell(*particle.get_surrounding_cell(),&(this->get_dof_handler()));

for (unsigned int mineral_i = 0; mineral_i < n_minerals; ++mineral_i)
{
Expand All @@ -316,6 +317,7 @@ namespace aspect
strain_rate_3d,
velocity_gradient_3d,
particle.get_location(),
cell,
temperature,
pressure,
velocity,
Expand Down Expand Up @@ -571,6 +573,7 @@ namespace aspect
const SymmetricTensor<2,3> &strain_rate_3d,
const Tensor<2,3> &velocity_gradient_tensor,
const Point<dim> &position,
const typename DoFHandler<dim>::active_cell_iterator &cell,
const double temperature,
const double pressure,
const Tensor<1,dim> &velocity,
Expand All @@ -592,6 +595,7 @@ namespace aspect

const DeformationType deformation_type = determine_deformation_type(deformation_type_selector[mineral_i],
position,
cell,
temperature,
pressure,
velocity,
Expand Down Expand Up @@ -836,6 +840,7 @@ namespace aspect
DeformationType
CrystalPreferredOrientation<dim>::determine_deformation_type(const DeformationTypeSelector deformation_type_selector,
const Point<dim> &position,
const typename DoFHandler<dim>::active_cell_iterator &cell,
const double temperature,
const double pressure,
const Tensor<1,dim> &velocity,
Expand Down Expand Up @@ -872,6 +877,7 @@ namespace aspect
material_model_inputs.velocity[0] = velocity;
material_model_inputs.composition[0] = compositions;
material_model_inputs.strain_rate[0] = strain_rate;
material_model_inputs.current_cell = cell;

MaterialModel::MaterialModelOutputs<dim> material_model_outputs(1,this->n_compositional_fields());
this->get_material_model().evaluate(material_model_inputs, material_model_outputs);
Expand Down
147 changes: 147 additions & 0 deletions tests/cpo_karato_simple_shearbox.prm
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -0,0 +1,147 @@
set Dimension = 3
set Pressure normalization = surface
set Surface pressure = 0
set Nonlinear solver scheme = single Advection, no Stokes
set End time = 1e5
set Use years in output instead of seconds = false
set Output directory = cpo_karato_simple_shearbox
set Maximum time step = 1e4

subsection Compositional fields
set Number of fields = 1
set Names of fields = water
set Compositional field methods = static
end

subsection Initial composition model
set Model name = function

subsection Function
set Variable names = x,y,z
set Function expression = 1000*x + 500
end
end

subsection Gravity model
set Model name = vertical

subsection Vertical
set Magnitude = 10
end
end

subsection Geometry model
set Model name = box

subsection Box
set X extent = 1
set Y extent = 1
set Z extent = 1
set Box origin X coordinate = -0.500
set Box origin Y coordinate = -0.500
set Box origin Z coordinate = -0.500
end
end

subsection Initial temperature model
set Model name = function

subsection Function
set Function expression = 1200 ## Annotation: Temperature function
end
end

subsection Boundary temperature model
set List of model names = box
set Fixed temperature boundary indicators = bottom, top, left, right, front, back

subsection Box
set Bottom temperature = 1200
set Left temperature = 1200
set Right temperature = 1200
set Top temperature = 1200
set Front temperature = 1200
set Back temperature = 1200
end
end

subsection Prescribed Stokes solution
set Model name = function

subsection Velocity function
set Variable names = x,y,z,t
set Function expression = z*1e-5;0;0 ## Annotation set velocity condition
end
end


# Set the viscosity to be 1.75e13, which results in a differential stress
# of 350 MPa when combined with the prescribed shear conditions.
subsection Material model
set Model name = visco plastic
subsection Visco Plastic
set Maximum viscosity = 1.75e13
set Minimum viscosity = 1.75e13
set Densities = 3300
end
end

# No need to refine the mesh in this case
subsection Mesh refinement
set Initial global refinement = 0
set Initial adaptive refinement = 0
set Time steps between mesh refinement = 10000
end



subsection Solver parameters
set Temperature solver tolerance = 1e-10
end

subsection Particles
set List of particle properties = crystal preferred orientation
set Particle generator name = reference cell

subsection Generator
subsection Reference cell
set Number of particles per cell per direction = 5
end
end

subsection Crystal Preferred Orientation
set Random number seed = 301
set Number of grains per particle = 5
set Property advection method = Backward Euler
set Property advection tolerance = 1e-15
set CPO derivatives algorithm = D-Rex 2004

subsection Initial grains
set Minerals = Olivine: Karato 2008
set Volume fractions minerals = 1.0
end

subsection D-Rex 2004
set Mobility = 125
set Stress exponents = 3.5
set Exponents p = 1.5
set Nucleation efficiency = 5
set Threshold GBS = 0.3
end
end

subsection CPO Bingham Average
set Random number seed = 200
end
end


subsection Postprocess
set List of postprocessors = particles

subsection Particles
set Time between data output = 1e5
set Data output format = gnuplot
set Exclude output properties = a_cosine_matrix, volume fraction, integrated strain invariant, crystal preferred orientation, cpo elastic tensor, cpo bingham average, elastic tensor decomposition
end
end
Loading