[cig-commits] [commit] inversion: Added a consistency inversion. Given an frac_pv and Fe fractions for each, see if you get the same answer out. (fcc314f)

[cig_noreply at geodynamics.org]

Repository : https://github.com/geodynamics/burnman

On branch  : inversion
Link       : https://github.com/geodynamics/burnman/compare/776551fd2a2dab635fe5b4124b8d3675d3a1b82a...fcc314ffd83374fe1d6a4d79f8bd2f8be734dcd1

>---------------------------------------------------------------

commit fcc314ffd83374fe1d6a4d79f8bd2f8be734dcd1
Author: Cayman Unterborn <kmanunterborn at gmail.com>
Date:   Sat Dec 13 17:11:11 2014 -0800

    Added a consistency inversion. Given an frac_pv and Fe fractions for each, see if you get the same answer out.


>---------------------------------------------------------------

fcc314ffd83374fe1d6a4d79f8bd2f8be734dcd1
 ...inversion_setup.py => knowndata_consistency.py} | 22 +++++++++++++---------
 1 file changed, 13 insertions(+), 9 deletions(-)

diff --git a/inversion/basic_inversion_setup.py b/inversion/knowndata_consistency.py
similarity index 93%
copy from inversion/basic_inversion_setup.py
copy to inversion/knowndata_consistency.py
index cfe8c4f..32fffb6 100644
--- a/inversion/basic_inversion_setup.py
+++ b/inversion/knowndata_consistency.py
@@ -25,16 +25,20 @@ number_of_points = 5 #set on how many depth slices the computations should be do
 
 
 # velocity constraints from seismology
+#seismic_model = burnman.seismic.PREM() # pick from .prem() .slow() .fast() (see code/seismic.py)
+fraction_pv = 0.5
+fe_pv = 0.3
+fe_pc = 0.3
+pv=minerals.other.mg_fe_perovskite(fe_pv)
+pc=minerals.other.ferropericlase(fe_pc)
+rock = burnman.Composite( [fraction_pv, 1.0-fraction_pv],[ pv,pc ] )
+
 seismic_model = burnman.seismic.PREM() # pick from .prem() .slow() .fast() (see code/seismic.py)
 depths = np.linspace(1000e3,2500e3, number_of_points)
-seis_p, seis_rho, seis_vp, seis_vs, seis_vphi = seismic_model.evaluate_all_at(depths)
-
-seis_G= seis_vs**2.*seis_rho
-seis_K= seis_vphi**2*seis_rho
-
-
-
+seis_p = seismic_model.pressure(depths)
 temperature = burnman.geotherm.brown_shankland(seis_p)
+seis_rho, seis_vp, seis_vs, seis_vphi, seis_K,seis_G = burnman.main.velocities_from_rock(rock,seis_p,temperature,burnman.averaging_schemes.VoigtReussHill())
+
 '''
 
 #velocities from known material (still named seis*)
@@ -147,14 +151,14 @@ if whattodo=="continue":
         print "*** run=%d/%d, # samples: %d" % (l, n_runs, db.trace('fraction_pv').stats()['n'] )
         S = pymc.MCMC(model, db=db)
         #S.sample(iter=100, burn=10, thin=1)
-        S.sample(iter=1000, burn=0, thin=10) # Search space for 100000 acceptable steps, forget first 1000 and save every 10.
+        S.sample(iter=100, burn=0, thin=10) # Search space for 100000 acceptable steps, forget first 1000 and save every 10.
         S.db.close()
 
 if whattodo=="plot":
 	files=sys.argv[2:]
 	print "files:",files
         
-        b=100#10000 # burn number
+        b=10#10000 # burn number
         i=1
 
         for t in things: