巖石側(cè)壓系數(shù)
Sheorey (1994) 基于一維層狀模型，考慮地球的曲率撒蟀、殼幔物質(zhì)屬性叙谨、溫度場和與熱相關(guān)的其他性質(zhì)，提出的一種計算地殼初始應(yīng)力的方法:
k = 0.25 + 7E(0.001 + 1/z)
下圖黑線是Tobias Hergert在2011年一篇論文中使用該公式畫出的曲線：
k-h
- k 為巖石的側(cè)壓系數(shù)保屯，是平均水平應(yīng)力(SH + Sh)/2 與垂向應(yīng)力SV 的比值(SH 和Sh 分別為最大和最小水平應(yīng)力)手负；
- E為巖石的楊氏模量（GPa），這里采用模型的平均巖石模量E=96 GPa, 該平均值依據(jù)CRUST1.0 計算得到配椭；
- z 為距離地表的深度虫溜，公式中采用以米為單位
GMT繪圖
當(dāng)E=95Ga，深度z以km為單位時股缸，公式可換算如下：
z = 0.665/(k-0.915)
#!/bin/bash
data=k-values.csv
gmt begin
gmt set FONT_ANNOT_PRIMARY=8p,4
gmt set FONT_LABEL=9p,4
gmt set MAP_TICK_LENGTH -2p
gmt set MAP_FRAME_PEN 0.5p
gmt set MAP_GRID_PEN_PRIMARY 0.3p,gray,-
gmt figure k-values pdf A0.2c,E800
gmt basemap -R0/3/0/8 -JX10c/-8c -Bxag0.5+l"k\075\050S@-H@-\053S@-h@-\051\0572S@-v@-" -Bya2fg1+l"depth(km)" -BWNrb
# 繪制實測應(yīng)力點
awk -F"," '{if(NR<=25) print $1,$2}' $data | gmt plot -Sc0.15c -W0.5p,blue -lKTB+jBR+gwhite+f8p+p0.5p,gray+o0.5p
awk -F"," '{if(NR>=26 && NR<=30)print $1,$2}' $data | gmt plot -Sd0.2c -W0.5p,red -lSAFOD
awk -F"," '{if(NR>=27) print $1, $2}' $data | gmt plot -Sx0.2c -W0.2p,gray -l"Global Dataset"
# 繪制測壓曲線
# k = 0.25 + 7E(0.001 + 1/z) E單位為GPa,Z單位為m;
# z = 0.665/(k-0.915) ; E=95GPa, z單位為km
gmt math -T0.915001/3/0.01 T 0.915 SUB INV 0.665 MUL = | gmt plot -W1p -l"E\07595GPa"
# k=1/3
gmt plot -W0.6p,- -l"k\0751/3" << EOF
0.333333 8
0.333333 0
EOF
gmt end
k-z
