use rayon::iter::{IndexedParallelIterator, IntoParallelIterator, IntoParallelRefIterator};
use super::utils::check_topology_state_sizes;
use crate::prelude::*;
use std::path::Path;
#[derive(Debug, Default)]
pub struct System {
pub(super) top: Topology,
pub(super) st: State,
}
impl System {
pub fn new(top: Topology, st: State) -> Result<Self, SelectionError> {
check_topology_state_sizes(&top, &st)?;
Ok(Self { top, st })
}
pub fn from_file(fname: impl AsRef<Path>) -> Result<Self, SelectionError> {
let mut fh = FileHandler::open(fname)?;
let (top, st) = fh.read()?;
Ok(Self::new(top, st)?)
}
pub fn sub_select(&self, ind: &impl IndexSliceProvider, def: impl SelectionDef) -> Result<Sel, SelectionError> {
unsafe {Ok(self.try_bind_sorted_slice(ind.get_index_slice())?.select(def)?)}
}
pub fn select_all(&self) -> Sel {
Sel(unsafe { SVec::from_sorted((0..self.len()).into_iter().collect()) })
}
pub fn select_all_bound(&self) -> SelOwnBound<'_> {
SelOwnBound {
sys: self,
index: unsafe { SVec::from_sorted((0..self.len()).into_iter().collect()) },
}
}
pub fn select_bound_mut(
&mut self,
def: impl SelectionDef,
) -> Result<SelOwnBoundMut<'_>, SelectionError> {
let index = def.into_sel_index(self, None)?;
Ok(SelOwnBoundMut { sys: self, index })
}
pub fn select_all_bound_mut(&mut self) -> SelOwnBoundMut<'_> {
let index = unsafe { SVec::from_sorted((0..self.len()).into_iter().collect()) };
SelOwnBoundMut { sys: self, index }
}
unsafe fn try_bind_sorted_slice<'a>(&'a self, sel: &'a [usize]) -> Result<SelBound<'a>, SelectionError> {
let last = unsafe { *sel.get_unchecked(sel.len() - 1) };
if last >= self.top.len() {
Err(SelectionError::IndexValidation(
*sel.first().unwrap(),
last,
self.top.len() - 1,
))
} else {
Ok(SelBound {
sys: self,
index: sel,
})
}
}
pub fn bind<'a>(&'a self, sel: &'a Sel) -> SelBound<'a> {
let last = unsafe { *sel.0.get_unchecked(sel.0.len() - 1) };
if last >= self.top.len() {
panic!("selection is out of bounds");
} else {
SelBound {
sys: self,
index: sel.0.as_slice(),
}
}
}
pub fn try_bind<'a>(&'a self, sel: &'a Sel) -> Result<SelBound<'a>, SelectionError> {
let last = unsafe { *sel.0.get_unchecked(sel.0.len() - 1) };
if last >= self.top.len() {
Err(SelectionError::IndexValidation(
*sel.0.first().unwrap(),
last,
self.top.len() - 1,
))
} else {
Ok(SelBound {
sys: self,
index: sel.0.as_slice(),
})
}
}
pub fn bind_mut<'a>(&'a mut self, sel: &'a Sel) -> SelBoundMut<'a> {
let last = unsafe { *sel.0.get_unchecked(sel.0.len() - 1) };
if last >= self.top.len() {
panic!("selection is out of bounds");
} else {
SelBoundMut {
sys: self,
index: sel.0.as_slice(),
}
}
}
pub fn try_bind_mut<'a>(&'a mut self, sel: &'a Sel) -> Result<SelBoundMut<'a>, SelectionError> {
let last = unsafe { *sel.0.get_unchecked(sel.0.len() - 1) };
if last >= self.top.len() {
Err(SelectionError::IndexValidation(
*sel.0.first().unwrap(),
last,
self.top.len() - 1,
))
} else {
Ok(SelBoundMut {
sys: self,
index: sel.0.as_slice(),
})
}
}
pub fn iter_par_split_mut<'a>(
&'a mut self,
par: &'a ParSplit,
) -> impl IndexedParallelIterator<Item = SelParMut<'a>> {
par.check_bounds(self);
par.selections
.par_iter()
.map(|sel| SelParMut::new(self, &sel.0))
}
pub fn iter_par_split<'a>(
&'a mut self,
par: &'a ParSplit,
) -> impl IndexedParallelIterator<Item = SelPar<'a>> {
par.check_bounds(self);
par.selections
.par_iter()
.map(|sel| SelPar::new(self, &sel.0))
}
pub fn split_bound<'a, RT, F>(&'a self, func: F) -> impl Iterator<Item = SelOwnBound<'a>>
where
RT: Default + std::cmp::PartialEq + 'a,
F: Fn(Particle) -> Option<RT> + 'a,
{
self.split(func).map(|sel| SelOwnBound {
sys: self,
index: sel.0,
})
}
pub fn split_resindex_bound(&self) -> impl Iterator<Item = SelOwnBound<'_>> {
self.split_bound(|p| Some(p.atom.resindex))
}
pub fn set_state(&mut self, st: State) -> Result<State, SelectionError> {
if self.len() != st.len() {
Err(SelectionError::IncompatibleState)
} else {
Ok(std::mem::replace(&mut self.st, st))
}
}
pub fn set_topology(&mut self, top: Topology) -> Result<Topology, SelectionError> {
if self.len() != top.len() {
Err(SelectionError::IncompatibleTopology)
} else {
Ok(std::mem::replace(&mut self.top, top))
}
}
pub fn release(self) -> (Topology, State) {
(self.top, self.st)
}
pub fn append_from_self(&mut self, sel: &Sel) -> Result<Sel, SelectionError> {
let old_last = self.len() - 1;
let last_ind = unsafe { *sel.0.get_unchecked(sel.0.len() - 1) };
if last_ind >= self.top.len() {
return Err(SelectionError::IndexValidation(
sel.0[0],
last_ind,
self.top.len(),
))?;
}
let pos: Vec<_> = sel.0.iter().map(|i| &self.st.coords[*i]).cloned().collect();
let atoms: Vec<_> = sel.0.iter().map(|i| &self.top.atoms[*i]).cloned().collect();
self.st.add_coords(pos.into_iter());
self.top.add_atoms(atoms.into_iter());
Ok(Sel::from_iter(old_last + 1..self.len())?)
}
pub fn append_atoms_coords<'a>(
&mut self,
atoms: impl AtomIterator<'a>,
coords: impl PosIterator<'a>,
) -> Result<Sel, SelectionError> {
let old_last = self.len() - 1;
self.st.add_coords(coords.cloned());
self.top.add_atoms(atoms.cloned());
check_topology_state_sizes(&self.top, &self.st)?;
Ok(Sel::from_iter(old_last + 1..self.len())?)
}
pub fn append(
&mut self,
data: &(impl AtomIterProvider + PosIterProvider),
) -> Result<Sel, SelectionError> {
let old_last = self.len() - 1;
self.st.add_coords(data.iter_pos().cloned());
self.top.add_atoms(data.iter_atoms().cloned());
check_topology_state_sizes(&self.top, &self.st)?;
Ok(Sel::from_iter(old_last + 1..self.len())?)
}
pub fn remove(
&mut self,
removed: impl Iterator<Item = usize> + Clone,
) -> Result<(), BuilderError> {
self.st.remove_coords(removed.clone())?;
self.top.remove_atoms(removed)?;
Ok(())
}
pub fn set_box_from(&mut self, src: &impl BoxProvider) {
self.st.pbox = src.get_box().cloned();
}
pub fn multiply_periodically(&mut self, nbox: [usize; 3]) -> Result<(), SelectionError> {
if self.get_box().is_none() {
return Err(PeriodicBoxError::NoPbc)?;
}
let m = self.require_box()?.get_matrix();
let all = self.select_all();
for x in 0..=nbox[0] {
for y in 0..=nbox[1] {
for z in 0..=nbox[2] {
if x == 0 && y == 0 && z == 0 {
continue;
}
let added = self.append_from_self(&all)?;
let shift =
m.column(0) * x as f32 + m.column(1) * y as f32 + m.column(2) * z as f32;
self.select_bound_mut(added)?.translate(&shift);
}
}
}
self.get_box_mut().unwrap().scale_vectors([
nbox[0] as f32,
nbox[1] as f32,
nbox[2] as f32,
])?;
self.top.assign_resindex();
Ok(())
}
pub fn assign_resindex(&mut self) {
self.top.assign_resindex();
}
}
impl Selectable for System {
fn select(&self, def: impl SelectionDef) -> Result<Sel, SelectionError> {
Ok(Sel(def.into_sel_index(self, None)?))
}
}
impl SelectableBound for System {
fn select_bound(&self, def: impl SelectionDef) -> Result<SelOwnBound<'_>, SelectionError> {
Ok(SelOwnBound {
index: def.into_sel_index(self, None)?,
sys: unsafe{&*self.get_system_ptr()},
})
}
}
impl SaveTopology for System {
fn iter_atoms_dyn(&self) -> Box<dyn Iterator<Item = &Atom> + '_> {
Box::new(self.iter_atoms())
}
}
impl SaveState for System {
fn iter_pos_dyn<'a>(&'a self) -> Box<dyn Iterator<Item = &'a Pos> + 'a> {
Box::new(self.iter_pos())
}
}
impl SaveTopologyState for System {}
impl SystemProvider for System {
fn get_system_ptr(&self) -> *const System {
self
}
}
impl SystemMutProvider for System {}
impl LenProvider for System {
fn len(&self) -> usize {
self.top.atoms.len()
}
}
impl IndexProvider for System {
unsafe fn get_index_unchecked(&self, i: usize) -> usize {
i
}
fn iter_index(&self) -> impl Iterator<Item = usize> {
(0..self.len()).into_iter()
}
}
impl IndexParProvider for System {
fn par_iter_index(&self) -> impl IndexedParallelIterator<Item = usize> {
(0..self.len()).into_par_iter()
}
}
impl<'a> std::ops::Shr<&'a System> for &'a Sel {
type Output = SelBound<'a>;
fn shr(self, rhs: &'a System) -> Self::Output {
rhs.bind(self)
}
}
impl<'a> std::ops::Shr<&'a mut System> for &'a Sel {
type Output = SelBoundMut<'a>;
fn shr(self, rhs: &'a mut System) -> Self::Output {
rhs.bind_mut(self)
}
}