[CIG-SHORT] Spherical 3d gravity initial stress
Demian Gomez
demiang at gmail.com
Wed Nov 9 08:55:12 PST 2016
Brad,
I'm analyzing the far-field (> 800 km from a thrust fault) co-seismic
displacements and the effect of the asthenospheric stiffness in
the lithospheric stress guide. It's mostly a "theory" analysis and I cannot
neglect the effect of sphericity since this would be highly criticized by
reviewers.
I am trying to implement a spatial grid db as you suggested, but I am
having a few problems. First for all, for some reason Pylith is crashing
with my spatial grid with the following error:
RuntimeError: Error occurred while reading spatial database file
'spatialdb/initial_stress-layer1.spatialdb'
I/O error while reading SimpleGridDB data.
I thought the problem could be the size of the spatialdb so I reduced the
number of
points but the problem persisted.
I then downloaded a traction db example from a benchmark (
https://github.com/geodynamics/pylith_benchmarks/tree/
master/dynamic/scecdynrup/tpv102) to compare it against my file. I couldn't
see any remarkable differences (except the coordinate system, number of
points, etc) so I used my header (modifying num-x, num-y and num-z) with
the data portion of the benchmark spatial db. Even though the number of
values (columns) in this spatial db do not agree with the 6 stress
components that I need (the benchmark uses 7 values) and the values itself
are completely meaningless in terms of stress, Pylith does now pass the
material initialization. It then crashes with the second layer where I left
the spatial db unmodified.
There seems to be some kind of problem with my data portion of the spatial
db, but I can't find it. I looked for CRLF/LF problems, field sizes, etc
just in case, but found no differences. Attached you will find the spatial
grid db that crashes Pylith. I'd appreciate if you could take a look at it
to see if there is anything obviously wrong.
Demián
--
*Dr. Demián D. Gómez*
Postdoctoral Researcher
The Ohio State University - School of Earth Sciences
275 Mendenhall Laboratory
125 South Oval Mall
Columbus, Ohio 43210
Cell: +1 (901) 900-7324
email: gomez.124 at osu.edu
On Wed, Nov 9, 2016 at 9:50 AM, Brad Aagaard <baagaard at usgs.gov> wrote:
> Demian,
>
> The code for initializing the initial stress and strain fields is in the
> method ElasticMaterial::_initializeInitialStress() in
> libsrc/pylith/materials/ElasticMaterial.cc.
>
> Before you head in that direction, I would seriously consider the size of
> the approximation in using a flat domain compared to a spherical domain
> relative to the other assumptions in your problem. In other words, do your
> have sufficient resolution in all of your other parameters that you know
> you need to use a domain with a spherical cap?
>
> Regards,
> Brad
>
>
>
> On 11/08/2016 04:08 PM, Demian Gomez wrote:
>
>> Hi Brad,
>>
>> Thanks for your answer.
>>
>> Would it be possible to modify Pylith's code to calculate the initial
>> tractions during runtime without having to specify them through a
>> spatial database? I can try to do it if you can provide me with some
>> orientation on where I should start looking. It seems that for a model
>> constructed in ECEF coordinates, calculating stress_xx, stress_yy and
>> stress_zz should be trivial if I have access to the cell centroids and
>> the material density.
>>
>> Thanks again,
>> Demián
>>
>> --
>> *Dr. Demián D. Gómez*
>> Postdoctoral Researcher
>> The Ohio State University - School of Earth Sciences
>> 275 Mendenhall Laboratory
>> 125 South Oval Mall
>> Columbus, Ohio 43210
>> Cell: +1 (901) 900-7324 <tel:%2B1%20%28901%29%20900-7324>
>> email: gomez.124 at osu.edu <mailto:gomez.124 at osu.edu>
>>
>>
>>
>> On Mon, Nov 7, 2016 at 12:42 PM, Brad Aagaard <baagaard at usgs.gov
>> <mailto:baagaard at usgs.gov>> wrote:
>>
>> Demian,
>>
>> I think you can get things to work with a 1-D spatial database with
>> some minor changes.
>>
>> WARNING: This will only work if stress_xx = stress_yy = stress_zz
>> and stress_xy = stress_yz = stress_xz = 0. This is due to the fact
>> that the orientation of the gravity vector is different in the
>> projected coordinate system than it is in the geocentric coordinate
>> system.
>>
>> You should construct your model in a geocentric local Cartesian
>> coordinate system. This is a geocentric coordinate system with the
>> origin at the surface.
>>
>> The initial tractions should be in a projected (flat) geographic
>> coordinate system, so that stress(z) = integral of rho*g from 0 to z.
>>
>> As an aside, the SimpleGridDB provides much faster interpolation
>> than SimpleDB for a logically Cartesian grid because it can find the
>> relevant points without a global search. The points need to conform
>> to a grid, but the x, y, and z coordinates do not have to be spaced
>> uniformly.
>>
>> Regards,
>> Brad
>>
>>
>>
>> On 11/06/2016 01:16 PM, Demian Gomez wrote:
>>
>> Dear all,
>>
>> I am trying to turn on gravity in a quasi-static problem that
>> uses a 3D
>> spherical cap model with layers (see attached figure for
>> geometry and
>> displacement field without gravity). The coordinate system of
>> the mesh
>> is ECEF and in pylithapp.cfg I've defined coordsys =
>> spatialdata.geocoords.CSGeo. My spatial databases for material
>> properties, etc are defined using:
>>
>> num-locs = 1
>> data-dim = 0
>> space-dim = 3
>> cs-data = geographic {
>> to-meters = 1.0
>> space-dim = 3
>> ellipsoid = sphere
>> datum-vert = ellipsoid
>> is-geocentric = true
>> } // cs-data
>>
>> Since the layers are uniform, I use db_properties.query_type =
>> nearest.
>>
>> To turn on gravity, I need to set the initial lithostatic stress =
>> g*rho*h for each element (I can't use a uniform database
>> anymore), but
>> because this is a large model (~15e6 elements), the sizes of the
>> initial
>> stress databases are > 300 MB which makes Pylith take a very
>> long time
>> to initialize the materials (so far I've waited almost 12 hours
>> and
>> Pylith was still initializing the first layer!).
>>
>> Following this thread:
>> http://lists.geodynamics.org/pipermail/cig-short/2014-July/0
>> 01818.html
>> <http://lists.geodynamics.org/pipermail/cig-short/2014-July/
>> 001818.html>
>>
>> I have now attempted to create an initial stress spatial
>> database that
>> is only 1D as follows:
>>
>> #SPATIAL.ascii 1
>> SimpleDB {
>> num-values = 6
>> value-names = stress-xx stress-yy stress-zz stress-xy
>> stress-yz stress-xz
>> value-units = Pa Pa Pa Pa Pa Pa // units
>> num-locs = 2
>> data-dim = 1
>> space-dim = 3
>> cs-data = geo-local-cartesian {
>> origin-lon = -72.6425
>> origin-lat = -35.6331
>> origin-elev = 0
>> to-meters = 1.0
>> ellipsoid = sphere
>> datum-vert = ellipsoid
>> } // cs-data
>> }
>> 0.0000 0.0000 0.0000 0.000
>> 0.000
>> 0.000 0.00 0.00 0.00
>> 0.0000 0.0000 -104911.1753 -2674950661.887
>> -2674950661.887
>> -2674950661 <tel:2674950661>.887 0.00 0.00 0.00
>>
>>
>> and setting db_initial_stress.query_type = linear.
>>
>> My idea was that Pylith would map the geo-local-cartesian from the
>> spatial database to ECEF coordinates in geocoords.CSGeo.
>> However, I got
>> the following error during execution:
>>
>> mpinemesis: spatialdb/SimpleDBQuery.cc:352: void
>> spatialdata::spatialdb::SimpleDBQuery::_findLinePt(std::vect
>> or<spatialdata::spatialdb::SimpleDBQuery::WtStruct>*):
>> Assertion `locIndexA >= 0' failed.
>>
>> My questions are:
>>
>> 1) Is it possible to mix coordinate systems the way I did?
>> 2) Is there another way to tell Pylith that the initial stresses
>> in my
>> model vary linearly with depth, even when I have an ECEF
>> coordinate
>> system? I would like to avoid having to assign the initial stress
>> to
>> every element one by one (which leads to a very long
>> initialization time).
>>
>> Any suggestions will be greatly appreciated. Let me know if you
>> need
>> more information.
>>
>> Cheers,
>> Demián
>>
>> --
>> *Dr. Demián D. Gómez*
>> Postdoctoral Researcher
>> The Ohio State University - School of Earth Sciences
>> 275 Mendenhall Laboratory
>> 125 South Oval Mall
>> Columbus, Ohio 43210
>> Cell: +1 (901) 900-7324 <tel:%2B1%20%28901%29%20900-7324>
>> email: gomez.124 at osu.edu <mailto:gomez.124 at osu.edu>
>> <mailto:gomez.124 at osu.edu <mailto:gomez.124 at osu.edu>>
>>
>>
>>
>> _______________________________________________
>> CIG-SHORT mailing list
>> CIG-SHORT at geodynamics.org <mailto:CIG-SHORT at geodynamics.org>
>> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/cig-short
>> <http://lists.geodynamics.org/cgi-bin/mailman/listinfo/cig-short>
>>
>>
>> _______________________________________________
>> CIG-SHORT mailing list
>> CIG-SHORT at geodynamics.org <mailto:CIG-SHORT at geodynamics.org>
>> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/cig-short
>> <http://lists.geodynamics.org/cgi-bin/mailman/listinfo/cig-short>
>>
>>
>>
>>
>> _______________________________________________
>> CIG-SHORT mailing list
>> CIG-SHORT at geodynamics.org
>> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/cig-short
>>
>>
> _______________________________________________
> CIG-SHORT mailing list
> CIG-SHORT at geodynamics.org
> http://lists.geodynamics.org/cgi-bin/mailman/listinfo/cig-short
>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.geodynamics.org/pipermail/cig-short/attachments/20161109/2e33c8ed/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: initial_stress-layer1.spatialdb
Type: application/octet-stream
Size: 76745 bytes
Desc: not available
URL: <http://lists.geodynamics.org/pipermail/cig-short/attachments/20161109/2e33c8ed/attachment-0001.obj>
More information about the CIG-SHORT
mailing list