Phase transition in ViscoPlastic model

Hi all,
I am trying to implement a phase transition in the ViscoPlastic model. However, I am having a problem with the prefactors for diffusion creep.
I originally wanted to directly set up the phase transition in the box without adding an additional field. However, for some reason, it couldn’t be achieved. So I added an additional field to achieve the phase transition.

  1. If I set it up like this:
set Prefactors for diffusion creep                          = background:1.25E-015,C_1:1.25E-015|6.12E-022|2.94E-022|5.4E-025

I got errors like this:
The string for property does not have
the expected format. Check that the string is either a comma separated
list of <double>' or : ||…, :
|… , … '. If the string looks correct, it is likely that
the length of the list of keys passed to parse_map_to_double_array
does not match the length of the comma separated property list.

  1. If I set it up like this:
set Prefactors for diffusion creep                          = background:1.25E-015|1.25E-015|1.25E-015|1.25E-015,C_1:1.25E-015|6.12E-022|2.94E-022|5.4E-025

I got errors like this:
The key in does not have
the expected number of values. It expects 1 or 1 values, but we found 4

  1. If I don’t set up the background prefactors for diffusion creep

set Prefactors for diffusion creep = C_1:1.25E-015|6.12E-022|2.94E-022|5.4E-025

I got errors like this:
The keyword in Prefactors for diffusion creep is not
listed, although it is expected. Check that you have at least one
value for each field id in your list (possibly plus background if a
background field is expected for this property).

Is there a way for me to set up the phase transition in the ViscoPlastic model? I would appreciate it if someone could help me! Thank you very much!


Hi Ziqi,

The syntax for setting up properties for different compositional fields and phase transitions is a bit tricky. But first it would be good to know why you need an additional compositional field. Having an additional field with different properties gives you a different physical behavior that adding a phase transition.

In addition, can you post the whole input file that does not work? Because for reading in material properties (like the Prefactors for diffusion creep) it is important how many compositional fields and how many phase transitions you have. You always need as many entries in your list as you have phases, in other words, one more than you have phase transitions. And you always need to list the properties for all compositional fields you have. So from your three examples, the second one looks most like what I would expect this to look like in the input file.

In your example, it seems like ASPECT only expects one value (and not the 4 values you’ve provided). So did you also add the lines that tell the material model that there should be three phase transitions?
This would be the parameters like Phase transition depths, Phase transition widths, Phase transition temperatures and Phase transition Clapeyron slopes.

You can see an example of how this looks like in the visco_plastic_phases_viscosity_max_min.prm test.


Hi Juliane,
Thank you very much for your reply! Your example is very helpful, and I have solved the problem.
Now I want to implement the 300 kg/m density jump at 410 and 400 kg/m at 660 by setting up Densities because I couldn’t find the phase transition density jump in the viscoplastic material model parameters. However, as a result, the plume couldn’t penetrate the 660 transition zone.

Below is part of my code:

subsection Compositional fields
  set Names of fields           = C_1
  set Number of fields          = 1
  set Compositional field methods = field

# Initial composition model
subsection Initial composition model
    set List of model names = function
    subsection Function
        set Coordinate system = cartesian
        set Function expression = 0.0

