Pteros  2.0
Molecular modeling library for human beings!
system.h
1 /*
2  *
3  * This source code is part of
4  * ******************
5  * *** Pteros ***
6  * ******************
7  * molecular modeling library
8  *
9  * Copyright (c) 2009-2013, Semen Yesylevskyy
10  *
11  * This program is free software; you can redistribute it and/or
12  * modify it under the terms of Artistic License:
13  *
14  * Please note, that Artistic License is slightly more restrictive
15  * then GPL license in terms of distributing the modified versions
16  * of this software (they should be approved first).
17  * Read http://www.opensource.org/licenses/artistic-license-2.0.php
18  * for details. Such license fits scientific software better then
19  * GPL because it prevents the distribution of bugged derivatives.
20  *
21 */
22 
23 #ifndef SYSTEM_H
24 #define SYSTEM_H
25 
26 #include <string>
27 #include <vector>
28 #include <functional>
29 #include <memory>
30 #include <Eigen/Core>
31 #include <Eigen/Dense>
32 #include "pteros/core/atom.h"
33 #include "pteros/core/force_field.h"
34 #include "pteros/core/periodic_box.h"
35 #include "pteros/core/typedefs.h"
36 
37 namespace pteros {
38 
42  float total;
44  float lj_14;
46  float q_14;
48  float lj_sr;
50  float q_sr;
51 
52  Energy_components(): total(0.0), lj_14(0.0), q_14(0.0), lj_sr(0.0), q_sr(0.0) {}
53 
56  std::string to_str();
57 
60  Energy_components& operator+=(const Energy_components& other);
61 };
62 
67 struct Frame {
69  std::vector<Eigen::Vector3f> coord;
73  float time;
74 
75  Frame(): time(0.0) {}
76 };
77 
78 //Forward declarations
79 class Selection;
80 class Atom_proxy;
81 class Mol_file;
82 typedef unsigned short Mol_file_content;
83 
95 class System {
96  // System and Selection are friends because they are closely integrated.
97  friend class Selection;
98  // Selection_parser must access internals of Selection
99  friend class Selection_parser;
100  // Mol_file needs an access too
101  friend class Mol_file;
102 public:
103  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
106 
108  System();
109 
111  System(std::string fname);
112 
114  System(const System& other);
115 
117  System& operator=(System other);
118 
120  virtual ~System();
121 
123 
127 
130  Selection append(const System& sys);
131 
135  Selection append(const Selection& sel);
136 
139  Selection append(const Atom& at, Vector3f_const_ref coord);
140 
152  Selection append(const Atom_proxy& at);
153 
157  void rearrange(const std::vector<std::string>& sel_strings);
158 
162  void rearrange(const std::vector<Selection>& sel_vec);
163 
165  void keep(const std::string& sel_str);
166 
168  void keep(const Selection& sel);
169 
171  void remove(const std::string& sel_str);
172 
175  void remove(Selection& sel);
176 
180  void distribute(const Selection sel, Vector3i_const_ref ncopies, Matrix3f_const_ref shift);
181 
183 
184 
185  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
188 
190  inline int num_frames() const { return traj.size(); }
191 
193  inline int num_atoms() const { return atoms.size(); }
195 
196  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
216 
218  Selection select(std::string str, int fr = 0);
219  Selection operator()(std::string str, int fr = 0);
220 
221  Selection select(int ind1, int ind2);
222  Selection operator()(int ind1, int ind2);
223 
224  Selection select(const std::vector<int>& ind);
225  Selection operator()(const std::vector<int>& ind);
226 
227  Selection select(std::vector<int>::iterator it1,
228  std::vector<int>::iterator it2);
229  Selection operator()(std::vector<int>::iterator it1,
230  std::vector<int>::iterator it2);
231 
232  Selection select(const std::function<void(const System&,int,std::vector<int>&)>& callback, int fr = 0);
233  Selection operator()(const std::function<void(const System&,int,std::vector<int>&)>& callback, int fr = 0);
235 
238  Selection select_all();
239  Selection operator()();
241 
242  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
245 
254  // Skip functionality suggested by Raul Mera
255  void load(std::string fname,
256  int b=0,
257  int e=-1,
258  int skip = 0,
259  std::function<bool(System*,int)> on_frame = 0);
260 
268  bool load(const std::unique_ptr<Mol_file> &handler,
269  Mol_file_content what,
270  std::function<bool(System*,int)> on_frame = 0);
271 
272 
274 
299  void set_filter(std::string str);
301  void set_filter(int ind1, int ind2);
302  void set_filter(const std::vector<int>& ind);
303  void set_filter(std::vector<int>::iterator it1,
304  std::vector<int>::iterator it2);
306 
308 
309 
310  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
313 
315  int frame_dup(int fr);
316 
318  void frame_append(const Frame& fr);
319 
321  void frame_copy(int fr1, int fr2);
322 
327  void frame_delete(int b = 0, int e = -1);
328 
330  void frame_swap(int fr1, int fr2);
332 
333 
334  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
337 
339  inline Periodic_box& Box(int fr){
340  return traj[fr].box;
341  }
342 
344  inline const Periodic_box& Box(int fr) const {
345  return traj[fr].box;
346  }
347 
349  inline float& Time(int fr){
350  return traj[fr].time;
351  }
352 
354  inline const float& Time(int fr) const {
355  return traj[fr].time;
356  }
357 
359  inline Eigen::Vector3f& XYZ(int ind, int fr){
360  return traj[fr].coord[ind];
361  }
362 
364  inline const Eigen::Vector3f& XYZ(int ind, int fr) const {
365  return traj[fr].coord[ind];
366  }
367 
369  inline Atom& Atom_data(int ind) {
370  return atoms[ind];
371  }
372 
374  inline const Atom& Atom_data(int ind) const {
375  return atoms[ind];
376  }
377 
379  inline Frame& Frame_data(int fr){
380  return traj[fr];
381  }
382 
384  inline const Frame& Frame_data(int fr) const {
385  return traj[fr];
386  }
387 
389 
390  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
393 
395  void dssp(std::string fname) const;
396 
398  void dssp(std::ostream& os) const;
399 
413  std::string dssp() const;
415 
416 
417  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
421 
423  Selection atoms_dup(const std::vector<int>& ind);
424 
426  Selection atoms_add(const std::vector<Atom>& atm,
427  const std::vector<Eigen::Vector3f>& crd);
428 
430  void atoms_delete(const std::vector<int>& ind);
432 
433 
434  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
437 
439  void wrap(int fr, Vector3i_const_ref dims_to_wrap = Eigen::Vector3i::Ones());
440 
442 
443  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
446 
448  float distance(int i, int j, int fr, bool is_periodic = true,
449  Vector3i_const_ref dims = Eigen::Vector3i::Ones()) const;
450 
452  float angle(int i, int j, int k, int fr, bool is_periodic = true,
453  Vector3i_const_ref dims = Eigen::Vector3i::Ones()) const;
454 
456  float dihedral(int i, int j, int k, int l, int fr, bool is_periodic = true,
457  Vector3i_const_ref dims = Eigen::Vector3i::Ones()) const;
458 
460 
461  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
464 
466  Energy_components non_bond_energy(int at1, int at2, int frame, bool is_periodic = true) const;
467 
469  Energy_components non_bond_energy(const std::vector<Eigen::Vector2i>& nlist, int fr, bool is_periodic=true) const;
470 
472 
473 
474  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
477 
479  void clear();
480 
483  return force_field.ready;
484  }
485 
489  if(force_field.ready)
490  return &force_field;
491  else
492  return nullptr;
493  }
494 
497  void assign_resindex(int start=0);
498 
501  void sort_by_resindex();
502 
504 
505 protected:
506 
507  // Holds all atom attributes except the coordinates
508  std::vector<Atom> atoms;
509 
510  // Coordinates for any number of frames
511  std::vector<Frame> traj;
512 
513  // Force field parameters
514  Force_field force_field;
515 
516  // Supplementary function to check if last added frame contains same number
517  // of atoms as topology
518  void check_num_atoms_in_last_frame();
519 
520  // Indexes for filtering
521  std::vector<int> filter;
522  // Filter selection text for text-based filters
523  std::string filter_text;
524 
525  void filter_atoms();
526  void filter_coord(int fr);
527 };
528 
529 }
530 #endif /* SYSTEM_H */
std::string to_str()
Writes all energy components to string in the following order: total, lj_sr, lj_14, q_sr, q_14.
Definition: system.cpp:855
Selection parser class.
Definition: selection_parser.h:145
Periodic_box & Box(int fr)
Read/write access for periodic box for given frame.
Definition: system.h:339
const Periodic_box & Box(int fr) const
Read only access for periodic box for given frame.
Definition: system.h:344
float total
Total energy.
Definition: system.h:42
Definition of single trajectory frame.
Definition: system.h:67
float lj_14
Lenard-Jones energy of 1-4 pairs.
Definition: system.h:44
Pteros namespace.
Definition: options.cpp:32
std::vector< Eigen::Vector3f > coord
Coordinates of atoms.
Definition: system.h:69
bool force_field_ready()
Returns true if the force field is set up properly and is able to compute energies.
Definition: system.h:482
Force field parameters of the system.
Definition: force_field.h:47
Auxilary type used to incapsulate the atom and its current coordinates Used internally in Selection::...
Definition: atom_proxy.h:1012
The system of atoms.
Definition: system.h:95
Selection class.
Definition: atom_proxy.h:57
Components of the non-bond energy.
Definition: system.h:40
const Eigen::Vector3f & XYZ(int ind, int fr) const
Read only access for given coordinate of given frame.
Definition: system.h:364
int num_atoms() const
Returns the number of atoms in selection.
Definition: system.h:193
Class which represents a single atom.
Definition: atom.h:33
const Atom & Atom_data(int ind) const
Read only access for given atom.
Definition: system.h:374
Atom & Atom_data(int ind)
Read/Write access for given atom.
Definition: system.h:369
int num_frames() const
Returns the number of frames in selection.
Definition: system.h:190
Eigen::Vector3f & XYZ(int ind, int fr)
Read/Write access for given coordinate of given frame.
Definition: system.h:359
Frame & Frame_data(int fr)
Get read/write reference for given frame.
Definition: system.h:379
float q_sr
Short-range Coloumb energy (within cut-off)
Definition: system.h:50
Force_field * get_force_field()
Returns pointer to internal Force_field object for direct manipulation or NULL if force field is not ...
Definition: system.h:488
Energy_components non_bond_energy(const Selection &sel1, const Selection &sel2, float cutoff=0.25, int fr=-1, bool periodic=true)
Definition: selection.cpp:1229
Generic API for reading and writing any molecule file formats.
Definition: mol_file.h:43
Periodic_box box
Periodic box.
Definition: system.h:71
float & Time(int fr)
Read/Write access to the time stamp of given frame.
Definition: system.h:349
const Frame & Frame_data(int fr) const
Get read only reference for given frame.
Definition: system.h:384
const float & Time(int fr) const
Read only access to the time stamp of given frame.
Definition: system.h:354
float time
Timestamp.
Definition: system.h:73
Energy_components operator+(const Energy_components &other)
Addition of energies.
Definition: system.cpp:863
float lj_sr
Short-range Lenard-Jones energy (within cut-off)
Definition: system.h:48
float q_14
Coloumb energy of 1-4 pairs.
Definition: system.h:46
Class encapsulating all operations with arbitrary triclinic periodic boxes This class stores the peri...
Definition: periodic_box.h:47