parent
a244c49256
commit
f436222852
1 changed files with 119 additions and 0 deletions
@ -0,0 +1,119 @@ |
||||
# Created: 20150528 |
||||
# Test module for stochastic fitting |
||||
|
||||
import numpy |
||||
from wpylib.params import struct |
||||
from wpylib.math.fitting.stochastic import StochasticFitting |
||||
|
||||
def setup_MC_TZ(): |
||||
""" |
||||
Internal identifier: Cr2_TZ_data_20140728uhf |
||||
Binding energy of Cr2, phaseless QMC/UHF, cc-pwCVTZ-DK, dt=0.01 |
||||
""" |
||||
global Cr2_TZ_data_20140728uhf |
||||
from numpy import asarray, matrix |
||||
Cr2_TZ_data_20140728uhf = asarray(matrix(""" |
||||
1.550 -1.613 0.025 ; |
||||
1.600 -1.897 0.036 ; |
||||
1.6788 -2.141 0.016 ; |
||||
1.720 -2.143 0.025 ; |
||||
1.800 -2.190 0.022 ; |
||||
1.900 -2.064 0.020 ; |
||||
2.000 -2.038 0.023 ; |
||||
2.400 -1.611 0.019 ; |
||||
3.000 -1.037 0.018 |
||||
""")) |
||||
|
||||
|
||||
|
||||
def test_fit_PEC_MC_TZ(show_samples=True, save_fig=False): |
||||
"""20150528 |
||||
PEC fitting, using `morse2` functional form. |
||||
|
||||
Options: |
||||
- show_samples : prints the sample data (recommended) |
||||
- save_fig : dumps the fit result for each stochastic snapshot |
||||
(not recommended unless you are debugging/diagnosing) |
||||
|
||||
Modeled after do_morse2_sfit routine in Cr2 analysis script. |
||||
|
||||
Example output: |
||||
|
||||
test_MC_TZ_PEC_fit:: |
||||
ansatz=morse2_fit_func |
||||
Raw data: (x, y, dy) |
||||
1.550 -1.613 0.025 |
||||
1.600 -1.897 0.036 |
||||
1.679 -2.141 0.016 |
||||
1.720 -2.143 0.025 |
||||
1.800 -2.190 0.022 |
||||
1.900 -2.064 0.020 |
||||
2.000 -2.038 0.023 |
||||
2.400 -1.611 0.019 |
||||
3.000 -1.037 0.018 |
||||
Using guess param from NLF: [-2.18625505 9.82497657 1.80455072 1.86192922] |
||||
Final parameters : -2.187(10) 9.85(61) 1.8044(66) 1.864(83) |
||||
Total execution time = 0.55 secs |
||||
All testings passed. |
||||
""" |
||||
from numpy import array |
||||
from wpylib.text_tools import matrix_str, str_indent |
||||
from wpylib.timer import block_timer as Timer |
||||
from wpylib.math.fitting.funcs_pec import morse2_fit_func |
||||
from wpylib.py import function_name |
||||
|
||||
global MC_TZ_PEC_fit |
||||
global Cr2_TZ_data_20140728uhf |
||||
|
||||
print("test_MC_TZ_PEC_fit::") |
||||
setup_MC_TZ() |
||||
|
||||
# parameters etc |
||||
ansatz = morse2_fit_func() |
||||
rng = dict(seed=378711, rng_class=numpy.random.RandomState) |
||||
rawdata = Cr2_TZ_data_20140728uhf |
||||
sfit = MC_TZ_PEC_fit = StochasticFitting() |
||||
|
||||
print("ansatz=%s" \ |
||||
% (function_name(ansatz),)) |
||||
|
||||
# This corresponds to fit_variant==1 in the subroutine |
||||
ansatz.fit_method = "leastsq" |
||||
ansatz.fit_opts['leastsq'] = dict(xtol=1e-8, epsfcn=1e-6) |
||||
|
||||
sfit.init_func(ansatz) |
||||
sfit.init_samples(x=rawdata[:,0], y=rawdata[:,1], dy=rawdata[:,2]) |
||||
sfit.init_rng(**rng) |
||||
if show_samples: |
||||
print("Raw data: (x, y, dy)") |
||||
print(str_indent(matrix_str(rawdata, fmt="%8.3f"))) |
||||
|
||||
with Timer() as tmr: |
||||
sfit.mcfit_loop_begin_() |
||||
sfit.mcfit_loop1_(num_iter=200, save_fig=save_fig) |
||||
sfit.mcfit_loop_end_() |
||||
sfit.mcfit_analysis_() |
||||
sfit.mcfit_report_final_params() |
||||
|
||||
nlf_fit_result = sfit.nlf_rec['xopt'] |
||||
|
||||
# Verify the results: |
||||
# * deterministic fit (only precise to these digits; more digits are discarded) |
||||
nlf_fit_ref = numpy.array([-2.18625505, 9.82497657, 1.80455072, 1.86192922]) |
||||
|
||||
assert numpy.allclose(nlf_fit_result, nlf_fit_ref, atol=0.6e-8, rtol=0) |
||||
|
||||
def match_mc_param(name, val_ref, err_ref, atol): |
||||
vv = sfit.final_mc_params[name] |
||||
assert numpy.allclose(vv.value(), val_ref, atol=atol, rtol=0) |
||||
assert numpy.allclose(vv.error(), err_ref, atol=atol, rtol=0) |
||||
|
||||
match_mc_param('E0', -2.187, 0.010, 0.6e-3) |
||||
match_mc_param('k', 9.85, 0.61, 0.6e-2) |
||||
match_mc_param('r0', 1.8044, 0.0066, 0.6e-4) |
||||
match_mc_param('a', 1.864, 0.083, 0.6e-3) |
||||
print("All testings passed.") |
||||
|
||||
|
||||
if __name__ == '__main__': |
||||
test_fit_PEC_MC_TZ() |
Loading…
Reference in new issue