[CIG-SHORT] Gravity on Large/Heterogeneous Meshes

[Charles Williams]

Scott,

Are you using linear interpolation for your initial stresses?

Charles


On 12/11/2012, at 4:43 PM, Scott Henderson wrote:

> Brad,
> 
> If I understand you correctly, I should bring the tolerances closer together, so instead of:
> ksp_rtol = 1.0e-10
> ksp_atol = 1.0e-20
> 
> I've tried (among other combinations):
> ksp_rtol = 1.0e-18
> ksp_atol = 1.0e-20
> ksp_max_it = 200
> 
> And I always converge to the same result (-0.017m surface subsidence).
> 
> Could this have something to do with the way the initial_stresses are getting assigned? Basically, I'm thinking that if the bottom elements are getting initial stresses slightly less than lithostatic there would be subsidence of the domain. For the heck of it I tried specifying the max stress at 99km instead of 100km. I thought this would result in an error, but Pylith runs to completion, and the displacements turn out to be +17m. 
> 
> Thanks,
> Scott
> 
> 
> 
> On Nov 11, 2012, at 7:41 PM, Brad Aagaard <baagaard at usgs.gov> wrote:
> 
>> Scott,
>> 
>> There are a couple of issues to consider:
>> (1) If you have nonuniform density, then the initial stresses as a 
>> function of depth are more complicated than a constant linear increase 
>> with depth. You need to calculate the lithostatic stress as a function 
>> of depth and use that for the initial stresses.
>> (2) Tighten the convergence tolerance (ksp_rtol and ksp_atol) to get a 
>> more accurate solution.
>> 
>> Regards,
>> Brad
>> 
>> 
>> On 11/11/12 4:25 PM, Scott Henderson wrote:
>>> Hello,
>>> 
>>> I've been modifying the step16.cfg tutorial for a larger mesh. If I
>>> increase the domain size from 6x6x4 to 200x200x100km, and change the
>>> initial_stress.spatialdb to correspond to the new lithostatic stress
>>> at the base of the mesh (-2451637500Pa). I'm finding that the elastic
>>> solution has ~2 cm of subsidence at the surface, but I expect zero
>>> displacement. This appears to be independent of element size.
>>> 
>>> Also, if I have non-uniform density in the mesh (for example an
>>> embedded weak layer or arbitrary shape), how might I set up the
>>> initial stresses to compensate for gravity? I'm guessing that the
>>> only way to go about this is run the solution without initial
>>> stresses, extract the stress values at a few points through each
>>> material, then set up separate spatial databases for each material
>>> for the next run?
>>> 
>>> Thanks, Scott
>>> 
>>> step16_1 = same settings as step16, but larger mesh step16_2 = just
>>> elastic solution
>>> 
>>> 
>>> 
>>> 
>>> 
>>> 
>>> 
>>> _______________________________________________ CIG-SHORT mailing
>>> list CIG-SHORT at geodynamics.org
>>> http://geodynamics.org/cgi-bin/mailman/listinfo/cig-short
>>> 
>> 
>> _______________________________________________
>> CIG-SHORT mailing list
>> CIG-SHORT at geodynamics.org
>> http://geodynamics.org/cgi-bin/mailman/listinfo/cig-short
> 
> _______________________________________________
> CIG-SHORT mailing list
> CIG-SHORT at geodynamics.org
> http://geodynamics.org/cgi-bin/mailman/listinfo/cig-short

Charles A. Williams
Scientist
GNS Science
1 Fairway Drive, Avalon
PO Box 30368
Lower Hutt  5040
New Zealand
ph (office): 0064-4570-4566
fax (office): 0064-4570-4600
C.Williams at gns.cri.nz

-------------- next part --------------
An HTML attachment was scrubbed...
URL: http://geodynamics.org/pipermail/cig-short/attachments/20121112/5d8b2aee/attachment.htm