[aspect-devel] Morency-Dorin with free surface

[Payman Janbakhsh]

Thanks Jonathan,
I have tried no-slip as well as tangential side boundary velocities but not much difference. I always keep surface vertical projection set to vertical when choosing top free surface. I even increased temperature polynomial degree to 3 to increase the accuracy.
Then I set thermal expansivity to zero to simplify the problem a bit. so it kept going for few time steps then froze again trying to solve Stokes systems.
I have noticed this happens when RMS velocities starts going up drastically . I have to figure out why velocities become so unphysical.
The whole purpose of using this material model was to incorporate a visco-plastic rheology to 3D problems. So need to find my way around to get there.
Any help is much appreciated from anyone.

Cheers

Payman

-----Original Message-----
From: Aspect-devel [mailto:aspect-devel-bounces at geodynamics.org] On Behalf Of Jonathan Perry-Houts
Sent: October 21, 2016 8:33 PM
To: aspect-devel at geodynamics.org
Subject: Re: [aspect-devel] Morency-Dorin with free surface

Hi Payman,

I haven't tried the Morency-Doin with a free surface, but it's caused me trouble in the past in relatively simple models with similar setups to the original Morency & Doin paper on which it's based. -- Point being that the material model is finnicky, and I'm not surprised it's having trouble.

A few thoughts:

Maybe no-slip on the sides is not a good choice, since the top is unconstrained? This shouldn't be a problem, but who knows. Maybe worth playing with?

Using:
>  subsection Free surface
>    set Surface velocity projection = vertical  end
sometimes helps stabilize free surfaces. Failure on the first time step makes this an unlikely cause for your problem, but maybe worth a try.

You mentioned that your model is 2D, so this isn't a problem for now, but keep in mind that the Morency-Doin model shouldn't be applied to 3D models. The formation of the second invariant of the strain rate tensor in the paper has never made sense to me, and so it's not generalized to three dimensions in the code (which reminds me, I should add an assertion that makes sure this is the case. Not sure how I missed that
before!) -- Any insight from someone else to what the authors meant by equation (4) in the paper (Morency and Doin, “Numerical Simulations of the Mantle Lithosphere Delamination.”) would be helpful! (e.g. is it a weird index notation I'm not familiar with?)

I think that equation, as implemented in the code, is almost certainly wrong. I included it for consistency with the paper after discussions with Timo and Wolfgang like two years ago. Perhaps try replacing it with the actual second invariant and try again?
(
https://github.com/geodynamics/aspect/blob/master/source/material_model/morency_doin.cc#L104
)

Beyond that I don't know what to tell you. Your parameters sound reasonable from what you described. Keep us updated on what works! I'd love to get that material model working for more than just a cookbook!

Cheers,
JPH

On 10/21/2016 01:48 PM, Payman Janbakhsh wrote:
> Hi John
> 
> 2D , delamination, geometry aspect ration of 2x1 , free top surface, 
> other boundaries zero velocity.
> 
> Linear temperature profile within the upper 100km and lower 100km and 
> a constant T in between.
> 
> 25km upper crust, 30km lower crust, a weak zone within the lower 
> crust( with a lower activation energy compare to the enclosing lower 
> crust, so can set off the delamination).
> 
> Have tried different reasonable densities and activation energies for 
> every layer except mantle.
> 
> Also have put Boundary composition model as : Initial composition  . 
> not sure if this is correct.
> 
> Running with IMPES or iterated Stokes made no difference
> 
> Also velocity polynomial set to 2 and other polynomials set to 1
> 
> also have turned off Local conservative discretization
> 
>  
> 
> cheers,
> 
>  
> 
> payman
> 
>  
> 
> *From:*Aspect-devel [mailto:aspect-devel-bounces at geodynamics.org] *On 
> Behalf Of *John Naliboff
> *Sent:* October 21, 2016 3:36 PM
> *To:* aspect-devel at geodynamics.org
> *Subject:* Re: [aspect-devel] Morency-Dorin with free surface
> 
>  
> 
> Hi Payman,
> 
> The Morency-Dorin material model should not have intrinsic problem 
> with handling a free surface.  The issue is likely related to some 
> other aspect of your model design.  What type of processes are you modeling?
> 
> Cheers,
> John
> 
> 
> *************************************************
> 
> John Naliboff
> 
> Assistant Project Scientist, CIG
> 
> Earth & Planetary Sciences Dept., UC Davis
> 
> On 10/21/2016 12:17 PM, Payman Janbakhsh wrote:
> 
>     Hi
> 
>     Does Morency-Dorin material model have difficulty with free surface?
> 
>     My models get stuck at solving Stokes system even at step zero when
>     I switch top boundary to free surface.
> 
>     I have set CFL to 0.1 and theta to 0.7 even chanh3d max time step to
>     1000 but still same issue. I have all other boundaries set to no
>     slip boundaries
> 
>      
> 
>     Thanks
> 
>     Payman
> 
>      
> 
>      
> 
>      
> 
>     Sent from my Samsung Galaxy smartphone.
> 
> 
> 
> 
>     _______________________________________________
> 
>     Aspect-devel mailing list
> 
>     Aspect-devel at geodynamics.org <mailto:Aspect-devel at geodynamics.org>
> 
>     http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
> 
>  
> 
> 
> 
> _______________________________________________
> Aspect-devel mailing list
> Aspect-devel at geodynamics.org
> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel
> 
_______________________________________________
Aspect-devel mailing list
Aspect-devel at geodynamics.org
http://lists.geodynamics.org/cgi-bin/mailman/listinfo/aspect-devel