OpenStructure: stage/include/ost/mol/mm/topology.hh Source File

00001 //------------------------------------------------------------------------------
00002 // This file is part of the OpenStructure project <www.openstructure.org>
00003 //
00004 // Copyright (C) 2008-2015 by the OpenStructure authors
00005 //
00006 // This library is free software; you can redistribute it and/or modify it under
00007 // the terms of the GNU Lesser General Public License as published by the Free
00008 // Software Foundation; either version 3.0 of the License, or (at your option)
00009 // any later version.
00010 // This library is distributed in the hope that it will be useful, but WITHOUT
00011 // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
00012 // FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public License for more
00013 // details.
00014 //
00015 // You should have received a copy of the GNU Lesser General Public License
00016 // along with this library; if not, write to the Free Software Foundation, Inc.,
00017 // 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
00018 //------------------------------------------------------------------------------
00019 
00020 #ifndef OST_MM_TOPOLOGY_HH
00021 #define OST_MM_TOPOLOGY_HH
00022 
00023 #include <vector>
00024 #include <map>
00025 #include <set>
00026 
00027 #include <boost/shared_ptr.hpp>
00028 
00029 #include <ost/message.hh>
00030 #include <ost/mol/atom_handle.hh>
00031 #include <ost/mol/residue_handle.hh>
00032 #include <ost/mol/entity_handle.hh>
00033 #include <ost/mol/mm/index.hh>
00034 #include <ost/mol/mm/buildingblock.hh>
00035 #include <ost/mol/mm/block_modifiers.hh>
00036 #include <ost/mol/mm/forcefield.hh>
00037 #include <ost/mol/mm/index.hh>
00038 #include <ost/mol/xcs_editor.hh>
00039 #include <ost/mol/bond_handle.hh>
00040 #include <ost/mol/residue_prop.hh>
00041 #include <ost/img/image_handle.hh>
00042 
00043 #include <time.h>
00044 
00045 
00046 
00047 namespace ost { namespace mol{ namespace mm{
00048 
00049 class Topology;
00050 typedef boost::shared_ptr<Topology> TopologyPtr;
00051 
00052 
00053 class Topology{
00054 
00055 public:
00056 
00057   Topology(const std::vector<Real>& masses);
00058 
00059   Topology() { } //should not be accessible from Python to avoid messing around
00060                  //with empty topology
00061 
00062   static TopologyPtr Load(const String& filename);
00063 
00064   void Save(const String& filename);
00065 
00066   uint AddHarmonicBond(uint index_one,
00067                        uint index_two, 
00068                        Real bond_length,
00069                        Real force_constant);
00070 
00071   uint AddHarmonicAngle(uint index_one,
00072                         uint index_two,
00073                         uint index_three,
00074                         Real angle,
00075                         Real force_constant);
00076 
00077   uint AddUreyBradleyAngle(uint index_one,
00078                            uint index_two,
00079                            uint index_three,
00080                            Real angle,
00081                            Real angle_force_constant,
00082                            Real bond_length,
00083                            Real bond_force_constant);
00084 
00085   uint AddPeriodicDihedral(uint index_one,
00086                            uint index_two,
00087                            uint index_three,
00088                            uint index_four,
00089                            int multiplicity,
00090                            Real phase,
00091                            Real force_constant);
00092 
00093   uint AddPeriodicImproper(uint index_one,
00094                            uint index_two,
00095                            uint index_three,
00096                            uint index_four,
00097                            int multiplicity,
00098                            Real phase,
00099                            Real force_constant);
00100 
00101   uint AddHarmonicImproper(uint index_one,
00102                            uint index_two,
00103                            uint index_three,
00104                            uint index_four,
00105                            Real angle,
00106                            Real force_constant);
00107 
00108   uint AddCMap(uint index_one,
00109                uint index_two,
00110                uint index_three,
00111                uint index_four,
00112                uint index_five,
00113                int dimension,
00114                std::vector<Real> values);
00115 
00116   uint AddLJPair(uint index_one, 
00117                  uint index_two,
00118                  Real sigma,
00119                  Real epsilon);
00120 
00121 
00122   uint AddDistanceConstraint(uint index_one, 
00123                              uint index_two,
00124                              Real distance);
00125 
00126   uint AddExclusion(uint index_one,
00127                     uint index_two);
00128 
00129   void AddPositionConstraint(uint index);
00130 
00131   void SetDensity(const ost::img::ImageHandle& d, Real r, Real s = 1.0);
00132 
00133   void DefineDensityBody(const std::vector<int>& indices);
00134 
00135   int GetNumDensityBodies() const { return density_bodies_.size(); }
00136 
00137   void ResetPositionConstraints() { position_constraints_.clear(); }
00138 
00139   void ResetExclusions() { exclusions_.clear(); }
00140 
00141   uint AddHarmonicPositionRestraint(uint index, const geom::Vec3& ref_position, Real k, 
00142                                     Real x_scale = 1.0, Real y_scale = 1.0, Real z_scale = 1.0);
00143 
00144   uint AddHarmonicDistanceRestraint(uint index_one, uint index_two, 
00145                                     Real length, Real force_constant);
00146 
00147   uint AddFGMDHBondDonor(uint index_one, uint index_two,
00148                          Real length, Real k_length, Real alpha, 
00149                          Real k_alpha, Real beta, Real k_beta);
00150 
00151   uint AddFGMDHBondAcceptor(uint index_one, uint index_two);
00152 
00153   //Single atom parameters are expected to be set at once...
00154   void SetSigmas(const std::vector<Real>& sigmas);
00155 
00156   void SetSigma(uint index, Real sigma);
00157 
00158   void SetEpsilons(const std::vector<Real>& epsilons);
00159 
00160   void SetEpsilon(uint index, Real epsilon);
00161 
00162   void SetGBSARadii(const std::vector<Real>& gbsa_radii);
00163 
00164   void SetGBSARadius(uint index, Real radius);
00165 
00166   void SetOBCScalings(const std::vector<Real>& obc_scaling);
00167 
00168   void SetOBCScaling(uint index, Real scaling);
00169 
00170   void SetCharges(const std::vector<Real>& charges);
00171 
00172   void SetCharge(uint index, Real charge);
00173 
00174   void SetMasses(const std::vector<Real>& masses);
00175 
00176   void SetMass(uint index, Real mass);
00177 
00178   void SetFudgeQQ(Real fudge) { fudge_qq_ = fudge; }
00179 
00180   void SetFudgeLJ(Real fudge) { fudge_lj_ = fudge;} 
00181 
00182 
00183   void GetHarmonicBondParameters(uint index, uint& index_one, uint& index_two,
00184                                  Real& bond_length, Real& force_constant) const;
00185 
00186   void GetHarmonicAngleParameters(uint index, uint& index_one, uint& index_two, uint& index_three,
00187                                   Real& angle, Real& force_constant) const;
00188 
00189   void GetUreyBradleyAngleParameters(uint index, uint& index_one, uint& index_two, uint& index_three,
00190                                      Real& angle, Real& angle_force_constant, Real& bond_length, Real& bond_force_constant) const;
00191 
00192   void GetPeriodicDihedralParameters(uint index, uint& index_one, uint& index_two, uint& index_three, uint& index_four,
00193                                      int& multiplicity, Real& phase, Real& force_constant) const;
00194 
00195   void GetPeriodicImproperParameters(uint index, uint& index_one, uint& index_two, uint& index_three, uint& index_four,
00196                                      int& multiplicity, Real& phase, Real& force_constant) const;
00197 
00198   void GetHarmonicImproperParameters(uint index, uint& index_one, uint& index_two, uint& index_three, uint& index_four,
00199                                      Real& angle, Real& force_constant) const;
00200 
00201   void GetCMapParameters(uint index, uint& index_one, uint& index_two, uint& index_three, uint& index_four, uint& index_five,
00202                          int& dimension, std::vector<Real>& map) const;
00203 
00204   void GetLJPairParameters(uint index, uint& index_one, uint& index_two,
00205                            Real& sigma, Real& epsilon) const;
00206 
00207   void GetDistanceConstraintParameters(uint index, uint& index_one, uint& index_two,
00208                                        Real& distance) const;
00209 
00210   void GetHarmonicPositionRestraintParameters(uint index, uint& atom_index, geom::Vec3& ref_position,
00211                                               Real& k, Real& x_scale, Real& y_scale, Real& z_scale) const;
00212 
00213   void GetHarmonicDistanceRestraintParameters(uint index, uint& atom_one, uint& atom_two, Real& length,
00214                                               Real& force_constant) const;
00215 
00216   void GetFGMDHBondDonorParameters(uint index, uint& index_one, uint& index_two,
00217                                    Real& length, Real& k_length, Real& alpha, 
00218                                    Real& k_alpha, Real& beta, Real& k_beta) const;
00219 
00220   void GetFGMDHBondAcceptorParameters(uint index, uint& index_one, uint& index_two) const;
00221 
00222   ost::img::ImageHandle GetDensity() { return density_; }
00223 
00224   Real GetDensityResolution() { return density_resolution_; }
00225 
00226   Real GetDensityScaling() { return density_scaling_; }
00227 
00228   void SetHarmonicBondParameters(uint index, const Real bond_length, const Real force_constant);
00229 
00230   void SetHarmonicAngleParameters(uint index, const Real angle, const Real force_constant);
00231 
00232   void SetUreyBradleyAngleParameters(uint index, const Real angle, const Real angle_force_constant, 
00233                                      const Real bond_length, const Real bond_force_constant);
00234 
00235   void SetPeriodicDihedralParameters(uint index, const int multiplicity, 
00236                                      const Real phase, const Real force_constant);
00237 
00238   void SetPeriodicImproperParameters(uint index, const int multiplicity, const Real phase, 
00239                                      const Real force_constant);
00240 
00241   void SetHarmonicImproperParameters(uint index, const Real angle, const Real force_constant);
00242 
00243   void SetCMapParameters(uint index, const int dimension, const std::vector<Real>& map);
00244 
00245   void SetLJPairParameters(uint index, const Real sigma, const Real epsilon);
00246 
00247   void SetDistanceConstraintParameters(uint index, const Real distance);
00248 
00249   void SetHarmonicPositionRestraintParameters(uint index, const geom::Vec3& ref_position, Real k, 
00250                                                Real x_scale = 1.0, Real y_scale = 1.0, Real z_scale = 1.0);
00251 
00252   void SetHarmonicDistanceRestraintParameters(uint index, Real length, Real force_constant);
00253 
00254   void SetFGMDHBondDonorParameters(uint index, Real length, Real k_length, 
00255                                    Real alpha, Real k_alpha, Real beta, 
00256                                    Real k_beta);
00257 
00258   const std::vector<std::pair<Index<2>, std::vector<Real> > >& GetHarmonicBonds() const { return harmonic_bonds_; }
00259  
00260   const std::vector<std::pair<Index<3>, std::vector<Real> > >& GetHarmonicAngles() const { return harmonic_angles_; }
00261  
00262   const std::vector<std::pair<Index<3>, std::vector<Real> > >& GetUreyBradleyAngles() const { return urey_bradley_angles_; }
00263  
00264   const std::vector<std::pair<Index<4>, std::vector<Real> > >& GetPeriodicDihedrals() const { return periodic_dihedrals_; }
00265  
00266   const std::vector<std::pair<Index<4>, std::vector<Real> > >& GetPeriodicImpropers() const { return periodic_impropers_; }
00267  
00268   const std::vector<std::pair<Index<4>, std::vector<Real> > >& GetHarmonicImpropers() const { return harmonic_impropers_; }
00269  
00270   const std::vector<std::pair<Index<5>, std::vector<Real> > >& GetCMaps() const { return cmaps_; }
00271  
00272   const std::vector<std::pair<Index<2>, std::vector<Real> > >& GetLJPairs() const { return lj_pairs_; }
00273  
00274   const std::vector<std::pair<Index<2>, std::vector<Real> > >& GetDistanceConstraints() const { return distance_constraints_; }
00275 
00276   const std::vector<std::pair<Index<1>, std::vector<Real> > >& GetHarmonicPositionRestraints() const { return harmonic_position_restraints_; }
00277  
00278   const std::vector<std::pair<Index<2>, std::vector<Real> > >& GetHarmonicDistanceRestraints() const{ return harmonic_distance_restraints_; }
00279 
00280   const std::vector<Index<2> >& GetExclusions() const { return exclusions_; }
00281  
00282   const std::vector<uint>& GetPositionConstraints() const { return position_constraints_; }
00283 
00284   const std::vector<std::pair<Index<2>, std::vector<Real> > >& GetFGMDHBondDonors() const { return fgmd_hbond_donors_; }
00285 
00286   const std::vector<Index<2> >& GetFGMDHBondAcceptors() const { return fgmd_hbond_acceptors_; }
00287 
00288   const std::vector<std::vector<int> >& GetDensityBodies() const { return density_bodies_; }
00289 
00290 
00291   std::vector<Real> GetSigmas() const { return sigmas_; }
00292 
00293   std::vector<Real> GetEpsilons() const { return epsilons_; }
00294  
00295   std::vector<Real> GetGBSARadii() const { return gbsa_radii_; }
00296 
00297   std::vector<Real> GetOBCScalings() const { return obc_scaling_; }
00298  
00299   std::vector<Real> GetCharges() const { return charges_; }
00300 
00301   Real GetCharge(uint index) const;
00302 
00303   Real GetMass(uint index) const;
00304 
00305   Real GetSigma(uint index) const;
00306 
00307   Real GetEpsilon(uint index) const;
00308 
00309   Real GetGBSARadius(uint index) const;
00310 
00311   Real GetOBCScaling(uint index) const;
00312 
00313   std::vector<Real> GetMasses() const { return atom_masses_; }
00314 
00315   Real GetFudgeQQ() const { return fudge_qq_; }
00316 
00317   Real GetFudgeLJ() const { return fudge_lj_; }
00318 
00319   std::vector<uint> GetHarmonicBondIndices(uint index_one,
00320                                             uint index_two) const;
00321 
00322   std::vector<uint> GetHarmonicAngleIndices(uint index_one,
00323                                              uint index_two,
00324                                              uint index_three) const;
00325 
00326   std::vector<uint> GetUreyBradleyAngleIndices(uint index_one,
00327                                                 uint index_two,
00328                                                 uint index_three) const;
00329 
00330   std::vector<uint> GetPeriodicDihedralIndices(uint index_one,
00331                                                 uint index_two,
00332                                                 uint index_three,
00333                                                 uint index_four) const;
00334 
00335   std::vector<uint> GetPeriodicImproperIndices(uint index_one,
00336                                                 uint index_two,
00337                                                 uint index_three,
00338                                                 uint index_four) const;
00339 
00340   std::vector<uint> GetHarmonicImproperIndices(uint index_one,
00341                                                 uint index_two,
00342                                                 uint index_three,
00343                                                 uint index_four) const;
00344 
00345   std::vector<uint> GetCMapIndices(uint index_one,
00346                                     uint index_two,
00347                                     uint index_three,
00348                                     uint index_four,
00349                                     uint index_five) const;
00350 
00351   int GetLJPairIndex(uint index_one,
00352                      uint index_two) const;
00353 
00354   int GetDistanceConstraintIndex(uint index_one,
00355                                  uint index_two) const;
00356 
00357 
00358 
00359   std::vector<uint> GetHarmonicDistanceRestraintIndices(uint index_one,
00360                                                         uint index_two) const;
00361 
00362   std::vector<uint> GetFGMDHBondDonorIndices(uint index_one,
00363                                              uint index_two) const;
00364 
00365   std::vector<uint> GetFGMDHBondAcceptorIndices(uint index_one,
00366                                                 uint index_two) const;
00367 
00368   std::vector<uint> GetHarmonicBondIndices(uint atom_index) const;
00369 
00370   std::vector<uint> GetHarmonicAngleIndices(uint atom_index) const;
00371 
00372   std::vector<uint> GetUreyBradleyAngleIndices(uint atom_index) const;
00373 
00374   std::vector<uint> GetPeriodicDihedralIndices(uint atom_index) const;
00375 
00376   std::vector<uint> GetPeriodicImproperIndices(uint atom_index) const;
00377 
00378   std::vector<uint> GetHarmonicImproperIndices(uint atom_index) const;
00379 
00380   std::vector<uint> GetCMapIndices(uint atom_index) const;
00381 
00382   std::vector<uint> GetLJPairIndices(uint atom_index) const;
00383 
00384   std::vector<uint> GetDistanceConstraintIndices(uint atom_index) const;
00385 
00386   std::vector<uint> GetHarmonicPositionRestraintIndices(uint atom_index) const;
00387 
00388   std::vector<uint> GetHarmonicDistanceRestraintIndices(uint atom_index) const;
00389 
00390   std::vector<uint> GetFGMDHBondDonorIndices(uint atom_index) const;
00391 
00392   std::vector<uint> GetFGMDHBondAcceptorIndices(uint atom_index) const;
00393 
00394   uint GetNumParticles() { return num_particles_; }
00395 
00396   uint GetNumHarmonicBonds() { return harmonic_bonds_.size(); }
00397 
00398   uint GetNumHarmonicAngles() { return harmonic_angles_.size(); }
00399 
00400   uint GetNumUreyBradleyAngles() { return urey_bradley_angles_.size(); }
00401 
00402   uint GetNumPeriodicDihedrals() { return periodic_dihedrals_.size(); }
00403 
00404   uint GetNumPeriodicImpropers() { return periodic_impropers_.size(); }
00405 
00406   uint GetNumHarmonicImpropers() { return harmonic_impropers_.size(); }
00407 
00408   uint GetNumCMaps() { return cmaps_.size(); }
00409 
00410   uint GetNumLJPairs() { return lj_pairs_.size(); }
00411 
00412   uint GetNumDistanceConstraints() { return distance_constraints_.size(); }
00413 
00414   uint GetNumPositionConstraints() { return position_constraints_.size(); }
00415 
00416   uint GetNumHarmonicPositionRestraints() { return harmonic_position_restraints_.size(); }
00417 
00418   uint GetNumHarmonicDistanceRestraints() { return harmonic_distance_restraints_.size();}
00419 
00420   uint GetNumExclusions() { return exclusions_.size(); }
00421 
00422   uint GetNumFGMDHBondDonors() { return fgmd_hbond_donors_.size(); }
00423 
00424   uint GetNumFGMDHBondAcceptors() { return fgmd_hbond_acceptors_.size(); }
00425 
00426   void Merge(ost::mol::EntityHandle& ent, TopologyPtr other, const ost::mol::EntityHandle& other_ent);
00427 
00428   void Merge(TopologyPtr other);
00429 
00430   template <typename DS>
00431   void Serialize(DS& ds){
00432 
00433 
00434     uint num_items = 0;
00435     Index<2> actual_index;
00436 
00437     ds & num_particles_;
00438     ds & fudge_qq_;
00439     ds & fudge_lj_;
00440 
00441     if(ds.IsSource()){
00442       ds & num_items;
00443       atom_masses_ = std::vector<Real>(num_items);
00444       ds & num_items;
00445       sigmas_ = std::vector<Real>(num_items);
00446       ds & num_items;
00447       epsilons_ = std::vector<Real>(num_items);
00448       ds & num_items;
00449       gbsa_radii_ = std::vector<Real>(num_items);
00450       ds & num_items;
00451       obc_scaling_ = std::vector<Real>(num_items);
00452       ds & num_items;
00453       charges_ = std::vector<Real>(num_items);
00454       ds & num_items;
00455       position_constraints_ = std::vector<uint>(num_items);
00456     }
00457     else{
00458       num_items = atom_masses_.size();
00459       ds & num_items;
00460       num_items = sigmas_.size();
00461       ds & num_items;
00462       num_items = epsilons_.size();
00463       ds & num_items;
00464       num_items = gbsa_radii_.size();
00465       ds & num_items;
00466       num_items = obc_scaling_.size();
00467       ds & num_items;
00468       num_items = charges_.size();
00469       ds & num_items;
00470       num_items = position_constraints_.size();
00471       ds & num_items;
00472     }
00473 
00474     for(std::vector<Real>::iterator i = atom_masses_.begin();
00475         i != atom_masses_.end(); ++i){
00476       ds & *i;
00477     }
00478 
00479     for(std::vector<Real>::iterator i = sigmas_.begin();
00480         i != sigmas_.end(); ++i){
00481       ds & *i;
00482     }
00483 
00484     for(std::vector<Real>::iterator i = epsilons_.begin();
00485         i != epsilons_.end(); ++i){
00486       ds & *i;
00487     }
00488 
00489     for(std::vector<Real>::iterator i = gbsa_radii_.begin();
00490         i != gbsa_radii_.end(); ++i){
00491       ds & *i;
00492     }
00493 
00494     for(std::vector<Real>::iterator i = obc_scaling_.begin();
00495         i != obc_scaling_.end(); ++i){
00496       ds & *i;
00497     }
00498 
00499     for(std::vector<Real>::iterator i = charges_.begin();
00500         i != charges_.end(); ++i){
00501       ds & *i;
00502     }
00503 
00504     for(std::vector<uint>::iterator i = position_constraints_.begin();
00505         i != position_constraints_.end(); ++i){
00506       ds & *i;
00507     }
00508 
00509     if(ds.IsSource()){
00510       ds & num_items;
00511       for(uint i = 0; i < num_items; ++i){
00512         harmonic_bonds_.push_back(std::make_pair(Index<2>(),std::vector<Real>(2)));
00513       }
00514 
00515       ds & num_items;
00516       for(uint i = 0; i < num_items; ++i){
00517         harmonic_angles_.push_back(std::make_pair(Index<3>(),std::vector<Real>(2)));
00518       }
00519 
00520       ds & num_items;
00521       for(uint i = 0; i < num_items; ++i){
00522         urey_bradley_angles_.push_back(std::make_pair(Index<3>(),std::vector<Real>(4)));
00523       }
00524 
00525       ds & num_items;
00526       for(uint i = 0; i < num_items; ++i){
00527         periodic_dihedrals_.push_back(std::make_pair(Index<4>(),std::vector<Real>(3)));
00528       }
00529 
00530       ds & num_items;
00531       for(uint i = 0; i < num_items; ++i){
00532         periodic_impropers_.push_back(std::make_pair(Index<4>(),std::vector<Real>(3)));
00533       }
00534 
00535       ds & num_items;
00536       for(uint i = 0; i < num_items; ++i){
00537         harmonic_impropers_.push_back(std::make_pair(Index<4>(),std::vector<Real>(2)));
00538       }
00539 
00540       ds & num_items;
00541       for(uint i = 0; i < num_items; ++i){
00542         uint num_values = 0;
00543         ds & num_values;
00544         cmaps_.push_back(std::make_pair(Index<5>(),std::vector<Real>(num_values)));
00545       }
00546 
00547       ds & num_items;
00548       for(uint i = 0; i < num_items; ++i){
00549         lj_pairs_.push_back(std::make_pair(Index<2>(),std::vector<Real>(2)));
00550       } 
00551 
00552       ds & num_items;
00553       for(uint i = 0; i < num_items; ++i){
00554         distance_constraints_.push_back(std::make_pair(Index<2>(),std::vector<Real>(1)));
00555       }
00556 
00557       ds & num_items;
00558       exclusions_ = std::vector<Index<2> >(num_items);
00559 
00560       ds & num_items;
00561       for(uint i = 0; i < num_items; ++i){
00562         harmonic_position_restraints_.push_back(std::make_pair(Index<1>(),std::vector<Real>(7)));
00563       }
00564 
00565       ds & num_items;
00566       for(uint i = 0; i < num_items; ++i){
00567         harmonic_distance_restraints_.push_back(std::make_pair(Index<2>(),std::vector<Real>(2)));
00568       }
00569 
00570       ds & num_items;
00571       for(uint i = 0; i < num_items; ++i){
00572         fgmd_hbond_donors_.push_back(std::make_pair(Index<2>(),std::vector<Real>(6)));
00573       }
00574 
00575       ds & num_items;
00576       for(uint i = 0; i < num_items; ++i){
00577         fgmd_hbond_acceptors_.push_back(Index<2>());
00578       }
00579     }
00580     else{
00581       num_items = harmonic_bonds_.size();
00582       ds & num_items;
00583       num_items = harmonic_angles_.size();
00584       ds & num_items;
00585       num_items = urey_bradley_angles_.size();
00586       ds & num_items;
00587       num_items = periodic_dihedrals_.size();
00588       ds & num_items;
00589       num_items = periodic_impropers_.size();
00590       ds & num_items;
00591       num_items = harmonic_impropers_.size();
00592       ds & num_items;
00593       num_items = cmaps_.size();
00594       ds & num_items;
00595       for(uint i = 0; i < cmaps_.size(); ++i){
00596         num_items = cmaps_[i].second.size();
00597         ds & num_items;
00598       }
00599       num_items = lj_pairs_.size();
00600       ds & num_items;
00601       num_items = distance_constraints_.size();
00602       ds & num_items;
00603       num_items = exclusions_.size();
00604       ds & num_items;
00605       num_items = harmonic_position_restraints_.size();
00606       ds & num_items;
00607       num_items = harmonic_distance_restraints_.size();
00608       ds & num_items;
00609       num_items = fgmd_hbond_donors_.size();
00610       ds & num_items;
00611       num_items = fgmd_hbond_acceptors_.size();
00612       ds & num_items;
00613     }
00614 
00615     for(std::vector<std::pair<Index<2>,std::vector<Real> > >::iterator i = harmonic_bonds_.begin();
00616         i != harmonic_bonds_.end(); ++i){
00617       ds & i->first;
00618       ds & i->second[0];
00619       ds & i->second[1];
00620     }
00621 
00622     for(std::vector<std::pair<Index<3>,std::vector<Real> > >::iterator i = harmonic_angles_.begin();
00623         i != harmonic_angles_.end(); ++i){
00624       ds & i->first;
00625       ds & i->second[0];
00626       ds & i->second[1];
00627     }
00628 
00629     for(std::vector<std::pair<Index<3>,std::vector<Real> > >::iterator i = urey_bradley_angles_.begin();
00630         i != urey_bradley_angles_.end(); ++i){
00631       ds & i->first;
00632       ds & i->second[0];
00633       ds & i->second[1];
00634       ds & i->second[2];
00635       ds & i->second[3];
00636     }
00637 
00638     for(std::vector<std::pair<Index<4>,std::vector<Real> > >::iterator i = periodic_dihedrals_.begin();
00639         i != periodic_dihedrals_.end(); ++i){
00640       ds & i->first;
00641       ds & i->second[0];
00642       ds & i->second[1];
00643       ds & i->second[2];
00644     }
00645 
00646     for(std::vector<std::pair<Index<4>,std::vector<Real> > >::iterator i = periodic_impropers_.begin();
00647         i != periodic_impropers_.end(); ++i){
00648       ds & i->first;
00649       ds & i->second[0];
00650       ds & i->second[1];
00651       ds & i->second[2];
00652     }
00653 
00654     for(std::vector<std::pair<Index<4>,std::vector<Real> > >::iterator i = harmonic_impropers_.begin();
00655         i != harmonic_impropers_.end(); ++i){
00656       ds & i->first;
00657       ds & i->second[0];
00658       ds & i->second[1];
00659     }
00660 
00661     for(std::vector<std::pair<Index<5>,std::vector<Real> > >::iterator i = cmaps_.begin();
00662         i != cmaps_.end(); ++i){
00663       ds & i->first;
00664       for(std::vector<Real>::iterator j = i->second.begin();
00665           j != i->second.end(); ++j){
00666         ds & (*j);
00667       }
00668     }
00669 
00670     for(std::vector<std::pair<Index<2>,std::vector<Real> > >::iterator i = lj_pairs_.begin();
00671         i != lj_pairs_.end(); ++i){
00672       ds & i->first;
00673       ds & i->second[0];
00674       ds & i->second[1];
00675     }
00676 
00677     for(std::vector<std::pair<Index<2>,std::vector<Real> > >::iterator i = distance_constraints_.begin();
00678         i != distance_constraints_.end(); ++i){
00679       ds & i->first;
00680       ds & i->second[0];
00681     }
00682 
00683     for(std::vector<Index<2> >::iterator i = exclusions_.begin();
00684         i != exclusions_.end(); ++i){
00685       ds & (*i);
00686     }
00687 
00688     for(std::vector<std::pair<Index<1>,std::vector<Real> > >::iterator i = harmonic_position_restraints_.begin();
00689         i != harmonic_position_restraints_.end(); ++i){
00690       ds & i->first;
00691       ds & i->second[0];
00692       ds & i->second[1];
00693       ds & i->second[2];
00694       ds & i->second[3];
00695       ds & i->second[4];
00696       ds & i->second[5];
00697       ds & i->second[6];
00698     }
00699 
00700     for(std::vector<std::pair<Index<2>,std::vector<Real> > >::iterator i = harmonic_distance_restraints_.begin();
00701         i != harmonic_distance_restraints_.end(); ++i){
00702       ds & i->first;
00703       ds & i->second[0];
00704       ds & i->second[1];
00705     }
00706 
00707     for(std::vector<std::pair<Index<2>,std::vector<Real> > >::iterator i = fgmd_hbond_donors_.begin();
00708         i != fgmd_hbond_donors_.end(); ++i){
00709       ds & i->first;
00710       ds & i->second[0];
00711       ds & i->second[1];
00712       ds & i->second[2];
00713       ds & i->second[3];
00714       ds & i->second[4];
00715       ds & i->second[5];
00716     }
00717 
00718     for(std::vector<Index<2> >::iterator i = fgmd_hbond_acceptors_.begin();
00719         i != fgmd_hbond_acceptors_.end(); ++i){
00720       ds & (*i);
00721     }
00722 
00723     if(ds.IsSource()){
00724       ds & num_items;
00725       for(uint i = 0; i < num_items; ++i){
00726         ds & actual_index;
00727         added_lj_pairs_.insert(actual_index);
00728       }
00729 
00730       ds & num_items;
00731       for(uint i = 0; i < num_items; ++i){
00732         ds & actual_index;
00733         added_distance_constraints_.insert(actual_index);
00734       }
00735 
00736       ds & num_items;
00737       for(uint i = 0; i < num_items; ++i){
00738         ds & actual_index;
00739         added_exclusions_.insert(actual_index);
00740       }
00741     }
00742     else{
00743       num_items = added_lj_pairs_.size();
00744       ds & num_items;
00745       for(std::set<Index<2> >::iterator i = added_lj_pairs_.begin(); 
00746           i != added_lj_pairs_.end(); ++i){
00747         actual_index = *i;
00748         ds & actual_index;
00749       }
00750       num_items = added_distance_constraints_.size();
00751       ds & num_items;
00752       for(std::set<Index<2> >::iterator i = added_distance_constraints_.begin(); 
00753           i != added_distance_constraints_.end(); ++i){
00754         actual_index = *i;
00755         ds & actual_index;
00756       }
00757       num_items = added_exclusions_.size();
00758       ds & num_items;
00759       for(std::set<Index<2> >::iterator i = added_exclusions_.begin(); 
00760           i != added_exclusions_.end(); ++i){
00761         actual_index = *i;
00762         ds & actual_index;
00763       }
00764     }
00765   }
00766 
00767 private:
00768 
00769   void MergeTop(TopologyPtr other);
00770 
00771   void MergeEnt(ost::mol::EntityHandle& ent, const ost::mol::EntityHandle& other_ent);
00772 
00773   void CheckTopToAdd(TopologyPtr other);
00774 
00775   void CheckEntToAdd(ost::mol::EntityHandle& ent, TopologyPtr other, 
00776                      const ost::mol::EntityHandle& other_ent);
00777 
00778   uint num_particles_;
00779 
00780   //fudge parameters for lj 1,4 pairs
00781   Real fudge_qq_;
00782   Real fudge_lj_;
00783 
00784   //single atom parameters
00785   std::vector<Real> atom_masses_;
00786   std::vector<Real> sigmas_;
00787   std::vector<Real> epsilons_;
00788   std::vector<Real> gbsa_radii_;
00789   std::vector<Real> obc_scaling_;
00790   std::vector<Real> charges_;
00791   std::vector<uint> position_constraints_;
00792 
00793   //interactions with multiple atoms involved
00794   std::vector<std::pair<Index<2>,std::vector<Real> > > harmonic_bonds_;
00795   std::vector<std::pair<Index<3>,std::vector<Real> > > harmonic_angles_;
00796   std::vector<std::pair<Index<3>,std::vector<Real> > > urey_bradley_angles_;
00797   std::vector<std::pair<Index<4>,std::vector<Real> > > periodic_dihedrals_;
00798   std::vector<std::pair<Index<4>,std::vector<Real> > > periodic_impropers_;
00799   std::vector<std::pair<Index<4>,std::vector<Real> > > harmonic_impropers_;
00800   std::vector<std::pair<Index<5>,std::vector<Real> > > cmaps_;
00801   std::vector<std::pair<Index<2>,std::vector<Real> > > lj_pairs_;
00802   std::vector<std::pair<Index<2>,std::vector<Real> > > distance_constraints_;
00803   std::vector<Index<2> > exclusions_;
00804   std::vector<std::pair<Index<1>,std::vector<Real> > > harmonic_position_restraints_;
00805   std::vector<std::pair<Index<2>,std::vector<Real> > > harmonic_distance_restraints_;
00806   std::vector<std::pair<Index<2>,std::vector<Real> > > fgmd_hbond_donors_;
00807   std::vector<Index<2> > fgmd_hbond_acceptors_;
00808 
00809   //the atoms of the interactions, that should be unique get tracked in here
00810   //note, that this is waste of memory, needs better implementation
00811   std::set<Index<2> > added_lj_pairs_;
00812   std::set<Index<2> > added_distance_constraints_;
00813   std::set<Index<2> > added_exclusions_;
00814 
00815   //density map
00816   ost::img::ImageHandle density_;
00817   Real density_resolution_;
00818   Real density_scaling_;
00819   std::vector<std::vector<int> > density_bodies_;
00820   std::vector<bool> in_density_body_;
00821 };
00822 
00823 
00824 }}} //ns
00825 
00826 #endif