subsection Material model
  set Model name = depth dependent
  subsection Depth dependent model
    set Base model                                    = visco plastic
    set Depth dependence method                       = File
    set Data directory                                = /home/mazq/software/aspect/aspect/data/material-model/depth-dependent/
    set Viscosity depth file                          = visc_depth.txt
    set Reference viscosity                           = 1e22

  subsection Visco Plastic
    set Maximum viscosity                             = 1e24
    set Minimum viscosity                             = 1e18
    set Thermal conductivities                        = 4

    set Phase transition depths                       = background:410e3|520e3|670e3, C_1: 410e3|520e3|670e3
    set Phase transition widths                       = background:5e3|5e3|5e3, C_1: 5e3|5e3|5e3
    set Phase transition temperatures                 = background:1662.0|1662.0|1662.0, C_1: 1662.0|1662.0|1662.0
    set Phase transition Clapeyron slopes             = background:4e6|4.1e6|-2.7e6, C_1: 4e6|4.1e6|-2.7e6

    set Prefactors for diffusion creep                = background:1.25e-15|6.12e-22|2.94e-22|5.4e-25, C_1:1.25e-15|6.12e-22|2.94e-22|5.4e-25
    set Activation energies for diffusion creep       = background:375e3, C_1:375e3|231e3|270e3|299e3
    set Activation volumes for diffusion creep        = background:6e-6, C_1:6e-6|6e-6|6e-6|1.5e-6

    set Prefactors for dislocation creep              = background: 8.33e-17, C_1:8.33e-17|2.05e-12|2.05e-19|1.e-40
    set Activation energies for dislocation creep     = background:530e3, C_1:530e3
    set Activation volumes for dislocation creep      = background:1.4e-5, C_1:1.4e-5|1.7e-5|1.7e-5|0

    set Densities                                     = background: 3400|3700|3700|4100, C_1:3400|3700|3700|4100
    set Heat capacities                               = background: 3300.0|3394.4|3453.2|3617.6, C_1: 3235.0|3372.3|3441.7|3680.8

Then I tried changing the density of the background to 3400, and it worked. The plume can penetrate the 660.

set Densities     = background: 3400, C_1:3400|3700|3700|4100

But I don’t understand why. Could you please help me make sure that I understand correctly how the different fields interact with each other when you are available? The reason why the phase transition in another field affects the plume dynamics in the background field is that another field means another particle field, am I right? And the properties carried by particles are interpolated back to the background field in each time step. Do I understand correctly? Thank you very much!


Hi Ziqi,

Based on what you’ve sent, it looks like your compositional field C_1 is zero everywhere (in your Initial composition model, you set the composition to zero). So your model just has the background properties everywhere. You can see that in the second plot you sent, the density does not have jumps at the phase transition anymore. That is why you change the model when you change the background field. If you want the plume to have a different composition, you need to set that in the initial and/or boundary composition model.

If your plume does not penetrate 660, that is just what happens for that parameter combination. You could either reduce the Clapeyron slope or density jump, make the plume bigger/more buoyant, or change some other parameter.

And no, you selected “field” under your Compositional field methods, so your field is not tracked on particles. ASPECT instead solves an equation for advecting that field.


Hi Juliane,
Thank you very much for your helpful reply! At this stage, we are trying to implement realistic (or at least as much as possible) composition for the mantle transition zone before moving on to whether to consider different plume compositions. I will try different parameter combinations, as you suggested until the plume can penetrate the 660.
Thank you for explaining the compositional field methods! I read the corresponding part of this parameter in the manual and have a clearer understanding of it.
I do have another question. When I set up the phase transition depth to be 410, 520 and 660 km, however, according to the specific heat result, the transitional phases seem to be at depths of 330, 440, 740 km (my model y extent is 1320 km).

set Phase transition depths = background:410e3|520e3|670e3, C_1: 410e3|520e3|670e3

I notice the 410 km is in fact at 330 km (80 km less), 520 km is in fact at 440 km (80 km less);the 670 is in fact 740km (70 km more). I am a bit confused about this.
Nonetheless, I changed my phase transition depth, increased the 410 km by 80 km to 490 km, increased the 520 km by 80 km to 600 km, and decreased the 670 km by 70 km to 600 km. Then I got the result that I wanted. I wonder why.

set Phase transition depths = background:490e3|600e3|600e3, C_1: 490e3|600e3|600e3

Sorry for asking so many questions! Thank you for going through all my questions so far!


Hi Ziqi,

Have you looked at the Phase transition temperatures parameter? If your temperature is different from the temperature you specify in that parameter, and your Clapeyron slope is not zero, your phase transition will shift up or down.


Hi Juliane,
Got it! Thank you very much!