Note that you DO have to work on redoubt
today!
Today you'll put stuff on the crust and simulate the displacement field. To do this, you will use a program I wrote a while ago and will hopefully finish one of these days. CrusDe implements a generic convolution that enables a scientific version of plug-n-play for Green's functions (i.e. Earth Models), load models, time history functions and so forth.
After playing a bit the tool, you will familiarize yourself a bit with elastic and viscous parameters of the Earth and how they impact the displacement field.
Log into redoubt
, go into your working directory, and make a subdirectory for
today's lab:
$> cd $GEOP572 $> mkdir lab14 $> cd lab14
Now get some of the testcases that I provide with CrusDe:
$> cp -r $CRUSDE_HOME/testcases ./ $> cd testcases $> ls alma disk hekla icelandic_rhythmics multi_load
These directories contain some sample data and model configuration files that can be used with CrusDe. You will use some of these later.
Let's get a little bit of background on how to run CrusDe first; simply running
crusde
on the command line will give you a very limited
help:
$> crusde [ CrusDe ] : Welcome! Usage: crusde| > This program runs an experiment that is parameterized by an XML file, or starts a gui for either the experiment manager or the plug-in manager. guiOptions: -P|p run Plug-in manager. -P|p install <plug-in.so> Install 'plug-in.so' from command line. -M|m run Experiment manager. Options: -h|? show this help. -Q|q quiet run, minimal output. -v|V show version information.
The main point to take away from this is that you run crusde
with an
XML input file. Change into the disk
directory and look at the
disk_elastic
file:
$> cd disk $> less disk_elastic.xml
You will find that the file is made up of several sections: some header commentary,
and an experiment
that consists of a few subsections. If might be
worthwhile to follow along on this schematic, too:
file
tag whos value
attribute specifies the name of the output file that will contain the resultsregion
tags that specify the east, west, south and north borders of the model region (here in meters)gridsize
that specifies, in the same units as the region, the discretization of the model domainkernel
definition. All you need to know about this is that CrusDe currently implements 2
fast 2d convolution
fast 3d convolution
greens_function
. The function
we use here goes by the name elastic halfspace (pinel)
, other ones are listed below. This utilizes
the green's functions for an elastic halfspace derived by Pinel et al. (2007).
For this we need to specify 3 parameters: acc. gravity, Poisson's ratio, Youngs modulus.
load_function
that is made up of a load
and can contain
load_history
functions. Here, we use a disk load, later you may use irregular
loads that are specified in external files. This specific example requires center coordinates
of the disk, a density and a height and radius.xy2r
, which converts horizontal displacements
into radial displacements. and the the output
part which can either be netcdf
or table
. netcdf
creates a gridfile that can be viewed with ncview
, table
creates ASCII ouput.
(and yes, ideally the filename was actually defined in there, that will happen at some point in the future.)
$> crusde disk_elastic.xml [ CrusDe ] : Welcome! [ CrusDe ] : found region of interest coordinates to be in METERS: north=20000.000000, south=-20000.000000, east=20000.000000, west=-20000.000000 [ CrusDe ] : Plugin::registerParameter() for: /usr/local/CrusDe/plugins/green/pinel_hs_elastic.so [ CrusDe ] : Plugin::registerParameter() for: disk load My id: 0 [ CrusDe ] : Plugin::registerParameter() for: modeling kernel [ CrusDe ] : Plugin::registerParameter() for: output [ CrusDe ] : Plugin::registerParameter() for: xy2r [ CrusDe ] : elastic halfspace (pinel) registers output fields... [ CrusDe ] : x: 0 [ CrusDe ] : y: 1 [ CrusDe ] : z: 2 [ CrusDe ] : (xy2r) registers output field... [ CrusDe ] : (xy2r) ... got 3 [ CrusDe ] : Iterating through requested_plugin_list .... [ CrusDe ] : Number of requested plug-ins: 0 [ CrusDe ] : (fast 2d convolution) planning FFT ... [ CrusDe ] : Initialization complete. [ CrusDe ] : -------------------------------------------------------------------------------- [ CrusDe ] : Let me start your experiment. [ CrusDe ] : Working on job: elastic halfspace (pinel) (no job name defined) [ CrusDe ] : Model step: 0 (Model time: 0) [ CrusDe ] : convolution operator starts ... [ CrusDe ] : (fast 2d convolution) Convolution of Green's function and load ... [ CrusDe ] : postprocessors start ... [ CrusDe ] : (xy2r) wrote xy2r conversion to position: 3 [ CrusDe ] : result handler starts ... [ CrusDe ] : writing to file: ./disk_elastic_high.nc ... data_varid: 4, nc_id: 65536, NDIR: 4 (4), NLAT: 41, NLON: 41 [ CrusDe ] : ------------------------- [ CrusDe ] : Starting the experiment manager ... [ CrusDe ] : run time= 0.000000 s [ CrusDe ] : Terminating normally! Bye.
There's a bunch of stuff that happens, but fundamentally you are interested in seeing the Terminating normally! Bye.
statement at the end. Your directory should now contain a file disk_elastic_high.nc
, which is the file that
we defined above in the experiment definition. Let's look at the results:
$> ncview disk_elastic_high.nc
Note that there is now a direction value, click on the box under Current
to skip through the displacement
fields for x, y, z, r
directions. This should reproduce the images shown here; note that the horizontal displacement field
will look different as ncview scales to the maximum value of all results.
An example for a time series simulation might be running (will take between 30s - 60s):
$> crusde disk_fourier_timeseries.xml
You will see that the experiment definition file has 2 additional sections, and we're using the fast 3d convulution kernel.:
1) There are now parameters to specify the number of time steps and the size of the steps:
<!-- MODEL TIME CONFIG --> <parameter name="timestep_size" value="1"/> <!-- each step is a day --> <parameter name="timesteps" value="365"/> <!-- 365 steps == 1 year -->
2) In the load function, we now have a load_history
:
<!-- load history: use low order terms of fourier series --> <load_history> <plugin name="fourier_series" /> <parameter name="period_length" value="365"/> <!-- 365 days --> <parameter name="bias" value="0" /> <!-- vertical offset --> <parameter name="trend" value="0" /> <!-- linear trend --> <parameter name="annual_sin" value="0.3" /> <!-- see name --> <parameter name="annual_cos" value="0.3" /> <!-- see name --> <parameter name="semiannual_sin" value="0.3" /> <!-- see name --> <parameter name="semiannual_cos" value="0.3" /> <!-- see name --> </load_history>which allows us to start simulating seasonality in a more sophisticated way than by just assuming a sine.
$> ncview disk_fourier_ts.nc
Set the current direction to "2" and click the double forward arrow to run the simulation (>>). This should look rather hypnotic. Good! Stop whenever you're ready; feel free to look at other directions, though.
First go into the icelandic_rhythmics
folder and run the existing testcase:
$> cd $GEOP572/lab14/testcases/icelandic_rhythmics $> crusde elastic.xml $> ncview elastic.nc
First questions:
Now, what happens when you change the Young's modulus to 10 GPa, what at 80 GPa? Rerun the simulations with these new values, and give minimum and maximum displacements / screen shots of profiles. Do the results confirm the definition of Young's modulus (look it up if you don't know what it is.) How would a Young's modulus for a Craton compare to that of young Midocean Ridge crust, say in Iceland? While we're at it; what's a good definition for the Poisson's ratio?
For your interest, the load file defines 4 loads that are all uniform in height, everything else is zero (units are meter of water equivalent snow load. Yes, Vatnajökull gets about 2.5 m weq snow each winter.)
$> awk '{print $3}' load/load_all_glaciers.dat | sort -u 0.000000 1.250000 1.500000 1.650000 2.500000
Now set the gridsize
in the final_relaxed.xml
case to 10000m and run it. This brings some viscoelastic
relaxation into the game (this takes longer as the Green's functions are more complex; approx. 2 -5 mins)
$> crusde final_relaxed.xml $> ncview experiment.nc
The major new parameters for the Green's function are:
rho_f
- the density of ductile material,H
- the effective thickness of the elastic material that's on top of the viscoelastic material
What happens:
Do you expect these parameters could trade off with respect to the properties of the
final relaxed displacement field? Or would you expect unique solutions
for each of these parameters, i.e. you can uniquely resolve elastic plate thickness
vs viscosity? If so, why, if not, why not? (A few sentences with your substantiated
opinion will be enough, feel free to run simulations for other values, if you wish.)
H
to 1000, 50000? (with rho_f=3100)rho_f
to 2500, 3500? (with H=10000)
$> cd $GEOP572/lab14/testcases/icelandic_rhythmics $> crusde elastic_timeseries.xml $> ncview elastic_timeseries.nc
Hit play - Happy Thanksgiving!
The names below are the exact names you will have to use for your plugin names.
crustal_decay exponential exponential_rate data_handler netcdf writer table writer green elastic halfspace (pinel) thick plate (pinel) final relaxed (pinel) inverse thick plate (pinel) elastic minus thickplate (pinel) load disk load irregular load load_history dirac_impulse boxcar boxcar_rate fourier_series fourier_series_rate heaviside ramp sinusoidal sinusoidal_rate kernel fast 2d convolution fast 3d convolution postprocess factor xy2r
rg <at> nmt <dot> edu | Last modified: November 25 2015 21:46.