ttau_plots.jl - seaice-experiments - sea ice experiments using Granular.jl
(HTM) git clone git://src.adamsgaard.dk/seaice-experiments
(DIR) Log
(DIR) Files
(DIR) Refs
(DIR) README
(DIR) LICENSE
---
ttau_plots.jl (6162B)
---
1 #!/usr/bin/env julia
2 import ArgParse
3 ENV["MPLBACKEND"] = "Agg"
4 using PyPlot
5
6 function plotTau(xvalues::Vector, tau::Vector, tau_s::Vector,
7 x_label::String, sim_id::String;
8 x_lim::Vector=[-1., -1.],
9 y_lim::Vector=[-1., -1.],
10 ycenter::Float64=-1,
11 shadex::Vector=[-1., -1.],
12 xvalues2::Vector=[-1.],
13 tau2::Vector=[-1.],
14 tau2_s::Vector=[-1.],
15 x_label2::String="",
16 x_lim2::Vector=[-1., -1.],
17 logx::Bool=false,
18 cohesive::Bool=false)
19
20 fig = figure(figsize=(5,4))
21 ax1 = gca()
22 if length(xvalues2) > 1
23 ax2 = ax1[:twiny]()
24
25 #new_position = [0.06;0.06;0.77;0.91] # Position Method 2
26 #ax2[:set_position](top=0.85)
27 end
28
29 if length(tau2) > 1
30 #semilogy(xvalues, tau, "xk", label="Frictional DEM")
31 #semilogy(xvalues, tau2, "+k", label="Cohesive DEM")
32 ax1[:errorbar](xvalues, tau, yerr=tau_s, fmt="xk",
33 label="Frictional DEM")
34
35 if length(xvalues2) > 1
36 ax2[:errorbar](xvalues2, tau2, yerr=tau2_s, fmt="+b",
37 label="Cohesive DEM")
38 ax1[:legend](fancybox="true", loc="lower right")
39 ax2[:legend](fancybox="true", loc="upper right")
40 else
41 ax1[:errorbar](xvalues, tau2, yerr=tau2_s, fmt="+b",
42 label="Cohesive DEM")
43 ax1[:legend](fancybox="true")
44 end
45 else
46 #semilogy(xvalues, tau, "ok")
47 if cohesive
48 errorbar(xvalues, tau, yerr=tau_s, fmt="+b")
49 else
50 errorbar(xvalues, tau, yerr=tau_s, fmt="xk")
51 end
52 end
53
54 yscale("log")
55 if logx
56 xscale("log")
57 end
58
59 if shadex[2] > 0.
60 ax1[:axvspan](shadex[1], shadex[2], color="gray", alpha=0.5, lw=0)
61
62 ax1[:text]((shadex[2] - shadex[1])/2. + shadex[1],
63 ycenter,
64 "No jamming observed",
65 horizontalalignment="center",
66 verticalalignment="center",
67 rotation="vertical")
68 end
69
70 for x in xvalues
71 ax1[:plot]([x, x], y_lim[1].*[0.9, 1.1], "-r", lw=4)
72 end
73 if length(xvalues2) > 1
74 for x in xvalues2
75 ax2[:plot]([x, x], y_lim[2].*[0.9, 1.1], "-r", lw=4)
76 end
77 end
78
79 #title(parsed_args["id"])
80 ax1[:set_xlabel](x_label)
81 ax1[:set_ylabel]("Time-scale for jamming \$T\$ [h]")
82 if x_lim[2] > 0.
83 ax1[:set_xlim](x_lim)
84 end
85 if y_lim[2] > 0.
86 ax1[:set_ylim](y_lim)
87 end
88
89 if length(xvalues2) > 1
90 font1 = Dict("color"=>"blue")
91 ax2[:set_xlabel](x_label2, fontdict=font1)
92 setp(ax2[:get_xticklabels](), color="blue")
93
94 if x_lim2[2] > 0.
95 ax2[:set_xlim](x_lim2)
96 end
97 end
98
99 tight_layout()
100
101 if length(xvalues2) > 1
102 savefig(sim_id * "_combined_tau.png")
103 savefig(sim_id * "_combined_tau.pdf")
104 else
105 savefig(sim_id * "_$(x_label)_" * "_tau.png")
106 savefig(sim_id * "_$(x_label)_" * "_tau.pdf")
107 end
108 end
109
110 # mu_sigmac
111
112 data = readdlm("mu_sigmac/mu_sigmac_fits.txt")
113 tau = data[:,2]
114 s = data[:,3]
115 plotTau([0.00, 0.10, 0.15, 0.20, 0.25, 0.30, 0.35, 0.40, 0.45],
116 tau[1:9],
117 s[1:9],
118 "Coulomb-friction coefficient \$\\mu\$ [-]", "mu_sigmac",
119 x_lim=[0.0, 0.5],
120 y_lim=[1e0, 2e2],
121 ycenter=1.05e1,
122 shadex=[0.0, 0.08])
123
124 plotTau(
125 [ 0., 10., 100., 200., 300., 400., 500., 600., 700., 800.],
126 vcat(tau[1], tau[10:end]),
127 vcat(s[1], s[10:end]),
128 "Tensile strength \$\\sigma_c\$ [kPa]", "mu_sigmac",
129 x_lim=[0, 1000.],
130 y_lim=[1e-1, 1e2],
131 ycenter=3e0,
132 shadex=[0., 150.],
133 cohesive=true)
134
135
136 ## combined mu_sigmac plot
137 plotTau([0.00, 0.10, 0.15, 0.20, 0.25, 0.30, 0.35, 0.40, 0.45],
138 tau[1:9],
139 s[1:9],
140 "Coulomb-friction coefficient \$\\mu\$ [-]",
141 "mu_sigmac",
142 #x_lim=[0.0, 0.5],
143 x_lim=[0.0, 0.9],
144 y_lim=[1e-1, 2e2],
145 ycenter=1.05e1,
146 shadex=[0.0, 0.09],
147 xvalues2=[ 0., 10., 100., 200., 300., 400., 500., 600., 700., 800.],
148 x_lim2=[50.0, 750.],
149 tau2=vcat(tau[1], tau[10:end]),
150 tau2_s=vcat(s[1], s[10:end]),
151 x_label2="Tensile strength \$\\sigma_c\$ [kPa]")
152
153
154 # psd_width
155 data = readdlm("psd_width/psd_width_fits.txt")
156 tau = data[:,2]
157 s = data[:,3]
158 plotTau(
159 [0., 1.15e3-8e2, 1.35e3-6e2, 1.55e3-4e2, 1.75e3-2e2]*2.,
160 tau[5:-1:1],
161 s[5:-1:1],
162 "PSD width \$d_{max} - d_{min}\$ [m]", "psd_width",
163 x_lim=[0, 3300.],
164 y_lim=[1e-0, 2e2],
165 ycenter=1e1,
166 shadex=[0., 500.],
167 tau2=tau[end:-1:6],
168 tau2_s=s[end:-1:6])
169
170 # stiffness
171 #=data = readdlm("stiffness/stiffness_fits.txt")
172 tau = data[:,2]
173 s = data[:,3]
174 plotTau([2e5, 2e6, 2e7, 2e8, 2e9],
175 tau[1:5],
176 s[1:5],
177 "Young's modulus \$E\$ [Pa]", "stiffness",
178 x_lim=[1e5, 1e10],
179 y_lim=[1e-1, 2e3],
180 ycenter=3e0,
181 shadex=[0., 0.],
182 tau2=tau[6:end],
183 tau2_s=s[6:end],
184 logx=true)=#
185
186 # thickness
187 data = readdlm("thickness/thickness_fits.txt")
188 tau = data[:,2]
189 s = data[:,3]
190 plotTau([1.0, 1.29, 1.67, 2.15, 2.78, 3.59, 4.64, 5.99, 7.74, 10.],
191 tau[1:10],
192 s[1:10],
193 "Ice-floe thickness \$h\$ [m]", "thickness",
194 x_lim=[0.9, 11.],
195 y_lim=[1e0, 1e3],
196 ycenter=3e0,
197 shadex=[0.1, 0.11],
198 tau2=tau[11:end],
199 tau2_s=s[11:end],
200 logx=true)
201
202 # width
203 data = readdlm("width/width_fits.txt")
204 tau = data[:,2]
205 s = data[:,3]
206 plotTau([5e3, 6e3, 7e3, 8e3, 9e3, 1e4],
207 tau[1:6],
208 s[1:6],
209 "Strait width \$W\$ [m]", "width",
210 x_lim=[4.5e3, 9e3],
211 y_lim=[1e0, 3e2],
212 ycenter=2e1,
213 shadex=[7.2e3, 9e3],
214 tau2=tau[7:end],
215 tau2_s=s[7:end],
216 logx=false)
217
218 # width_monodisperse
219 data = readdlm("width_monodisperse/width_monodisperse_fits.txt")
220 tau = data[:,2]
221 s = data[:,3]
222 plotTau([6e3, 7e3, 8e3, 9e3, 1e4],
223 tau[1:5],
224 s[1:5],
225 "Strait width \$W\$ [m]", "width_monodisperse",
226 x_lim=[4e3, 1.1e4],
227 y_lim=[1e-1, 3e2],
228 ycenter=3e0,
229 shadex=[0.1, 0.11],
230 tau2=tau[6:end],
231 tau2_s=s[6:end],
232 logx=false)
233