Skip to content

Commit 8bebbc8

Browse files
bobmyhillbobmyhill
committed
updated tests
1 parent c14e7bd commit 8bebbc8

File tree

7 files changed

+387
-9
lines changed

7 files changed

+387
-9
lines changed

‎Readme.md‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,6 @@
11
[![DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.14238360.svg)](https://doi.org/10.5281/zenodo.14238360)
22
[![DOI](https://joss.theoj.org/papers/10.21105/joss.05389/status.svg)](https://doi.org/10.21105/joss.05389)
3-
[![coverage](https://img.shields.io/badge/coverage-%3E86%25-green")](https://github.com/geodynamics/burnman/actions/workflows/coverage.yml)
3+
[![coverage](https://img.shields.io/badge/coverage-%3E90%25-green")](https://github.com/geodynamics/burnman/actions/workflows/coverage.yml)
44

55
# BurnMan - a Python toolkit for planetary geophysics, geochemistry and thermodynamics
66

‎tests/test_eos.py‎

Lines changed: 13 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -477,7 +477,7 @@ def test_debye(self):
477477
self.assertAlmostEqual(Cv, Cv2)
478478
self.assertAlmostEqual(S, S2)
479479

480-
def test_pressure_finding(self):
480+
def test_pressure_finding_SLB(self):
481481
with warnings.catch_warnings(record=True) as w:
482482
warnings.simplefilter("always")
483483
m = alpha_quartz()
@@ -488,6 +488,18 @@ def test_pressure_finding(self):
488488
m.set_state_with_volume(V, T)
489489
self.assertAlmostEqual(m.pressure / 1.0e5, P / 1.0e5, places=3)
490490

491+
def test_pressure_finding_HP(self):
492+
with warnings.catch_warnings(record=True) as w:
493+
warnings.simplefilter("always")
494+
m = minerals.HP_2011_ds62.mt()
495+
P = 1.0e9
496+
T = 1000.0
497+
m.set_state(P, T)
498+
V = m.V
499+
m.set_state(P * 1.01, T)
500+
m.set_state_with_volume(V, T)
501+
self.assertAlmostEqual(m.pressure / 1.0e5, P / 1.0e5, places=3)
502+
491503

492504
if __name__ == "__main__":
493505
unittest.main()

‎tests/test_eos_consistency.py‎

Lines changed: 39 additions & 7 deletions
Original file line numberDiff line numberDiff line change
@@ -11,9 +11,41 @@ class EosConsistency(BurnManTest):
1111
def test_HP(self):
1212
P = 10.0e9
1313
T = 3000.0
14+
per = burnman.minerals.HP_2011_ds62.per()
1415
self.assertEqual(
1516
check_eos_consistency(
16-
burnman.minerals.HP_2011_ds62.per(),
17+
per,
18+
P,
19+
T,
20+
including_shear_properties=False,
21+
),
22+
True,
23+
)
24+
25+
def test_HP98(self):
26+
P = 10.0e9
27+
T = 3000.0
28+
29+
params = burnman.minerals.HP_2011_ds62.per().params
30+
params["equation_of_state"] = "hp98"
31+
params["dKdT_0"] = -1.0e7
32+
per = burnman.Mineral(params)
33+
self.assertEqual(
34+
check_eos_consistency(
35+
per,
36+
P,
37+
T,
38+
including_shear_properties=False,
39+
),
40+
True,
41+
)
42+
43+
def test_HP_liquid(self):
44+
P = 10.0e9
45+
T = 3000.0
46+
self.assertEqual(
47+
check_eos_consistency(
48+
burnman.minerals.HP_2011_ds62.foL(),
1749
P,
1850
T,
1951
including_shear_properties=False,
@@ -64,13 +96,13 @@ def test_DKS_liquid(self):
6496
def test_DKS_solid(self):
6597
P = 10.0e9
6698
T = 3000.0
99+
params = burnman.minerals.DKS_2013_solids.periclase().params
100+
params2 = burnman.minerals.SLB_2011.periclase().params
101+
for p in ["G_0", "Gprime_0", "eta_s_0"]:
102+
params[p] = params2[p]
103+
m = burnman.Mineral(params)
67104
self.assertEqual(
68-
check_eos_consistency(
69-
burnman.minerals.DKS_2013_solids.periclase(),
70-
P,
71-
T,
72-
including_shear_properties=False,
73-
),
105+
check_eos_consistency(m, P, T, including_shear_properties=True),
74106
True,
75107
)
76108

‎tests/test_fitting.py‎

Lines changed: 139 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -2,9 +2,13 @@
22
import unittest
33
from util import BurnManTest
44
import numpy as np
5+
from numpy import random
6+
import matplotlib.pyplot as plt
57

68
import burnman
9+
from burnman.optimize.eos_fitting import fit_XPTp_data
710
from burnman.optimize.nonlinear_fitting import nonlinear_least_squares_fit
11+
from burnman.utils.misc import attribute_function, pretty_string_values
812

913
path = os.path.dirname(os.path.abspath(__file__))
1014

@@ -141,8 +145,143 @@ def test_fit_bounded_PVT_data(self):
141145
fo, params, PTV, bounds=bounds, verbose=False
142146
)
143147

148+
cp_bands = burnman.nonlinear_fitting.confidence_prediction_bands(
149+
model=fitted_eos,
150+
x_array=PTV,
151+
confidence_interval=0.95,
152+
f=attribute_function(fo, "V"),
153+
flag="V",
154+
)
155+
self.assertEqual(len(cp_bands[0]), len(PTV))
156+
self.assertEqual(len(cp_bands), 4)
157+
144158
self.assertFloatEqual(3.0, fitted_eos.popt[2])
145159

160+
s = pretty_string_values(
161+
fitted_eos.popt,
162+
fitted_eos.pcov,
163+
extra_decimal_places=1,
164+
combine_value_and_sigma=False,
165+
)
166+
167+
self.assertEqual(len(s), 3)
168+
self.assertEqual(len(s[0]), 3)
169+
self.assertEqual(len(s[1]), 3)
170+
self.assertEqual(len(s[2]), 3)
171+
172+
s = pretty_string_values(
173+
fitted_eos.popt,
174+
fitted_eos.pcov,
175+
extra_decimal_places=1,
176+
combine_value_and_sigma=True,
177+
)
178+
179+
self.assertEqual(len(s), 3)
180+
self.assertEqual(len(s[0]), 3)
181+
self.assertEqual(len(s[1]), 3)
182+
self.assertEqual(len(s[2]), 3)
183+
184+
def test_bounded_solution_fitting(self):
185+
solution = burnman.minerals.SLB_2011.mg_fe_olivine()
186+
solution.set_state(1.0e5, 300.0)
187+
fit_params = [["V_0", 0], ["V_0", 1], ["V", 0, 1]]
188+
189+
n_data = 5
190+
data = []
191+
data_covariances = []
192+
flags = []
193+
194+
f_Verror = 1.0e-3
195+
196+
# Choose a specific seed so that the test is reproducible.
197+
random.seed(10)
198+
for i in range(n_data):
199+
x_fa = random.random()
200+
P = random.random() * 1.0e10
201+
T = random.random() * 1000.0 + 300.0
202+
X = [1.0 - x_fa, x_fa]
203+
solution.set_composition(X)
204+
solution.set_state(P, T)
205+
f = 1.0 + (random.normal() - 0.5) * f_Verror
206+
V = solution.V * f
207+
208+
data.append([1.0 - x_fa, x_fa, P, T, V])
209+
data_covariances.append(np.zeros((5, 5)))
210+
data_covariances[-1][4, 4] = np.power(solution.V * f_Verror, 2.0)
211+
212+
flags = ["V"] * 5
213+
214+
n_data = 2
215+
f_Vperror = 1.0e-2
216+
217+
for i in range(n_data):
218+
x_fa = random.random()
219+
P = random.random() * 1.0e10
220+
T = random.random() * 1000.0 + 300.0
221+
X = [1.0 - x_fa, x_fa]
222+
solution.set_composition(X)
223+
solution.set_state(P, T)
224+
f = 1.0 + (random.normal() - 0.5) * f_Vperror
225+
Vp = solution.p_wave_velocity * f
226+
227+
data.append([1.0 - x_fa, x_fa, P, T, Vp])
228+
data_covariances.append(np.zeros((5, 5)))
229+
data_covariances[-1][4, 4] = np.power(
230+
solution.p_wave_velocity * f_Vperror, 2.0
231+
)
232+
flags.append("p_wave_velocity")
233+
234+
data = np.array(data)
235+
data_covariances = np.array(data_covariances)
236+
flags = np.array(flags)
237+
delta_params = np.array([1.0e-8, 1.0e-8, 1.0e-8])
238+
bounds = np.array([[0, np.inf], [0, np.inf], [-np.inf, np.inf]])
239+
240+
fitted_eos = fit_XPTp_data(
241+
solution=solution,
242+
flags=flags,
243+
fit_params=fit_params,
244+
data=data,
245+
data_covariances=data_covariances,
246+
delta_params=delta_params,
247+
bounds=bounds,
248+
param_tolerance=1.0e-5,
249+
verbose=False,
250+
)
251+
252+
self.assertEqual(len(fitted_eos.popt), 3)
253+
254+
cp_bands = burnman.nonlinear_fitting.confidence_prediction_bands(
255+
model=fitted_eos,
256+
x_array=data,
257+
confidence_interval=0.95,
258+
f=attribute_function(solution, "V"),
259+
flag="V",
260+
)
261+
self.assertEqual(len(cp_bands[0]), len(data))
262+
self.assertEqual(len(cp_bands), 4)
263+
264+
good_data_confidence_interval = 0.9
265+
_, indices, probabilities = burnman.nonlinear_fitting.extreme_values(
266+
fitted_eos.weighted_residuals, good_data_confidence_interval
267+
)
268+
self.assertEqual(len(indices), 0)
269+
self.assertEqual(len(probabilities), 0)
270+
271+
# Just check plotting doesn't return an error
272+
fig, ax = plt.subplots()
273+
burnman.nonlinear_fitting.plot_residuals(
274+
ax=ax,
275+
weighted_residuals=fitted_eos.weighted_residuals,
276+
flags=fitted_eos.flags,
277+
)
278+
fig, ax = plt.subplots()
279+
burnman.nonlinear_fitting.weighted_residual_plot(ax, fitted_eos)
280+
281+
fig, ax = burnman.nonlinear_fitting.corner_plot(
282+
fitted_eos.popt, fitted_eos.pcov
283+
)
284+
146285

147286
if __name__ == "__main__":
148287
unittest.main()

‎tests/test_seismic.py‎

Lines changed: 29 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,8 +1,10 @@
11
from __future__ import absolute_import
22
import unittest
33
from util import BurnManTest
4+
import numpy as np
45

56
import burnman
7+
from burnman.tools import output_seismo
68

79

810
class test_seismic(BurnManTest):
@@ -84,6 +86,33 @@ def test_evaluate(self):
8486
# print "'%s': %s," % (name, result)
8587
self.assertArraysAlmostEqual(result, ref[name])
8688

89+
def test_output(self):
90+
mg_fe_perovskite = burnman.minerals.SLB_2011.mg_fe_perovskite()
91+
mg_fe_perovskite.set_composition([0.9, 0.1, 0]) # frac_mg, frac_fe, frac_al
92+
rock = burnman.Composite([mg_fe_perovskite], [1.0])
93+
depths = np.linspace(2890e3, 670e3, 20)
94+
lower_mantle = burnman.Layer(name="Pv LM", radii=6371.0e3 - depths)
95+
lower_mantle.set_material(rock)
96+
lower_mantle.set_temperature_mode(
97+
temperature_mode="adiabatic", temperature_top=1900.0
98+
)
99+
lower_mantle.set_pressure_mode(
100+
pressure_mode="self-consistent", pressure_top=2.4e10, gravity_bottom=10.0
101+
)
102+
lower_mantle.make()
103+
104+
# Make sure that all the writing functions return a value
105+
tvel = output_seismo.tvel_formatted_data_and_header(
106+
lower_mantle, background_model=burnman.seismic.PREM()
107+
)
108+
axisem = output_seismo.axisem_formatted_data_and_reference([lower_mantle])
109+
mineos = output_seismo.mineos_formatted_data_and_reference([lower_mantle])
110+
111+
# Write to nothing
112+
output_seismo.write_tvel_file(lower_mantle, [], burnman.seismic.PREM())
113+
output_seismo.write_axisem_input([lower_mantle], [], verbose=False)
114+
output_seismo.write_mineos_input([lower_mantle], [], verbose=False)
115+
87116

88117
if __name__ == "__main__":
89118
unittest.main()

‎tests/test_solidsolution.py‎

Lines changed: 46 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -511,6 +511,16 @@ def setup_ol_ss2(self):
511511
ol_ss.set_state(P, T)
512512
return fo, ol_ss
513513

514+
def test_ol_ss_formulae(self):
515+
P = 1.0e5
516+
T = 1000.0
517+
518+
ol_ss = olivine_ss()
519+
ol_ss.set_composition([0.4, 0.6])
520+
ol_ss.set_state(P, T)
521+
522+
self.assertEqual(ol_ss.site_formula(), "[Mg0.40Fe0.60]2SiO4")
523+
514524
def test_1_gibbs(self):
515525
fo, fo_ss = self.setup_1min_ss()
516526
with warnings.catch_warnings(record=True) as w:
@@ -1066,6 +1076,9 @@ def test_formula(self):
10661076
formula_to_string(sum_formulae(ol.endmember_formulae, [0.5, 0.5])),
10671077
"MgFeSiO4",
10681078
)
1079+
s = np.array([[2, 0, 1, 4], [0, 2, 1, 4]])
1080+
self.assertArraysAlmostEqual(ol.stoichiometric_array[0], s[0])
1081+
self.assertArraysAlmostEqual(ol.stoichiometric_array[1], s[1])
10691082

10701083
def test_stoichiometric_matrix_binary_solution(self):
10711084
olivine = burnman.minerals.SLB_2011.mg_fe_olivine()
@@ -1262,6 +1275,39 @@ def test_relaxed_solution(self):
12621275
# Give 5 J tolerance.
12631276
self.assertTrue(np.abs(mu) < 5.0)
12641277

1278+
V = ropx.V
1279+
ropx.set_state(1.0e10, T)
1280+
1281+
with warnings.catch_warnings(record=True) as w:
1282+
warnings.simplefilter("always")
1283+
ropx.set_state_with_volume(V * 1.000001, T)
1284+
self.assertTrue(ropx.pressure < 1.0e9)
1285+
ropx.set_state_with_volume(V * 0.999999, T)
1286+
self.assertTrue(ropx.pressure > 1.0e9)
1287+
1288+
self.assertAlmostEqual(
1289+
ropx.isothermal_compressibility_reuss,
1290+
1.0 / ropx.isothermal_bulk_modulus_reuss,
1291+
)
1292+
self.assertAlmostEqual(
1293+
ropx.isentropic_compressibility_reuss,
1294+
1.0 / ropx.isentropic_bulk_modulus_reuss,
1295+
)
1296+
self.assertTrue(
1297+
ropx.isothermal_compressibility_reuss
1298+
> ropx.isentropic_compressibility_reuss
1299+
)
1300+
self.assertTrue(ropx.molar_heat_capacity_p > ropx.molar_heat_capacity_v)
1301+
1302+
self.assertTrue(ropx.thermal_expansivity > ropx.unrelaxed.thermal_expansivity)
1303+
self.assertTrue(
1304+
ropx.molar_heat_capacity_p > ropx.unrelaxed.molar_heat_capacity_p
1305+
)
1306+
self.assertTrue(
1307+
ropx.isothermal_compressibility_reuss
1308+
> ropx.unrelaxed.isothermal_compressibility_reuss
1309+
)
1310+
12651311

12661312
if __name__ == "__main__":
12671313
unittest.main()

0 commit comments

Comments
 (0)