mirror of
https://github.com/dnlbauer/WHAM.git
synced 2026-09-10 14:15:31 +00:00
clean code
This commit is contained in:
7
Cargo.lock
generated
7
Cargo.lock
generated
@@ -20,6 +20,12 @@ dependencies = [
|
|||||||
"winapi",
|
"winapi",
|
||||||
]
|
]
|
||||||
|
|
||||||
|
[[package]]
|
||||||
|
name = "assert_approx_eq"
|
||||||
|
version = "1.1.0"
|
||||||
|
source = "registry+https://github.com/rust-lang/crates.io-index"
|
||||||
|
checksum = "3c07dab4369547dbe5114677b33fbbf724971019f3818172d59a97a61c774ffd"
|
||||||
|
|
||||||
[[package]]
|
[[package]]
|
||||||
name = "atty"
|
name = "atty"
|
||||||
version = "0.2.14"
|
version = "0.2.14"
|
||||||
@@ -355,6 +361,7 @@ name = "wham"
|
|||||||
version = "0.9.6"
|
version = "0.9.6"
|
||||||
dependencies = [
|
dependencies = [
|
||||||
"GSL",
|
"GSL",
|
||||||
|
"assert_approx_eq",
|
||||||
"clap",
|
"clap",
|
||||||
"error-chain",
|
"error-chain",
|
||||||
"rand",
|
"rand",
|
||||||
|
|||||||
@@ -13,6 +13,9 @@ rand = "0.5.5"
|
|||||||
GSL = "1.0.0"
|
GSL = "1.0.0"
|
||||||
rayon = "1.0.3"
|
rayon = "1.0.3"
|
||||||
|
|
||||||
|
[dev-dependencies]
|
||||||
|
assert_approx_eq = "1.1.0"
|
||||||
|
|
||||||
[profile.release]
|
[profile.release]
|
||||||
opt-level = 2
|
opt-level = 2
|
||||||
|
|
||||||
|
|||||||
@@ -239,10 +239,10 @@ mod tests {
|
|||||||
fn get_x_for_bin() {
|
fn get_x_for_bin() {
|
||||||
let ds = build_hist_set();
|
let ds = build_hist_set();
|
||||||
let expected: Vec<f64> = vec![0,1,2,3,4,5,6,7,8].iter()
|
let expected: Vec<f64> = vec![0,1,2,3,4,5,6,7,8].iter()
|
||||||
.map(|x| *x as f64 + 0.5).collect();
|
.map(|x| *x as f64 + 0.5).collect();
|
||||||
for i in 0..9 {
|
expected.iter().enumerate().for_each(|(i, exp)| {
|
||||||
assert_eq!(expected[i], ds.get_coords_for_bin(i)[0]);
|
assert_approx_eq!(exp, &ds.get_coords_for_bin(i)[0]);
|
||||||
}
|
})
|
||||||
}
|
}
|
||||||
|
|
||||||
#[test]
|
#[test]
|
||||||
|
|||||||
27
src/io.rs
27
src/io.rs
@@ -35,7 +35,7 @@ pub fn read_data(cfg: &Config) -> Result<Dataset> {
|
|||||||
let bin_width: Vec<f64> = (0..cfg.dimens).map(|idx| {
|
let bin_width: Vec<f64> = (0..cfg.dimens).map(|idx| {
|
||||||
(cfg.hist_max[idx] - cfg.hist_min[idx])/(cfg.num_bins[idx] as f64)
|
(cfg.hist_max[idx] - cfg.hist_min[idx])/(cfg.num_bins[idx] as f64)
|
||||||
}).collect();
|
}).collect();
|
||||||
let num_bins = cfg.num_bins.iter().fold(1, |state, &bins| state*bins);
|
let num_bins = cfg.num_bins.iter().product();
|
||||||
let dimens_length = cfg.num_bins.clone();
|
let dimens_length = cfg.num_bins.clone();
|
||||||
|
|
||||||
let f = File::open(&cfg.metadata_file).chain_err(|| "Failed to open metadata file")?;
|
let f = File::open(&cfg.metadata_file).chain_err(|| "Failed to open metadata file")?;
|
||||||
@@ -118,7 +118,7 @@ fn is_in_time_boundaries(time: f64, cfg: &Config) -> bool {
|
|||||||
false
|
false
|
||||||
}
|
}
|
||||||
|
|
||||||
// parse a timeseries file into a histogram
|
// parse a time series file into a histogram
|
||||||
fn read_window_file(window_file: &str, cfg: &Config) -> Result<Histogram> {
|
fn read_window_file(window_file: &str, cfg: &Config) -> Result<Histogram> {
|
||||||
let f = File::open(window_file)
|
let f = File::open(window_file)
|
||||||
.chain_err(|| format!("Failed to open sample data file {}.", window_file))?;
|
.chain_err(|| format!("Failed to open sample data file {}.", window_file))?;
|
||||||
@@ -139,7 +139,7 @@ fn read_window_file(window_file: &str, cfg: &Config) -> Result<Histogram> {
|
|||||||
while buf.read_line(&mut line).chain_err(|| "Failed to read line")? > 0 {
|
while buf.read_line(&mut line).chain_err(|| "Failed to read line")? > 0 {
|
||||||
linecount += 1;
|
linecount += 1;
|
||||||
// skip comments and empty lines
|
// skip comments and empty lines
|
||||||
if line.starts_with("#") || line.starts_with("@") || line.len() == 0 {
|
if line.starts_with('#') || line.starts_with('@') || line.is_empty() {
|
||||||
line.clear();
|
line.clear();
|
||||||
continue;
|
continue;
|
||||||
}
|
}
|
||||||
@@ -182,14 +182,14 @@ pub fn write_results(out_file: &str, ds: &Dataset, free: &[f64],
|
|||||||
|
|
||||||
let header: String = (0..ds.dimens_lengths.len()).map(|d| format!("coord{}", d+1))
|
let header: String = (0..ds.dimens_lengths.len()).map(|d| format!("coord{}", d+1))
|
||||||
.collect::<Vec<String>>().join(" ");
|
.collect::<Vec<String>>().join(" ");
|
||||||
writeln!(buf, "#{} {} {} {} {}",
|
writeln!(buf, "#{} Free Energy +/- Probability +/-", header).unwrap();
|
||||||
header, "Free Energy", "+/-", "Probability", "+/-").unwrap();
|
|
||||||
|
|
||||||
for bin in 0..free.len() {
|
for bin in 0..free.len() {
|
||||||
let coords = ds.get_coords_for_bin(bin);
|
let coords = ds.get_coords_for_bin(bin);
|
||||||
let coords_str: String = coords.iter().map(|c| {format!("{:8.6} ", c)})
|
let coords_str: String = coords.iter().map(|c| {format!("{:8.6} ", c)})
|
||||||
.collect::<Vec<String>>().join("\t");
|
.collect::<Vec<String>>().join("\t");
|
||||||
writeln!(buf, "{}{:8.6} {:8.6} {:8.6} {:8.6}", coords_str, free[bin], free_std[bin], prob[bin], prob_std[bin])
|
writeln!(buf, "{}{:8.6} {:8.6} {:8.6} {:8.6}", coords_str,
|
||||||
|
free[bin], free_std[bin], prob[bin], prob_std[bin])
|
||||||
.chain_err(|| "Failed to write to file.")?;
|
.chain_err(|| "Failed to write to file.")?;
|
||||||
}
|
}
|
||||||
Ok(())
|
Ok(())
|
||||||
@@ -198,6 +198,7 @@ pub fn write_results(out_file: &str, ds: &Dataset, free: &[f64],
|
|||||||
#[cfg(test)]
|
#[cfg(test)]
|
||||||
mod tests {
|
mod tests {
|
||||||
use super::*;
|
use super::*;
|
||||||
|
use assert_approx_eq::assert_approx_eq;
|
||||||
|
|
||||||
fn cfg() -> Config {
|
fn cfg() -> Config {
|
||||||
Config {
|
Config {
|
||||||
@@ -227,12 +228,12 @@ mod tests {
|
|||||||
let h = super::read_window_file(&f, &cfg).unwrap();
|
let h = super::read_window_file(&f, &cfg).unwrap();
|
||||||
println!("{:?}", h);
|
println!("{:?}", h);
|
||||||
assert_eq!(5000, h.num_points);
|
assert_eq!(5000, h.num_points);
|
||||||
assert_eq!(0.0, h.bins[2]);
|
assert_approx_eq!(0.0, h.bins[2]);
|
||||||
assert_eq!(11.0, h.bins[3]);
|
assert_approx_eq!(11.0, h.bins[3]);
|
||||||
assert_eq!(2236.0, h.bins[4]);
|
assert_approx_eq!(2236.0, h.bins[4]);
|
||||||
assert_eq!(2714.0, h.bins[5]);
|
assert_approx_eq!(2714.0, h.bins[5]);
|
||||||
assert_eq!(39.0, h.bins[6]);
|
assert_approx_eq!(39.0, h.bins[6]);
|
||||||
assert_eq!(0.0, h.bins[7]);
|
assert_approx_eq!(0.0, h.bins[7]);
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|
||||||
@@ -244,7 +245,7 @@ mod tests {
|
|||||||
assert_eq!(25, ds.num_windows);
|
assert_eq!(25, ds.num_windows);
|
||||||
assert_eq!(cfg.num_bins.len(), ds.dimens_lengths.len());
|
assert_eq!(cfg.num_bins.len(), ds.dimens_lengths.len());
|
||||||
assert_eq!(cfg.num_bins[0], ds.dimens_lengths[0]);
|
assert_eq!(cfg.num_bins[0], ds.dimens_lengths[0]);
|
||||||
assert_eq!(cfg.temperature * k_B, ds.kT);
|
assert_approx_eq!(cfg.temperature * k_B, ds.kT);
|
||||||
assert_eq!(25, ds.histograms.len())
|
assert_eq!(25, ds.histograms.len())
|
||||||
}
|
}
|
||||||
|
|
||||||
|
|||||||
37
src/lib.rs
37
src/lib.rs
@@ -5,6 +5,10 @@ extern crate error_chain;
|
|||||||
extern crate rand;
|
extern crate rand;
|
||||||
extern crate rgsl;
|
extern crate rgsl;
|
||||||
extern crate rayon;
|
extern crate rayon;
|
||||||
|
#[cfg(test)]
|
||||||
|
#[macro_use]
|
||||||
|
extern crate assert_approx_eq;
|
||||||
|
|
||||||
|
|
||||||
pub mod io;
|
pub mod io;
|
||||||
pub mod histogram;
|
pub mod histogram;
|
||||||
@@ -21,7 +25,7 @@ pub mod errors { error_chain!{} }
|
|||||||
use errors::*;
|
use errors::*;
|
||||||
|
|
||||||
#[allow(non_upper_case_globals)]
|
#[allow(non_upper_case_globals)]
|
||||||
static k_B: f64 = 0.0083144621; // kJ/mol*K
|
static k_B: f64 = 0.008_314_462_1; // kJ/mol*K
|
||||||
|
|
||||||
// Application config
|
// Application config
|
||||||
#[derive(Debug)]
|
#[derive(Debug)]
|
||||||
@@ -180,16 +184,11 @@ pub fn run(cfg: &Config) -> Result<()>{
|
|||||||
|
|
||||||
let (P, F, F_prev) = perform_wham(&cfg, &dataset)?;
|
let (P, F, F_prev) = perform_wham(&cfg, &dataset)?;
|
||||||
|
|
||||||
let P_std: Vec<f64>;
|
let (P_std, free_energy_std) = if cfg.bootstrap > 0 {
|
||||||
let free_energy_std: Vec<f64>;
|
error_analysis::run_bootstrap(&cfg, dataset.clone(), &P, cfg.bootstrap)
|
||||||
if cfg.bootstrap > 0 {
|
|
||||||
let error_est = error_analysis::run_bootstrap(&cfg, dataset.clone(), &P, cfg.bootstrap);
|
|
||||||
P_std = error_est.0;
|
|
||||||
free_energy_std = error_est.1;
|
|
||||||
} else {
|
} else {
|
||||||
P_std = vec![0.0; P.len()];
|
(vec![0.0; P.len()], vec![0.0; P.len()])
|
||||||
free_energy_std = vec![0.0; P.len()];
|
};
|
||||||
}
|
|
||||||
|
|
||||||
// calculate free energy and dump state
|
// calculate free energy and dump state
|
||||||
println!("Finished. Dumping final PMF");
|
println!("Finished. Dumping final PMF");
|
||||||
@@ -292,21 +291,21 @@ mod tests {
|
|||||||
let F = vec![1.0; dataset.num_bins] ;
|
let F = vec![1.0; dataset.num_bins] ;
|
||||||
let expected = vec!(0.0, 0.082_529_668_703_131_6, 40.923_558_470_974_93,
|
let expected = vec!(0.0, 0.082_529_668_703_131_6, 40.923_558_470_974_93,
|
||||||
124_226.700_033_77, 2_308_526_035.528_374_7);
|
124_226.700_033_77, 2_308_526_035.528_374_7);
|
||||||
for b in 0..dataset.num_bins {
|
expected.iter().enumerate().for_each(|(i, exp)| {
|
||||||
let p = super::calc_bin_probability(b, &dataset, &F);
|
let p = super::calc_bin_probability(i, &dataset, &F);
|
||||||
assert_delta!(expected[b], p, 0.000_000_1);
|
assert_delta!(exp, p, 0.000_000_1);
|
||||||
}
|
})
|
||||||
}
|
}
|
||||||
|
|
||||||
#[test]
|
#[test]
|
||||||
fn calc_bias_offset() {
|
fn calc_bias_offset() {
|
||||||
let dataset = create_test_dataset();
|
let dataset = create_test_dataset();
|
||||||
let probability = vec!(0.0, 0.1, 0.2, 0.3, 0.4);
|
let probability = vec!(0.0, 0.1, 0.2, 0.3, 0.4);
|
||||||
let expected = vec!(15.927_477_169_990_633, 15.927_477_169_990_633);
|
let expected = vec!(15.927_477_169_990_633, 15.927_477_169_990_633);
|
||||||
for window in 0..dataset.num_windows {
|
expected.iter().enumerate().for_each(|(i, exp)| {
|
||||||
let F = super::calc_window_F(window, &dataset, &probability);
|
let F = super::calc_window_F(i, &dataset, &probability);
|
||||||
assert_delta!(expected[window], F, 0.000_000_1);
|
assert_delta!(exp, F, 0.000_000_1);
|
||||||
}
|
})
|
||||||
}
|
}
|
||||||
|
|
||||||
#[test]
|
#[test]
|
||||||
|
|||||||
Reference in New Issue
Block a user