Question on mantle melting beneath mid-ocean ridges with a constant melt gradient

ASPECT
1

Good morning all,

I have a simulation adapted from ’mid_ocean_ridge.prm’ in cookbooks using Aspect 2.5.0. This numerical model based on the framework of Ribe (1985), designed to simulate one-dimensional mantle upwelling and melting beneath a mid-ocean ridge.

To simplify the model and isolate the mechanical processes, I eliminated thermal effects and adopted a constant melt gradient. I altered the boundary velocity and the initial temperature/material models, disabled thermal conduction and thermal diffusion, and the temperature is held constant at 0°C throughout the domain to isolate the decompression melting process. The goal is to have material entering from the bottom and moving upward, generating melt with a linear fraction gradient.

However, I am encountering severe numerical instabilities right from the first timestep:

  • Widespread negative pressure (P < 0) across the domain.

  • Unphysical upward velocity in the upper matrix (which should be fixed or have near-zero velocity depending on my boundary conditions).

  • Excessively high porosity values that exceed physically plausible bounds.

  • In later timesteps, the porosity field develops sharp, oscillatory oscillations instead of a smooth, increasing gradient.

I have meticulously reviewed my code and model setup multiple times but cannot pinpoint the root cause of these anomalies.I have also attached here the abnormal screenshot and both the original input file ‘mid_ocean_ridge.prm’ and my modified version ‘melt test.prm‘ for reference. Any help in comparing them or identifying the source of the instability would be invaluable.

porosity
porosity1268×1388 22.1 KB

​​I would be extremely grateful if someone could shed some light on why this is happening and how I might adjust my model to achieve a stable, physical solution.​​

Best regards,

Cedric Sage

mid_ocean_ridge.prm (10.0 KB)

melt test.prm (4.1 KB)

2

Hi @Cedric_Sage,

Welcome and thank you for posting to the forum.

A quick question to start - if you remove the melting components from melt test.prm (i.e., only include solid dynamics), do you still observe negative pressures and unphysical velocities?

Your boundary conditions looked reasonable, but perhaps it would be good start with your solid model setup and use something the simple material model to verify it is behaving as expected.

In you have not seen it, the 1D water solubility benchmark has a setup similar to what you are trying to achieve. That may be another starting point to work from.

John

3

Hi John@jbnaliboff,
Thank you so much for your helpful guidance.

Following your advice, I removed the melting components and ran the simulation with only solid dynamics. The results were very encouraging: ​​all parameters, including pressure and velocity, behaved normally and physically.​​ This suggests that the root cause of the anomalies is indeed isolated to the melting-related setup or processes.

Also, thank you for pointing me to the ​​1D water solubility benchmark​​. I have tried to adapt my boundary conditions and other settings based on its example. However, despite these adjustments, the problem persists. I am still encountering the same frustrating issues with widespread negative pressures and unphysical porosity oscillations when melting is enabled.

This leads me to my next question: I am specifically using a ​​constant melting gradient​​ to simplify the thermal mechanics. Could this specific approach be interacting with the melt migration physics in an unexpected way that might cause these numerical instabilities? For instance, could the rapid, prescribed generation of melt be creating conditions that are difficult for the solver to resolve?

Do you have any deeper insights or experiences regarding the implementation of such a constant melt production rate? Any further pointers on how to stabilize this would be immensely appreciated.

Thank you again for your time and suggestions.

Best regards,

Cedric_Sage

4

Hi @Cedric_Sage:

Thank you for performing the tests without melting and glad to hear things are behaving as expected in those models.

For instance, could the rapid, prescribed generation of melt be creating conditions that are difficult for the solver to resolve?

This certainly could be case.

Do you have any deeper insights or experiences regarding the implementation of such a constant melt production rate? Any further pointers on how to stabilize this would be immensely appreciated.

I have a few suggestions on how to proceed:

  1. Try using very small time steps and and/or a very small prescribed melting rate to see if this improves the behavior. Slowly ramping up the melting rate through time may also help.
  2. Try modifying the setup so that the melting rate is driven strictly by the P-T conditions, which give the desired melting rate following a melting parameterization of choice. Note that you can setup an ASPECT model such that the temperature is held constant through time.

If none of these suggestions help, I think we should dive further into the details of what is happening on the first problematic time steps.

Cheers,

John