B bB@s.dZddlmZmZddlZddlZddlZddlmZddl m Z ddl Z ddl m Zddl mZdd lmZmZdd lmZdd lmZmZmZmZdd lmZmZmZmZmZm Z m!Z!m"Z"m#Z#m$Z$m%Z%m&Z&m'Z'm(Z(m)Z)m*Z*m+Z+m,Z,m-Z-m.Z.m/Z/m0Z0m1Z1m2Z2m3Z3m4Z4m5Z5m6Z6m7Z7dd l8m9Z9m:Z:m;Z;ddlm?Z?m@Z@mAZAddlBmCZCmDZDddlEmFZFy(ddlGmHZHddlImJZKddlLmMZMWneNk rdZHZKYnXeOePZQddZRddZSddZTddZUddZVe=dd ZWGd!d"d"eXZYGd#d$d$eXZZGd%d&d&eXZ[d'd(Z\d)d*Z]dS)+z This module contains the core base class for all of the chemical structures with various topological and force field features. )absolute_importdivisionN) defaultdict)copy)unit)residue) DEG_TO_RADSMALL)ParameterError)STANDARD_BOND_LENGTHS_SQUARED!box_lengths_and_angles_to_vectors!box_vectors_to_lengths_and_angles distance2) AcceptorDonorAngleAtomAtomListBond ChiralFrameCmapDihedral DihedralTypeDihedralTypeList ExtraPointGroupImproperMultipoleFrameNonbondedException NoUreyBradleyOutOfPlaneBendOutOfPlaneExtraPointFrame PiTorsion ResidueList StretchBendThreeParticleExtraPointFrameTorsionTorsion TrackedList TrigonalAngleTwoParticleExtraPointFrameUnassignedAtomType UreyBradleyLink)PYPYfind_atom_pairs tag_molecules) needs_openmm) integer_types iteritems string_types)rangezip)Vec3)app)openmm)reducePeriodicBoxVectorsc Csg}x|D]}t|rN|jtjr:||tjq||tjq t|t rbt dq yt |t |}Wnt k rYnX||q W|S)NzUnit cell cannot have strings) u is_quantityrZ is_compatibledegreeappend value_in_unit angstroms isinstancer3 TypeErroriter_strip_box_units)argsZnew_argsargrF/lib/python3.7/site-packages/parmed/structure.pyrC0s     rCcCsX|jdkrdS|jdkrdS|jdkr*dS|jdkr8dS|jd krFd S|jd krTd Sd S)Ng333333?rg333333?g?g@g?g?g?) atomic_number)atomrFrFrG_bondiFs      rOcCs@|jdkr8|j}|djdkr"dS|djdkr4dSdSt|S)Nrr)rHrIg?)rLg?g333333?)rM bond_partnersrO)rNZbondedsrFrFrG_mbondiOs rRcCs*|jdkr"|jdjdkrdSdSt|S)NrrrIg?g333333?)rMrQrO)rNrFrFrG_mbondi2Ys  rScCsh|jjdkr*|jds$|jdrRdSn(|jjdkrR|jdsN|jdrRdS|jd kr`dSt|S) N)ZGLUZASPZGL4ZAS4ZOEZODgffffff?ZARGZHHZHEgQ?ZOXT)rname startswithrS)rNrFrFrG_mbondi3`s   rVcCs|tjtjfkrt|}nt|}|jdkr4|dfS|jdkrF|dfS|jdkrX|dfS|jdkrj|dfS|jdkr||d fS|jd kr|d fS|jd kr|d fS|dfS)a Gets the default GB parameters for a given atom according to a specific Generalized Born model Parameters ---------- atom : :class:`Atom` The atom to get the default GB parameters for model : ``app.HCT, app.OBC1, or app.OBC2`` The GB model to get the default parameters for (app.GBn and app.GBn2 are already handled in Structure._get_gb_parameters) Returns ------- radius, screen [,alpha, beta, gamma] : ``float, float [,float, float, float]`` The intrinsic radius of the atom and the screening factor of the atom. If the model is GBn2, alpha, beta, and gamma parameters are also returned rg333333?rHg ףp= ?rIgHzG?rP g)\(?rKgQ?rLgQ?g?)r7OBC1OBC2rSrRrM)rNZmodelZradrFrFrG_gb_rad_screenks$        rZc@seZdZdZdZdZdZdZdZdZ dZ d Z d Z d Z d Zd ZdZddZddZdddZddZddZdddZddZdd Zd!d"Zd#d$Zd%d&Zd'd(Zd)d*Zd+d,Zd-d.Zd/d0Z d1d2Z!e"d3d4Z#d5d6Z$dd8d9Z%d:d;Z&e"dd=Z'e"e)d?d@Z*e"dAdBZ+e+j(dCdBZ+e"dDdEZ,e,j(dFdEZ,ddHdIZ-e"dJdKZ.e.j(dLdKZ.ddMdNZ/e"dOdPZ0e0j(dQdPZ0dRdSZ1e"dTdUZ2e2j(dVdUZ2e)d7dWe3j4dXe3j4d7dYd7d7dXe3j5e3j6dZe3j7d[d\ddYd7d]dYddfd^d_Z8e)d`daZ9e)dbdcZ:e)ddddeZ;e)ddfdgZe)dldmZ?e)dndoZ@e)dpdqZAe)d7d e3j4de3j4d]d\fdrdsZBdtduZCdvdwZDe)d7dxe3j4d[d\d7dXe3j5e3jEdZe3j7dYfdydzZFe)d{d|ZGe)d}d~ZHe)ddZIe)ddZJe)ddZKddZLddZMdddZNddZOddZPddZQddZRddZSddZTddZUddZVddZWddZXddZYddZZe)ddZ[ddZ\ddZ]ddZ^ddZ_e_Z`dddZaddZbebZcedddZeddZfddZgd7S) Structurea A chemical structure composed of atoms, bonds, angles, torsions, and other topological features Attributes ---------- atoms : :class:`AtomList` List of all atoms in the structure residues : :class:`ResidueList` List of all residues in the structure bonds : :class:`TrackedList` (:class:`Bond`) List of all bonds in the structure angles : :class:`TrackedList` (:class:`Angle`) List of all angles in the structure dihedrals : :class:`TrackedList` (:class:`Dihedral`) List of all dihedrals in the structure rb_torsions : :class:`TrackedList` (:class:`Dihedral`) List of all Ryckaert-Bellemans torsions in the structure urey_bradleys : :class:`TrackedList` (:class:`UreyBradley`) List of all Urey-Bradley angle bends in the structure impropers : :class:`TrackedList` (:class:`Improper`) List of all CHARMM-style improper torsions in the structure cmaps : :class:`TrackedList` (:class:`Cmap`) List of all CMAP objects in the structure trigonal_angles : :class:`TrackedList` (:class:`TrigonalAngle`) List of all AMOEBA-style trigonal angles in the structure out_of_plane_bends : :class:`TrackedList` (:class:`OutOfPlaneBends`) List of all AMOEBA-style out-of-plane bending angles pi_torsions : :class:`TrackedList` (:class:`PiTorsion`) List of all AMOEBA-style pi-torsion angles stretch_bends : :class:`TrackedList` (:class:`StretchBend`) List of all AMOEBA-style stretch-bend compound bond/angle terms torsion_torsions : :class:`TrackedList` (:class:`TorsionTorsion`) List of all AMOEBA-style coupled torsion-torsion terms chiral_frames : :class:`TrackedList` (:class:`ChiralFrame`) List of all AMOEBA-style chiral frames defined in the structure multipole_frames : :class:`TrackedList` (:class:`MultipoleFrame`) List of all AMOEBA-style multipole frames defined in the structure adjusts : :class:`TrackedList` (:class:`NonbondedException`) List of all nonbonded pair-exception rules acceptors : :class:`TrackedList` (:class:`AcceptorDonor`) List of all H-bond acceptors, if that information is present donors : :class:`TrackedList` (:class:`AcceptorDonor`) List of all H-bond donors, if that information is present groups : :class:`TrackedList` (:class:`Group`) List of all CHARMM-style GROUP objects (whatever those are used for) box : ``list of 6 floats`` Box dimensions (a, b, c, alpha, beta, gamma) for the unit cell. If no box is defined, `box` is set to `None` space_group : ``str`` The space group of the structure (default is "P 1") nrexcl : ``int`` The number of bonds away that an atom can be in order to be excluded from the direct nonbonded summation title : ``str`` Cosmetic only, it is an arbitrary title assigned to the system. Default value is an empty string positions : u.Quantity(list(Vec3), u.angstroms) Unit-bearing atomic coordinates. If not all atoms have coordinates, this property is None coordinates : np.ndarray of shape (nframes, natom, 3) If no coordinates are set, this is set to None. The first frame will match the coordinates present on the atoms. symmetry : :class:`Symmetry` if no symmetry is set, this is set to None. links : :class:`TrackedList` (:class:`Link`) The list of Link definitions for this Structure Notes ----- This class also has a handful of type lists for each of the attributes above (excluding `atoms`, `residues`, `chiral_frames`, and `multipole_frames`). They are all TrackedList instances that are designed to hold the relevant parameter type. The list is: bond_types, angle_types, dihedral_types, urey_bradley_types, improper_types, cmap_types, trigonal_angle_types, out_of_plane_bend_types, pi_torsion_types, stretch_bend_types, torsion_torsion_types, adjust_types dihedral_types _may_ be a list of :class:`DihedralType` instances, since torsion profiles are often represented by a Fourier series with multiple terms rrrHrIrPrW cCsDt|_t|_t|_t|_t|_t|_t|_ t|_ t|_ t|_ t|_ t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_t|_ t|_!t|_"t|_#t|_$d|_%d|_&d|_'d|_(d|_)d|_*d|_+d|_,dS)NzP 1Fr]lorentz)-ratomsr#residuesr'bondsangles dihedrals rb_torsions urey_bradleys improperscmapstrigonal_anglesout_of_plane_bends pi_torsions stretch_bendstorsion_torsions chiral_framesmultipole_framesadjusts acceptorsdonorsgroups bond_types angle_typesdihedral_typesurey_bradley_typesimproper_typesrb_torsion_types cmap_typestrigonal_angle_typesout_of_plane_bend_typespi_torsion_typesstretch_bend_typestorsion_torsion_types adjust_typeslinks_box _coordinates space_groupunknown_functionalnrexcltitle_combining_rulesymmetry)selfrFrFrG__init__sTzStructure.__init__cCs&t|j}t|j}tdd|jD}dt|j|fg}|dkrR|d||d|t|j}|d||jdk r|dt |jd d t kst |jd d t kst |jd d t kr|d n |dt|jdkr|jdjdk r|dn |dd |S)NcSsg|]}t|tqSrF)r@r).0arFrFrG ,sz&Structure.__repr__..z <%s %d atomsrz [%d EPs]z ; %d residuesz ; %d bondsz; PBCr]Zr^r_z (triclinic)z (orthogonal)z; parametrized>z; NOT parametrized>rc) lenrerfsumtype__name__r=rgboxabsr join)rnatomnresnextraretstrnbondrFrFrG__repr__)s&     ,  "  zStructure.__repr__rccCs&|j|||||||j|dS)a Adds a new atom to the Structure, adding a new residue to `residues` if it has a different name or number as the last residue added and adding it to the `atoms` list. Parameters ---------- atom : :class:`Atom` The atom to add to this residue list resname : ``str`` The name of the residue this atom belongs to resnum : ``int`` The number of the residue this atom belongs to chain : ``str`` The chain ID character for this residue inscode : ``str`` The insertion code ID character for this residue (it is stripped) segid : ``str`` The segment identifier for this residue (it is stripped) Notes ----- If the residue name and number differ from the last residue in this list, a new residue is added and the atom is added to that residue N)rfadd_atomrer=)rrNZresnameZresnumchainZinscodesegidrFrFrGrDszStructure.add_atomcCsb|j|jk rtd|jd}|||jr<||jdkrJ|j|n|j|jd|dS)a Adds a new atom to the Structure at the end if the given residue Parameters ---------- atom : :class:`Atom` The atom to add to the system residue : :class:`Residue` The residue to which to add this atom. It MUST be part of this Structure instance already or a ValueError is raised Notes ----- This atom is added at the end of the residue and is inserted into the `atoms` list in such a way that all residues are composed of atoms contiguous in the atoms list. For large systems, this may be a relatively expensive operation z$Residue is not part of the structurerN)listrf ValueErrorrerr=insertidx)rrNrZ last_atomrFrFrGadd_atom_to_residuecs   zStructure.add_atom_to_residuecCs|t|S)z A deep copy of the Structure )rr)rrFrFrG__copy__szStructure.__copy__Fc* Cs |}x:|jD]0}|j}t|}|||j|j|j|j|jqWx|j D]}|j t|qJW|j x|j D]}|j t|qtW|j |r,d} i} xt |jD]v\} } t| drxb| D].} |j t| | | g | | d7} qWq| | g | | d7} |j t| qWn"x |jD]} |j t| q4W|j x |jD]}|j t|q`W|j x |jD]}|j t|qW|j x |jD]}|j t|qW|j x |jD]}|j t|qW|j x |jD]}|j t|qW|j x |jD]}|j t|q||#jj||#jj||#j$j|#j o|j|#j jqLWxd|j?D]Z} |j? t@|| jj|| jj|| j$j|| j*j|| j5j| j o|j| j jqWx6|jAD],}$|jA tB||$jj||$jj|$jC qWx8|jDD].}%|jD tE||%jFj|%jG|%jH|%jI|%jJ q|jSD]4}'|jS tT||'jj||'jj|'jU|'jV|'jW qXWt|jX|_Xt|jY|_Y|jZ|_Zx$t[|j\|j\D]\}(})|)j]|(_] qW|S)a Makes a copy of the current structure as an instance of a specified subclass Parameters ---------- cls : Structure subclass The returned object is a copy of this structure as a `cls` instance split_dihedrals : ``bool`` If True, then the Dihedral entries will be split up so that each one is paired with a single DihedralType (rather than a DihedralTypeList) Returns ------- *cls* instance The instance of the Structure subclass `cls` with a copy of the current Structure's topology information r__iter__rr)improper ignore_endrN)r)^rerrrrTnumberrinsertion_coderryr=claimrz enumerater{hasattr setdefaultr|r}r~rrrrrrrrgratom1ratom2rfunctrhratom3rir4rrratom4rZ_functrkrr+rlrrjrmratom5rnr(ror rpr"atom6rqr$rrr&rsr chiralityrtrrN frame_pt_numvectailvecheadnvecrurrvrrwrxrmoverr,lengthZ symmetry_op1Z symmetry_op2rrcombining_ruler5rfter)*rclsZsplit_dihedralscrNresrZbtatZndtZmapdtZidtZdttutitZrtctZtaZotZptstttrebdiZieZtitypubrcmopsZchmglZr1Zr2rFrFrGrsh "                                        $                      " " &   zStructure.copycCsddlm}||S)aM Generates a DataFrame from the current Structure's atomic properties Returns ------- df : DataFrame DataFrame with all atomic properties See Also -------- :func:`parmed.utils.pandautils.create_dataframe` for full documentation of the generated DataFrame r)create_dataframe)parmed.utils.pandautilsr)rrrFrFrG to_dataframeqs zStructure.to_dataframecCsddlm}|||S)av Loads atomic properties from an input DataFrame Parameters ---------- df : pandas.DataFrame A pandas DataFrame with atomic properties that will be used to set the properties on the current list of atoms See Also -------- :func:`parmed.utils.pandautils.load_dataframe` for full documentation r)load_dataframe)rr)rdfrrFrFrGrs zStructure.load_dataframecCsH|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|j jpF|j jpF|j jpF|j jpF|j jpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|jjpF|j jpF|jjpF|jjpF|jjpF|jjpF|j jS)zB Determines if any of the topology has changed for this structure )!rechangedrfrgrnrirkrlrmrhrorprqrrrsrtrurvrwrxryrzr{r|rr}rrrrrrjr~)rrFrFrG is_changeds&    zStructure.is_changedcCsd|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j _d|j _d|j _d|j _d|j _d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_d|j_dS)z< Toggles all lists so that they do not indicate any changes FN) rerrfrgrhrirkrlrmrnrorprqrrrsrtrurvrwrxryrzr{r|r}rrrrrrr)rrFrFrGunchanges>zStructure.unchangecCs||||||||||| | | | | |dS)z Looks through all of the topological lists and gets rid of terms in which at least one of the atoms is None or has an `idx` attribute set to -1 (indicating that it has been removed from the `atoms` atom list) N)_prune_empty_bonds_prune_empty_angles_prune_empty_dihedrals_prune_empty_rb_torsions_prune_empty_ureys_prune_empty_impropers_prune_empty_cmaps_prune_empty_trigonal_angles_prune_empty_out_of_plane_bends_prune_empty_pi_torsions_prune_empty_stretch_bends_prune_empty_torsion_torsions_prune_empty_chiral_frames_prune_empty_multipole_frames_prune_empty_adjusts)rrFrFrGprune_empty_termsszStructure.prune_empty_termscCst}g}x|jD]}|jrq|j|jjks:|j|jjkrBd|_q|jj|jjf|kr^d|_qt|j t r|j j dkr| |q| |jj|jjf| |jj|jjfqWx&|D]}|jj|jjf|krd|_qWdS)a2 Nonbonded exclusions and exceptions have the following priority: bond -> angle -> dihedral Since bonds and angles are completely excluded, any ring systems in which two atoms are attached by a bond or angle as well as a dihedral should be completely excluded as the bond and angle exclusion rules take precedence. If a Bond or Angle was _added_ to the structure between a pair of atoms previously connected only by a dihedral term, it's possible that those two atoms have both an exclusion *and* an exception defined. The result of this scenario is that sander and pmemd will happily compute an energy, _including_ the 1-4 nonbonded terms between atoms now connected by a bond or an Angle. OpenMM, on the other hand, will complain about an exception specified multiple times. This method scans through all of the dihedrals in which `ignore_end` is `False` and turns it to `True` if the two end atoms are in the bond or angle partners arrays TrN)setrirrrrQZangle_partnersrr@rrperr=add)rZset14Zdeferred_dihedralsZdihedralrFrFrGupdate_dihedral_exclusionss"   z$Structure.update_dihedral_exclusionscs6ddlm}t||r2|j|k r(td|nht|trL|||nNy t|Wn&tk r~tdt|j YnXt t |j krt dt ddtD}xt|D] }|j |=qW||j||jdk r2trtfd d|jD|_n|jddtdkf|_dS) a  Deletes a subset of the atoms corresponding to an atom-based selection. Parameters ---------- selection : :class:`AmberMask`, ``str``, or ``iterable`` This is the selection of atoms that will be deleted from this structure. If it is a string, it will be interpreted as an AmberMask. If it is an AmberMask, it will be converted to a selection of atoms. If it is an iterable, it must be the same length as the `atoms` list. r) AmberMaskz(passed mask does not belong to Structurez#Selection not a supported type [%s]zSelection iterable wrong lengthcSsg|]\}}|r|qSrFrF)rrrrFrFrGr5sz#Structure.strip..Ncs"g|]}fddt|DqS)cs,g|]$\}\}}}|dkr|||gqS)rrF)rrxyz)selrFrGr@s z.Structure.strip...)r)rZcrd)rrFrGr@s) parmed.amberrr@parmrA Selectionr3rrrrrersortedrreversedrrfZprunerrr-nparray)r selectionrZatomlistrrF)rrGstrips4          zStructure.stripc Gsddlm}t|}x|D]}||qWt}t}t}x|jD]}t||} | dkrt||j||qFdd|jD} x|jD]} | j | j kr|| q| j | j j | _ xX| j | j j D]F} | j | kr| | j | j kr| | | j j kr|j t| | | j qWqWqFWxZt|jddD]B\} }t||} |js0|jrp|j|j| djkrpq0t|j| d|}| dkr|dkrq0| dkrH|jdkrq0x:|j| djD]}|j |jj krPqWtdq0xP|jD]F} t| j |j f}t| ||kr| |j kr<|j t| |PqWq0| jdkrXq0|dkrx&|jD]}|j | jj krjPqjWq0xZ|j| djD]F} t| j |j f}t| ||kr| |j kr|j t| |PqWq0|jdkrq0x&|jD]}|j | jj krPqWq0x0|j| djD]}|j |jj kr4Pq4Wq0||j kr0|j t||q0Wx|jD]}t|j d kstj|j sqtj|j jd krqxB|jD]8} | j d krt| j d kr|| || PqWqW|jdkrdSt !t"t#}t$|||}x|D]}xR||j%D]D}t|j |j f}t|||krP||j krP|j t||qPW|j|ks@||krq@x^|jjD]R}||krΐqt|j |j f}t|||kr||j kr|j t||qWq@WdS) a Assigns bonds to all atoms based on the provided residue template libraries. Atoms whose names are *not* in the templates, as well as those residues for whom no template is found, is assigned to bonds based on distances. Parameters ---------- reslibs : dict{str: ResidueTemplate} Any number of residue template libraries. By default, assign_bonds knows about the standard amino acid, RNA, and DNA residues. r)StandardBiomolecularResiduesNcSsi|] }||jqSrF)rT)rrrFrFrG gsz*Structure.assign_bonds..rrzCould not find the head atom of the next template! Bond pattern may be wrong, which could lead to extra TER cards in a PDB filer]ZCYSZSGr\)&Zparmed.modellerrrupdaterrf_res_in_templlibrerrTmaprMrQrgr=rrrrheadLOGGERZwarningr rtailrrAminoAcidResiduehasgetabbr coordinatesmathsqrtmaxvaluesr.r)rZreslibsrZ all_residueslibZunassigned_residuesZunassigned_atomsZ cysteine_sgrZtemplZresatomsrZbprZntemplr Zmaxdistr ZmindistZpairsrNpartnerrFrFrG assign_bondsIs         $  (                      zStructure.assign_bondscOs.|jdkrtdddlm}||f||S)aUse nglview for visualization. This only works with Jupyter notebook and require to install `nglview` Examples -------- >>> import parmed as pmd >>> parm = pmd.download_PDB('1tsu') >>> parm.visualize() Parameters ---------- args and kwargs : positional and keyword arguments given to nglview, optional Nzcoordinates must not be Noner) show_parmed)rrZnglviewr)rrDkwargsrrFrFrG visualizes  zStructure.visualizec s^t|tr|j|S||}|dkr0t|St|tr>|St|}|dkrXt|S|dgx0tdt|D]} |d||qrWddt |Dt|x^t |jD]P\}}||sq|j }|j dkr|j}n|j }t||j||j|j|jqWfdd}|jj|j|jdd g|jj|j|jdd d g|jj|j|jdd d d d d g|jj|j|jdd d d d d g|jj|j|jdd g|jj |j|j dd d d g|j!j"|j!|j"dd d d dg|j#j$|j#|j$dd d d g|j%j&|j%|j&dd d d g|j'j(|j'|j(dd d d ddg|j)j*|j)|j*dd d g|j+j,|j+|j,dd d d dg|j-g|j-gdd dg|j.g|j.gdddddg|j/j0|j/|j0dd g|j1g|j1gdd g|j2g|j2gdd g|j3g|j3gdddg|j4_4|j5_5|j6_6S)aY Allows extracting a single atom from the structure or a slice of atoms as a new Structure instance. The following syntaxes are allowed: - struct[str] : str is interpreted as an Amber selection mask - struct[[sel,[sel,]],sel] : sel can be a list of indices (or strings for chain selection) an integer, or a slice Parameters ---------- selection : str, or slice|iter|str, slice|iter|int, slice|iter|int Which atoms to select for the new structure Returns ------- struct : :class:`Structure` If more than one atom is selected, the resulting return value is a new structure containing all selected atoms and all valence terms (and types) in which all involved atoms are present or atom : :class:`Atom` If the selection is a single integer, a tuple of two integers, or a tuple of a string and two integers. This is equivalent to selecting a particular atom, a particular atom from a particular residue, or a particular atom from a particular residue from a particular chain, respectively. All other selections return a Structure instance, even if that selection happens to only select a single atom. Notes ----- When selecting more than 1 atom, this is a costly operation. It is currently implemented by making a full copy of the object and then stripping the unused atoms. Raises ------ ``TypeError`` : if selection (or any part of the selection) is not an allowed type ``ValueError`` : if the selection is a boolean-like list and its length is not the same as the number of atoms in the system NrrcSsg|]\}}||qSrFrF)rrrrFrFrGr?sz)Structure.__getitem__..csJdd|D}dd|D}x|D]fdd|D}fdd|D}t|sTq"t} tdrxjtkrxt| d<n*|rjdk r|jj| d<d|jj<x6t|D]*\} } t| trj | jd || <qW| t|| td r"j |d _ q"Wx(t ||D]\} } | r| | qWt|d rF| dS) z Copies the valence terms from one list to another; oval=Other VALence; otyp=Other TYPe; sval=Self VALence; styp=Self TYPe; attrlist=ATTRibute LIST (atom1, atom2, ...) cSsg|] }t|qSrF)r)rrrFrFrGrOszEStructure.__getitem__..copy_valence_terms..cSsg|]}dqS)FrF)rrrFrFrGrPscsg|]}t|qSrF)getattr)rattr)valrFrGrRscs g|]}t|tr|jqSrF)r@rr)rr)scanrFrGrTsrNTrrrr)alldictrrrrrr@rrer=rr5r)ovalotypsvalstypattrlistotypcpZ used_typesatsindiceskwsrrZusedr)r struct)rrGcopy_valence_termsJs0      z1Structure.__getitem__..copy_valence_termsrrrrrrrrrrNrrrrrr)7r@r1re_get_selection_arrayrrrr4rr=r5rrrrrrrTrrrrgryrhrzrir{rjr~rkr|rlr}rmrrnrrorrprrqrrrrrsrtrurrwrvrxrrr)rrsumselrrNrZnumr-rF)r r,rG __getitem__s-         $        zStructure.__getitem__csddlm}t|tr,|||}|}nt|trvdd|jD}x&ttt |j|D] }d||<q^W|}nft|t r*t |dkr*t |dkr|\}}t|t rt|t r|j ||Sd}nJt |d kr|\} }}t t} x|j D]} | | j| qWt| trnt|t rnt|t rnt| } y| | ||Stk rltd | YnXt| tr| | krt| gnd Sn^t| tr|j djg} x,|j D]"} | j| d kr| | jqWt| | nt| xD]}|| krPqWd Sd }t|tr@tttt |j |n,t|tsXt|t rdt|gnt|t|trtttt |j|n,t|tst|t rt|gnt||rz.r)r\r]r\Fr]zNo chain %s in StructureNrTcsPg|]H}|jjkoJ|jjks,|jjkoJ|jkpJ|j|jdjkqS)r)rrrTr)rr)atomsetchainsetressetrFrGrscsDg|]<}|jks(|j|jdjko>|jjkp>|jjkqS)r)rTrr)rr)r1r3rFrGrscSsg|]}dqS)rrF)rrNrFrFrGrszSelection iterable is too longzSelected atom out of rangecSsg|]}tt|qSrF)intbool)rrrFrFrGrs)rrr@r3rslicererr4rtupler1rfrr'rr=r"KeyError IndexErrorrr2rrr)rrrmaskrrZresselZatomselZ has_chainZchainselZchainmaprZchainsrrZ chain_namerrrF)r1r2r3rGr.s                         zStructure._get_selection_arraycCst|S)z Returns an indexable object that can be indexed like a standard Structure, but returns a *view* rather than a copy See Also -------- Structure.__getitem__ )_StructureViewerCreator)rrFrFrGviews zStructure.viewcslt|dd|jD}t|}g}g}t}ddt|D}x"|jD]}||jd|qHWxt|D]|}tdd|D} t| dkr$|dj } t dd| D} t d d| D} t d d| D} t d d| D}| j t| | | | |f}||kr$||| qpd }xt |D]\}}t|jt|kr2xt|j|D]\}}d |j |j fks~td|j j |j j krP|j |j krPd|jd|jkrPd|jd|jkrPd|jd|jkr\Pq\W|| d}Pq2W|spt| dkr(t|||<|fdd|jD}|||tgqpWtt||S)a Split the current Structure into separate Structure instances for each unique molecule. A molecule is defined as all atoms connected by a graph of covalent bonds. Returns ------- [structs, counts] : list of (:class:`Structure`, list) tuples List of all molecules in the order that they appear first in the parent structure accompanied by the list of the molecule numbers in which that molecule appears in the Structure cSsg|] }|jqSrF)marked)rrNrFrFrGr2sz#Structure.split..cSsg|]}gqSrFrF)rrrFrFrGr8srcss|]}|jjVqdS)N)rr)rrNrFrFrG =sz"Structure.split..rcss|] }|jVqdS)N)rT)rrrFrFrGr>Escss|]}d|jVqdS)z%.6fN)charge)rrrFrFrGr>Fscss|]}d|jVqdS)z%.6fN)rmin)rrrFrFrGr>Gscss|]}d|jVqdS)z%.6fN)epsilon)rrrFrFrGr>HsFNzResidues must all be setz%.6fTcsg|]}|jdkqS)r)r=)rrN)rrFrGr_s)r/rerr"r4r=r=rrrr7rTrrr5AssertionErrorr?r@rAr)rZmollistZnmolZstructsZcountsZ res_moleculesZmolatomsrNrZinvolved_residuesrnamesZchargesZrminsZepsilonsZ oneres_keyZ is_duplicatejr,a1a2ZmolrF)rrGsplit$s^      zStructure.splitNcKsddlm}m}m}m}ddddddddd d d d d d d dd} t|dsftj|rv|svt d|n|dkrvt dd} x6|j D],} t | j tr| jtk rt| j| _ qd} qWz|dk r|}nXtj|\} } | dkrtj| d} y | | }Wn"tk rtd|YnX|dkr<|j|f|nf|dkrX|j|f|nJ|dkrx|jj||f|n*|dkr|j|}|j|f|n|d kr|jj||f|n|d kr|jj||f|n|d kr|jj||fddi|n|dkr0|j|}|j|f|nr|d krP|j j||f|nR|dkrF|j!s|j"s|j#s|j$s|j%s|j&s|j'r|j(|}|j)|f|n|j*s|j+s|j,r|j-|}|j)|f|nby|j.|}WnBt/k r4}z"dt|kr"|j-|}nWdd}~XYnX|j)|f|n\|dkr|j0fdt1|j i|}|j2|_2|j3|_4|j5|_5|j||dkdn td|Wd| rx|j D]} t6| j| _ qWXdS)a Saves the current Structure in the requested file format. Supported formats can be specified explicitly or determined by file-name extension. The following formats are supported, with the recognized suffix and ``format`` keyword shown in parentheses: - PDB (.pdb, pdb) - PDBx/mmCIF (.cif, cif) - PQR (.pqr, pqr) - Amber topology file (.prmtop/.parm7, amber) - CHARMM PSF file (.psf, psf) - CHARMM coordinate file (.crd, charmmcrd) - Gromacs topology file (.top, gromacs) - Gromacs GRO file (.gro, gro) - Mol2 file (.mol2, mol2) - Mol3 file (.mol3, mol3) - Amber ASCII restart (.rst7/.inpcrd/.restrt, rst7) - Amber NetCDF restart (.ncrst, ncrst) Parameters ---------- fname : str or file-like object Name of the file or file-like object to save. If ``format`` is ``None`` (see below), the file type will be determined based on the filename extension. If ``fname`` is file-like object, ``format`` must be provided. If the type cannot be determined, a ValueError is raised. format : str, optional The case-insensitive keyword specifying what type of file ``fname`` should be saved as. If ``None`` (default), the file type will be determined from filename extension of ``fname`` overwrite : bool, optional If True, allow the target file to be overwritten. Otherwise, an IOError is raised if the file exists. Default is False kwargs : keyword-arguments Remaining arguments are passed on to the file writing routines that are called by this function Raises ------ ValueError if either filename extension or ``format`` are not recognized TypeError if the structure cannot be converted to the desired format for whatever reason IOError if the file cannot be written either because it exists and ``overwrite`` is ``False``, the filesystem is read-only, or write permissions are not granted for the user r)ambercharmmformatsgromacsZPDBZPQRZCIFZAMBERZPSFZGROMACSZGROZMOL2ZMOL3Z CHARMMCRDRST7NCRST)z.pdbz.pqrz.cifz.pdbxz.parm7z.prmtopz.psfz.topz.groz.mol2z.mol3z.crdz.rst7z.inpcrdz.restrtz.ncrstwritez%s exists; not overwritingNz7Must provide supported format if using file-like objectTF)z.bz2z.gzrz#Could not determine file type of %sZmol3zCannot translate exceptions)rLrMr)ZnetcdfzNo file type matching %s)7ZparmedrHrIrJrKrospathexistsIOError RuntimeErrorrer@rr1 atom_typer*struppersplitextr8rZ write_pdbZ write_cifZPQRFilerNZ CharmmPsfFileZfrom_structureZ write_psfZGromacsGroFileZMol2FileZGromacsTopologyFileZ CharmmCrdFilernrorrrprqrsrtZ AmoebaParmZ write_parmrkrlrmZ ChamberParmZ AmberParmrAZRst7rr velocitiesZvelsrr4)rfnameformatZ overwriterrHrIrJrKZextmapZall_intsrNbaseZextreZrst7rFrFrGsavefs/                    zStructure.savecCs|js dStdd|jDr"dStdd|jDr:dSt}t}t}xt|jD]\}}|j|jkr|j|j |j |jf}n|j|j |j |jf}||kr|| |j | |qXt ||<}| |j | |||_ qXW||_x&t|D]}|j||j|=qWdS)a5 Joins multi-term torsions into a single term and makes all of the parameters DihedralTypeList instances. If any dihedrals are *already* DihedralTypeList instances, or any are not parametrized, or there are no dihedral_types, this method returns without doing anything Ncss|]}t|tVqdS)N)r@r)rrrFrFrGr>sz+Structure.join_dihedrals..css|]}|jdkVqdS)N)r)rrrFrFrGr>s)r{anyrirr"r'rrrrrr=rrrdelete)rZdihedrals_to_deleteZdihedrals_processedZnew_dihedral_typesrrkeyZdtlrFrFrGjoin_dihedralss0     zStructure.join_dihedralscCs|jS)N)r)rrFrFrGrszStructure.combining_rulecCs|dkrtd||_dS)N)rd geometricz/combining_rule must be 'lorentz' or 'geometric')rr)rthingrFrFrGr!sc Csn|s&y|jStk r"Yq6Xn||t|_}|}y|jdj }d}d}Wnt k rz|jSXxt |j D]\}}|j |k r||j j kr|j j }|}|j }||j j|}ytjj|j}Wntk rd}YnX||j||qWt| } x,|jD]"} || | jj| | jjq$W|jdk rj|tt|j|S)a" The OpenMM Topology object. Cached when possible, but any changes to the Structure instance results in the topology being deleted and rebuilt Notes ----- This function calls ``prune_empty_terms`` if any topology lists have changed. rN) rZ _topologyAttributeErrorrrr7ZTopologyZaddChainrfrr9rrerZ addResiduerTelementZElementZgetByAtomicNumberrMr8ZaddAtomrrgaddBondrrrrZsetPeriodicBoxVectorsr9r ) rtoprZ last_chainZ last_residueZlast_omm_residuerrNelemrebondrFrFrGtopology)sF       " zStructure.topologycCs0ydd|jDtjStk r*dSXdS)a, A list of 3-element Quantity tuples of dimension length representing the atomic positions for every atom in the system. If set with unitless numbers, those numbers are assumed to be in angstroms. If any atoms do not have coordinates, this is simply ``None``. cSsg|]}t|j|j|jqSrF)r6xxxyxz)rrrFrFrGrnsz'Structure.positions..N)rer:r?rd)rrFrFrG positionseszStructure.positionscCst|r|tj}t|t|jkrVxt|jD]\}}||\|_|_|_ q4Wn\t|dt|jkrxDt|jD],\}}|d}|||d\|_|_|_ qxWnt ddS)a A list of 3-element Quantity tuples of dimension length representing the atomic positions for every atom in the system. If set with unitless numbers, those numbers are assumed to be in angstroms. Raises ------ ValueError if the positions are not either a (natom, 3)-shape iterable of iterables or a (natom*3)-length iterable. r]zWrong shape for position arrayN) r:r;r>r?rrerrkrlrmr)rvaluerrNZi3rFrFrGrnrs  "cCs4ydd|jD}Wntk r(dSXt|S)NcSsg|]}|j|j|jgqSrF)rkrlrm)rrrFrFrGrsz)Structure.coordinates..)rerdrr)rcoordsrFrFrGrs zStructure.coordinatesc Cs|dkrFd|_x|jD]*}y|`|`|`Wqtk r>YqXqWnt|r\|tj }t |}t j |t j ddd}|dt|jdf}t|dkrx,t|j|dD]\}}|\|_|_|_qW||_n>x6|jD],}y|`|`|`Wqtk rYqXqWd|_dS)z? Setting coordinates will also set xx, xy, and xz on the atoms NFT)dtypersubokrr]r)rrerkrlrmrdr:r;r>r?rrrfloat64reshaperr5)rrorNrprxyzrFrFrGrs,       r!cCsydd|jD}Wntk r,d}Yn Xt|}|jdk r|dkrRd|_n@|j|jjddkrpd|_n"t||jdtkrd|_|dkr|jdk r|jS|dk r| dt |jdfSdS|jdkr|dkr|dk r|St d|j|S) a In some cases, multiple conformations may be stored in the Structure. This function retrieves a particular frame's coordinates Parameters ---------- frame : int or 'all', optional The frame number whose coordinates should be retrieved. Default is 'all' Returns ------- coords : np.ndarray, shape([#,] natom, 3) or None If frame is 'all', all coordinates are returned with shape (#, natom, 3). Otherwise the requested frame is returned with shape (natom, 3). If no coordinates exist and 'all' is requested, None is returned Raises ------ IndexError if there are fewer than ``frame`` coordinates cSsg|]}|j|j|jgqSrF)rkrlrm)rrrFrFrGrsz-Structure.get_coordinates..Nrrr!r]zNo coordinate frames present) rerdrrrshaperrr rtrr9)rframerprFrFrGget_coordinatess.     zStructure.get_coordinatescCs|jdkrdS|jdS)Nr)r)rrFrFrGrs z Structure.boxcCs|dkrd|_nlt|tjr"|}n tt|}tj|tjddd}|jdkrpt |jdksh|jddkrpt d| d |_dS) NFT)rqrrr)rHr\rrHz Box information must be 6 floats)rrH) rr@rZndarrayrCrrrsrvrrrt)rrorrFrFrGrs   cCs:|jdkr|dkrtddS|dkr,|jS|j|SdS)a In some cases, multiple conformations may be stored in the Structure. This function retrieves a particular frame's unit cell (box) dimensions Parameters ---------- frame : int or 'all', optional The frame number whose unit cell should be retrieved. Default is 'all' Returns ------- box : np.ndarray, shape([#,] 6) or None If frame is 'all', all unit cells are returned with shape (#, 6). Otherwise the requested frame is returned with shape (6,). If no unit cell exist and 'all' is requested, None is returned Raises ------ IndexError if there are fewer than ``frame`` unit cell dimensions Nr!zNo unit cell frames present)rr9)rrwrFrFrGget_boxs zStructure.get_boxcCs0ytdd|jDStk r*dSXdS)z A (natom, 3)-shape numpy array with atomic velocities for every atom in the system (in units of angstrom/picosecond), or None if there are no velocities cSsg|]}|j|j|jgqSrF)vxvyvz)rrrFrFrGr$sz(Structure.velocities..N)rrrerd)rrFrFrGrXszStructure.velocitiesc Cst|r|tjtj}|dkr\x|jD]*}y|`|`|`Wq,t k rTYq,Xq,WnNt j |dd dt |jdf}x,t|j|dD]\}}|\|_|_|_qWdS)z A list of 3-element Quantity tuples of dimension length representing the atomic velocities for every atom in the system NF)rrr]r)r:r;r>r?Z picosecondsrerzr{r|rdrrrtrr5)rrorNrurFrFrGrX(s    cCsht}x*|jD] }|jtkrq|j|t|j<qWx0t|D]$\}}x|jD]}||krLdSqLWqr:r?Zdegreesrrrsr)rrorrrABGrFrFrGr_s      g @gTgfffffr@g?gS@gMb@?cCs |jrtd|dkrtj}t}|tjtjtjfkrJ|j dkrJt ddd|j D}|dk rt |rx|t j}|dkrt dxf|j D]\}|jdkrd}x|jD]}|jdkr|}PqW|dk r|||j<||j||j8<qWx|D]}||qW|||||r0||td ||||||td |||||td ||||td |||td ||| td|||!td|||"td|||#td|||$td|||%td|||&td|||'td|dk rd}nd}|||(||||||dk rtd|||)|||| | ||| | | r|*t+|j dk r|j,t-|j.|/||S)a$ Construct an OpenMM System representing the topology described by the prmtop file. Parameters ---------- nonbondedMethod : cutoff method This is the cutoff method. It can be either the NoCutoff, CutoffNonPeriodic, CutoffPeriodic, PME, or Ewald objects from the simtk.openmm.app namespace nonbondedCutoff : float or distance Quantity The nonbonded cutoff must be either a floating point number (interpreted as nanometers) or a Quantity with attached units. This is ignored if nonbondedMethod is NoCutoff. switchDistance : float or distance Quantity The distance at which the switching function is turned on for van der Waals interactions. This is ignored when no cutoff is used, and no switch is used if switchDistance is 0, negative, or greater than the cutoff constraints : None, app.HBonds, app.HAngles, or app.AllBonds Which type of constraints to add to the system (e.g., SHAKE). None means no bonds are constrained. HBonds means bonds with hydrogen are constrained rigidWater : bool=True If True, water is kept rigid regardless of the value of constraints. A value of False is ignored if constraints is not None. implicitSolvent : None, app.HCT, app.OBC1, app.OBC2, app.GBn, app.GBn2 The Generalized Born implicit solvent model to use. implicitSolventKappa : float or 1/distance Quantity = None This is the Debye kappa property related to modeling saltwater conditions in GB. It should have units of 1/distance (1/nanometers is assumed if no units present). A value of None means that kappa will be calculated from implicitSolventSaltConc (below) implicitSolventSaltConc : float or amount/volume Quantity=0 moles/liter If implicitSolventKappa is None, the kappa will be computed from the salt concentration. It should have units compatible with mol/L temperature : float or temperature Quantity = 298.15 kelvin This is only used to compute kappa from implicitSolventSaltConc soluteDielectric : float=1.0 The dielectric constant of the protein interior used in GB solventDielectric : float=78.5 The dielectric constant of the water used in GB useSASA : bool=False If True, use the ACE non-polar solvation model. Otherwise, use no SASA-based nonpolar solvation model. removeCMMotion : bool=True If True, the center-of-mass motion will be removed periodically during the simulation. If False, it will not. hydrogenMass : float or mass quantity = None If not None, hydrogen masses will be changed to this mass and the difference subtracted from the attached heavy atom (hydrogen mass repartitioning) ewaldErrorTolerance : float=0.0005 When using PME or Ewald, the Ewald parameters will be calculated from this value flexibleConstraints : bool=True If False, the energies and forces from the constrained degrees of freedom will NOT be computed. If True, they will (but those degrees of freedom will *still* be constrained). verbose : bool=False If True, the progress of this subroutine will be printed to stdout splitDihedrals : bool=False If True, the dihedrals will be split into two forces -- proper and impropers. This is primarily useful for debugging torsion parameter assignments. Notes ----- This function calls prune_empty_terms if any Topology lists have changed z'Cannot createSystem: unknown functionalNz(No periodic boundary conditions detectedcSsg|] }|jqSrF)mass)rrNrFrFrGrsz*Structure.createSystem..rzHydrogen mass must be positiverzAdding bonds...zAdding angles...zAdding dihedrals...z%Adding Ryckaert-Bellemans torsions...zAdding Urey-Bradleys...zAdding improper torsions...zAdding CMAP torsions...zAdding trigonal angle terms...zAdding out-of-plane bends...zAdding pi-torsions...zAdding stretch-bends...zAdding torsion-torsions...zAdding Nonbonded force...g?gS@zAdding GB force...)0rr r7NoCutoffmmZSystemCutoffPeriodicPMEEwaldrrrer:r;r>ZdaltonrerQrr addParticleomm_add_constraintsrrrr info_add_force_to_systemomm_bond_forceomm_angle_forceomm_dihedral_forceomm_rb_torsion_forceomm_urey_bradley_forceomm_improper_forceomm_cmap_forceomm_trigonal_angle_forceomm_out_of_plane_bend_forceomm_pi_torsion_forceomm_stretch_bend_forceomm_torsion_torsion_forceomm_nonbonded_forceomm_gbsa_forceaddForceZCMMotionRemoverZsetDefaultPeriodicBoxVectorsr9romm_set_virtual_sites)rnonbondedMethodnonbondedCutoffswitchDistance constraints rigidWaterimplicitSolventimplicitSolventKappaimplicitSolventSaltConc temperaturesoluteDielectricsolventDielectricuseSASAZremoveCMMotionZ hydrogenMassewaldErrorToleranceflexibleConstraintsverboseZsplitDihedralssystemZmassesrNZ heavy_atomrFrZrf_dielcrFrFrG createSystempsY                                zStructure.createSystemcCs|dkr|sdS|dtjtjtjfkr2td|tjtj}t }t }t |}|tjksh|tjkrx|j D]:}t |j trqpt |jtrqp|t|j j|jjfqpWnj|tjkrx\|j D]R}t |j trqt |jtrq|j jdks|jjdkr|t|j j|jjfqW|rx^|j D]T}t |j tr>q(t |jtrPq(||j jjr(|t|j j|jjfq(WxH|j D]>}t|j j|jjf|kr||j j|jj|jj|qW|tjkrVx~|jD]t} d} | j jdkr| d7} | jjdkr| d7} | dks4| dkr| jjdkr|| j j| jj| jjfqW|rx<|jD]2} || j jjrd|| j j| jj| jjfqdWx|jD]} | j j| jj| jjf|krt| j j| jjf|krqd} } xX| jj D]L}|| kr| j |kr|jj|} n"|| kr| j|kr|jj|} qW| dks| dkr^qt| jjt} t| | | | d| | | }|| j j| jj|qWdS)a Adds constraints to a given system Parameters ---------- system : mm.System The OpenMM system for which constraints should be added constraints : None, app.HBonds, app.AllBonds, or app.HAngles Which kind of constraints should be used rigidWater : bool If True, water bonds are constrained regardless of whether constrains is None Nz$Unrecognized constraints option (%s)rrr\rP)r7HBondsAllBondsHAnglesrr:angstromconversion_factor_to nanometerr_settlerrgr@rrrr frozensetrrerZ addConstraintrreqrhrrZcostheteqr r)rrrr length_convZconstraint_bond_setZconstraint_angle_setis_waterriangleZnumHl1l2ZcostrrFrFrGrs            $  "   &zStructure.omm_add_constraintsc Cs|t|jkrtdx|jD]}t|ts2q"|j}|}|}t|t r|\}}|\}} | |j t |j |j || q"t|tr|\}}} |\}} } | |j t |j |j | j || | q"t|tr"|\}}} |\}} } | |j t |j |j | j || | q"WdS)an Sets the virtual sites in a given OpenMM `System` object from the extra points defined in this system Parameters ---------- system : mm.System The system for which the virtual sites will be set. All particles must have already been added to this System before calling this method z.OpenMM System does not correspond to StructureN)getNumParticlesrrerr@rZ frame_typeZ get_weightsZ get_atomsr)ZsetVirtualSiterrZTwoParticleAverageSiter%ZThreeParticleAverageSiter!ZOutOfPlaneSite) rrrNrZweightsZrefatomsrErFZw1Zw2Za3Zw3rFrFrGrts,        "   zStructure.omm_set_virtual_sitesc Cs|s|tjtjfks|jsdStjtj}tjtj d}tj tj d}||}t |j drt |j drt}||j jd|||j jd|dnt}||j|rt|} ndd|jD} x|jD]} d | jj| jjfkr|s|tjkrq|s(| | jjjr(q| jdkr.rz Cannot find necessary parametersr)$r7rrrgr:r?r nanometerskilocalorie_per_molerkilojoule_per_molerrryrZAmoebaBondForceZsetAmoebaGlobalBondCubicrZsetAmoebaGlobalBondQuarticHarmonicBondForce setForceGroupBOND_FORCE_GROUPrrfrrerrrrrr rfrkZ getNumBonds) rrrrr_ambfrc_ommfrcfrc_convforcerrirFrFrGrs:      zStructure.omm_bond_forcec Cs<|js dStjtj}t|jdrzt|jdrz|jj}t }| |d| |d| |d| |dnt}||jx|jD]}|jjdk|jjdk}|tjkr|d ks|dkr|jjd kr|sq|jdkrtd ||jj|jj|jj|jjtd |jj|qW|d kr8dS|S) a Creates an OpenMM HarmonicAngleForce object (or AmoebaAngleForce if the angles are for an Amoeba-parametrized system) Parameters ---------- constraints : None, app.HBonds, app.AllBonds, or app.HAngles The types of constraints that are on the system. If flexibleConstraints is False, then the constrained bonds will not be added to the resulting Force flexibleConstraints : bool=True If True, all bonds are added to the force regardless of `constraints` Returns ------- force HarmonicAngleForce (or AmoebaAngleForce if this is an Amoeba system), or None if there are no angles to add Nr<rr]r^r_rHrr\rPzCannot find angle parametersr)rhr: kilocaloriesr kilojoulesrrzrrZAmoebaAngleForceZsetAmoebaGlobalAngleCubicZsetAmoebaGlobalAngleQuarticZsetAmoebaGlobalAnglePenticZsetAmoebaGlobalAngleSexticZHarmonicAngleForcerANGLE_FORCE_GROUPrrerr7rrrr addAnglerrr rZ getNumAngles)rrrrrrrZnum_hrFrFrGrs2      zStructure.omm_angle_forcec Cs6|js dStjtj}t}t}||j||j x|jD]}|j dkr^t d|rn|j rn|}n|}t |j trx|j D]<}||jj|jj|jj|jjt|j|jt|j|qWqH||jj|jj|jj|jjt|j j|j jt|j j|qHW|dkr|S|dkr.|S||fS)a- Creates the OpenMM PeriodicTorsionForce modeling dihedrals Parameters ---------- split : bool, optional, default=False If True, separate PeriodicTorsionForce instances with the propers in the first and impropers in the second return item. If no impropers or propers are present, the instances with zero terms are not returned. Returns ------- PeriodicTorsionForce[, PeriodicTorsionForce] Or None if no torsions are present in this system NzCannot find torsion parametersr)rir:rrrrZPeriodicTorsionForcerDIHEDRAL_FORCE_GROUPIMPROPER_FORCE_GROUPrr rr@r addTorsionrrrrrr4rZphaser phi_kZgetNumTorsions)rrGrZproperrtorrrrFrFrGr s2       "zStructure.omm_dihedral_forcecCs|js dStjtj}t}||jx||jD]r}|j dkrJt d| |j j |jj |jj |jj |j j||j j||j j||j j||j j||j j| q4W|S)z Creates the OpenMM RBTorsionForce for Ryckaert-Bellemans torsions Returns ------- RBTorsionForce Or None if no torsions are present in this system Nz"Cannot find R-B torsion parameters)rjr:rrrrZRBTorsionForcerRB_TORSION_FORCE_GROUPrr rrrrrrZc0c1c2Zc3Zc4Zc5)rZconvrrrFrFrGr= s    &zStructure.omm_rb_torsion_forcecCs|js dStjtj}tjtjd}tjtjd}||}t }| |j xL|jD]B}|j dkrttd||jj|jj|j j|d|j j|q^W|S)z Creates the OpenMM Urey-Bradley force Returns ------- HarmonicBondForce Or None, if no urey-bradleys are present Nr\z#Cannot find urey-bradley parameters)rkr:r?rrrrrrrrrUREY_BRADLEY_FORCE_GROUPrr rfrrrrr)rrrrrrZureyrFrFrGrT s     z Structure.omm_urey_bradley_forcec Cs|js dStjtj}d}|d7}|d7}|dtj7}t|}| d| d| |j xV|jD]L}|j dkrt d||jj|jj|jj|jj|j j||j jtfqlW|S) z Creates the OpenMM improper torsion force (quadratic bias) Returns ------- CustomTorsionForce With the formula k*(phi-phi0)^2, or None if there are no impropers Nzk*dtheta_torus^2;z8dtheta_torus = dtheta - floor(dtheta/(2*pi)+0.5)*(2*pi);zdtheta = theta - theta0;zpi = %f;rZtheta0z'Cannot find improper torsion parameters)rlr:rrrrpirZCustomTorsionForceZaddPerTorsionParameterrrrr rrrrrrZpsi_kZpsi_eqr )rrZenergy_functionrimprFrFrGrm s$       zStructure.omm_improper_forcec s|js dStjtjt}||jg}t }xv|jD]l}|j dkrTt dt |j |kr>|j }| t ||jj}||jfdd|D}||t |<q>WxR|jD]H}||t |j |jj|jj|jj|jj|jj|jj|jj|jj qW|S)z Creates the OpenMM CMAP torsion force Returns ------- CMAPTorsionForce Or None, if no CMAP terms are present Nz#Cannot find CMAP torsion parameterscsg|] }|qSrFrF)rr)rrFrGr sz,Structure.omm_cmap_force..)rmr:rrrrZCMAPTorsionForcerCMAP_FORCE_GROUPr"rr idr=gridZ switch_rangeTZaddMap resolutionrrrrrrr)rrZcmap_type_listZcmap_mapcmaprrrrF)rrGr s*      zStructure.omm_cmap_forcec s"|js dStjtj}tjtj}t |j t |rP| tj }|dksb|tjkrrtjjn|tjkrtjj|n|tjkrtjj|nj|tjkrtjj||n<|tjkrtjj||n td|||dd}x0|jD]&} | j| j|t| j|qBWfddx |jD]} | | d|jqWd|}|jdkrd d } n|jd krd d } |j sؐx|j!|j"D]} | j#rqt$| j%t&r|d} } d}xJ| dks(| dkr\|t'| j%kr\| j%|j(} | j%|j)} |d 7}qW| dksr| dkrtdn| j%j(} | j%j)} y"| j*j+j,t-| j.j+\}}}}WnZt/t0fk rt| j*j1| j.j1}t23||| }| | j*j4| j.j4}||9}YnX||| }|||}| j*j| j.j| }5| j*j6| j.j6|||dxz| j*j7D]n}t|j1| j.j1}t23||| }| |j4| j.j4}||9}|j| j.j| }5|j6| j.j6|||dq\Wxz| j.j7D]n}t|j1| j*j1}t23||| }| |j4| j*j4}||9}|j| j*j| }5|j6| j*j6|||dqWx| j*j7D]z}xr| j8j7D]f}t|j1|j1}t23||| }| |j4|j4}||9}|j|j| }5|j6|j6|||dqbWqTWqWxT|j D]J}|j*j|j8j|j%j9}5|j*j6|j8j6||j%j||j%j|dqWx:|jD]0} x(| j:D]}5| j6|j6ddddqBWq6W|r|tjk rt |r| tj }d|kr|krnn;d<||jd kr|=r|>||?|fS|?|fS|=r|>|fSS)a Creates the OpenMM NonbondedForce instance Parameters ---------- nonbondedMethod : cutoff method This is the cutoff method. It can be either the NoCutoff, CutoffNonPeriodic, CutoffPeriodic, PME, or Ewald objects from the simtk.openmm.app namespace nonbondedCutoff : float or distance Quantity The nonbonded cutoff must be either a floating point number (interpreted as nanometers) or a Quantity with attached units. This is ignored if nonbondedMethod is NoCutoff. switchDistance : float or distance Quantity The distance at which the switching function is turned on for van der Waals interactions. This is ignored when no cutoff is used, and no switch is used if switchDistance is 0, negative, or greater than the cutoff ewaldErrorTolerance : float=0.0005 When using PME or Ewald, the Ewald parameters will be calculated from this value reactionFieldDielectric : float=78.5 If the nonbondedMethod is CutoffPeriodic or CutoffNonPeriodic, the region beyond the cutoff is treated using a reaction field method with this dielectric constant. It should be set to 1 if another implicit solvent model is being used (e.g., GB) Returns ------- NonbondedForce This just implements the very basic NonbondedForce with the typical charge-charge and 12-6 Lennard-Jones interactions with the Lorentz-Berthelot combining rules. Notes ----- Subclasses of Structure for which this nonbonded treatment is inadequate should override this method to implement what is needed. If nrexcl is set to 3 and no exception parameters are stored in the adjusts list, the 1-4 interactions are determined from the list of dihedrals Nz!Unrecognized nonbondedMethod (%s)r\g)N>?c s||kr dSx|jD]}||kr"q||k rB|j|jddddx&|jD]}|j|jddddqJWx&|jD]}|j|jddddqrWx6|jD],}x&|jD]}|j|jddddqWqW|||d|qWdS)Ngg?Tr)rQ addExceptionrchildren)originrNlevelendrchildZchild2) exclude_torrFrGr s      z1Structure.omm_nonbonded_force..exclude_torrdcSs d||S)Ng?rF)sig1sig2rFrFrG z/Structure.omm_nonbonded_force..rbcSst||S)N)rr)rrrFrFrGr rrzADetected scaling constants of 0 for dihedral containing 1-4 info!Tgg?)@rer:rrrrrrZNonbondedForcerNONBONDED_FORCE_GROUPr;r>rr7rsetNonbondedMethodCutoffNonPeriodicsetCutoffDistancerrZsetEwaldErrorTolerancerrZsetReactionFieldDielectricrr?sigmarrArrrurirjrr@rrrsceescnbrrTr}rUrr8rdZ epsilon_14rrZsigma_14rrrrZchgscaleZexclusion_partnerssetUseSwitchingFunctionsetSwitchingDistancer~_omm_nbfixed_force_omm_geometric_force)rrrrrZreactionFieldDielectricrene_convZ sigma_scalerNZcomb_sigZdihrrrrijwdijrij14wdij14ZepsprodZsigprodZchgprodrrrZpairrFrF)rrrGr s.                 &    &   "    *   $        zStructure.omm_nonbonded_forcec! Cstjtj}tjtj}dd|jD}gg}}d}g} xt|jD]\} } | j} || rbqJ|d7}||| <| j t | j f} | | | | j | t | j xnt | dt|jD]V}|j|j}||dkrq|| kr|||<q| js|j t |j f}| |kr|||<qWqJWddt ||D}|dd}xt |D]} xt |D]}| |j}y| | j|\}}}}WnDtk r|| |||}t|| |||}YnX||9}||9}|d}t|||| ||<d|||| ||<qVWqHWtd }|d t||||d t||||d ||j|tjks|tjks|tjkr|tjjnD|tj kr|tjj n(|tj!kr|tjj!n t"d |x|D]} |#| dfqWx6t |$D]&} |%| \}}}|&| |ddq Wx6t |'D]&}|(|\} }}}} |)| |qDW|*d|tj kr|tjj nN|tj!kr|tjj!n2|tjtjtjfkr|tjjn t"d||+|,|-r|.d|/|0|S)a Private method for creating a CustomNonbondedForce with a lookup table. This should not be called by users -- you have been warned. Parameters ---------- nonbfrc : NonbondedForce NonbondedForce for the "standard" nonbonded interactions. This will be modified (specifically, L-J ixns will be zeroed) nonbondedMethod : Nonbonded Method (e.g., NoCutoff, PME, etc.) The nonbonded method to apply here. Ewald and PME will be interpreted as CutoffPeriodic for the CustomNonbondedForce. Returns ------- force : CustomNonbondedForce The L-J force with NBFIX implemented as a lookup table cSsg|]}dqS)rrF)rrNrFrFrGr sz0Structure._omm_nbfixed_force..rrcSsg|]}dqS)rrF)rrrFrFrGr sNrHr\zP(a/r6)^2-b/r6; r6=r2*r2*r2; r2=r^2; a=acoef(type1, type2); b=bcoef(type1, type2)acoefbcoefrz"Unrecognized nonbonded method [%s]g?gTzUnsupported nonbonded method %s)1r:r?rrrrrerrTr@rrAr=r4rr}rTr8rrrCustomNonbondedForceZaddTabulatedFunctionZDiscrete2DFunctionaddPerParticleParameterrrr7rrrrrrrBrrgetParticleParameterssetParticleParametersgetNumExceptionsgetExceptionParameters addExclusionsetUseLongRangeCorrectionrgetCutoffDistancegetUseSwitchingFunctionrrgetSwitchingDistance)!rnonbfrcrrrZ lj_idx_listZlj_radiiZ lj_depthsZ num_lj_typesZ lj_type_listrrNZatypeZljtyperDZatype2Zljtype2rrZnamejrrrrZrij6rchgsigepsiiqqsseerFrFrGrm s         $              zStructure._omm_nbfixed_forcecCstjtj}tjtj}td}|d|d| |j |t j ksd|t j ksd|t jkrt|tjjn@|t jkr|tjjn&|t jkr|tjjn td|x@|jD]6}t|j|d}t|j|}|||fqWx6t|D]&} || \} }}|| | ddqWx6t|D]&} || \} } } }}|| | q@c Csddlm} m} m} m} m}yddlm}Wntk rHdd}YnX|dkrVdS| r`d}nd}|dkrrtj }|tj tj tj tj tjfkrtd|||t|}|dkr t|r|tjtj}|}t|r|tj}d t|||}|d 9}nt|r&|tjd }|tj kr8d}nt|rR|tj}n|}|tj krv| |||||d }n|tj kr| |||||d }nh|tj kr| |||||d }nH|tj kr| |||||d }n(|tjkr||||||d }ntd x0t|j|D] \}}||j gt!|q Wy |"Wnt#k rRYnX|tj krp|$t%j&j n>|tj'kr|$t%j&j'|(|n|$t%j&j)|(||*|j+|S)a Creates a Generalized Born force for running implicit solvent calculations Parameters ---------- implicitSolvent : app.HCT, app.OBC1, app.OBC2, app.GBn, app.GBn2 The Generalized Born implicit solvent model to use. nonbondedMethod : cutoff method This is the cutoff method. It can be either the NoCutoff, CutoffNonPeriodic, CutoffPeriodic, PME, or Ewald objects from the simtk.openmm.app namespace. Default is NoCutoff nonbondedCutoff : float or distance Quantity The nonbonded cutoff must be either a floating opint number (interpreted as nanometers) or a Quantity with attached units. This is ignored if nonbondedMethod is NoCutoff implicitSolventKappa : float or 1/distance Quantity = None This is the Debye kappa property related to modeling saltwater conditions in GB. It should have units of 1/distance (1/nanometers is assumed if no units present). A value of None means that kappa will be calculated from implicitSolventSaltConc (below) implicitSolventSaltConc : float or amount/volume Quantity=0 moles/liter If implicitSolventKappa is None, the kappa will be computed from the salt concentration. It should have units compatible with mol/L temperature : float or temperature Quantity = 298.15 kelvin This is only used to compute kappa from implicitSolventSaltConc soluteDielectric : float=1.0 The dielectric constant of the protein interior used in GB solventDielectric : float=78.5 The dielectric constant of the water used in GB r) GBSAHCTForce GBSAOBC1Force GBSAOBC2Force GBSAGBnForce GBSAGBn2Force)convertParameterscSs|S)NrF)ZparamsZchoicerFrFrGrM rz*Structure.omm_gbsa_force..NZACEz#Unrecognized implicit solvent modelgX2ı*I@g333333@r)Zkappaz7Unexpected implicit solvent model... should not be here),Z(simtk.openmm.app.internal.customgbforcesrrrrrr ImportErrorr7rZHCTrXrYGBnGBn2r_get_gb_parametersrUr:r;r>molesliterkelvinrrrrrBr5rerr?rfinalizerdrrZ CustomGBForcerrrrr)rrrrrrrrrrrrrrrrZsasaZgb_parmsZsccutoffrrNZparmsrFrFrGr s~)                         zStructure.omm_gbsa_forcec Cs|js dStjtj}t|jdr0t|jds8tdt }|jj }| |d| |d| |d||d||jxP|jD]F}|jdkrtd ||jj|jj|jj|jj|jj|jj|qW|S) z Creates the Amoeba trigonal-angle force Returns ------- AmoebaInPlaneAngleForce The trigonal in-plane Angle force Nr<rz5Do not have the trigonal angle force table parametersr]r^r_rHz!Missing trigonal angle parameters)rnr:rrrrrr rZAmoebaInPlaneAngleForcerZ setAmoebaGlobalInPlaneAngleCubicZ"setAmoebaGlobalInPlaneAngleQuarticZ!setAmoebaGlobalInPlaneAnglePenticZ!setAmoebaGlobalInPlaneAngleSexticrTRIGONAL_ANGLE_FORCE_GROUPrrrrrrrrr)rrrrangrFrFrGr s&      z"Structure.omm_trigonal_angle_forcec Cs|js dStjtj}t|jdr0t|jds8tdt }|jj }| |d| |d| |d||d||jxN|jD]D}|jdkrtd ||jj|jj|jj|jjd |jj|qW|S) z Creates the Amoeba out-of-plane bend force Returns ------- AmoebaOutOfPlaneBendForce The out-of-plane bend Angle force Nr<rz5Do not have the trigonal angle force table parametersr]r^r_rHz$Missing out-of-plane bend parametersr\)ror:rrrrrr rZAmoebaOutOfPlaneBendForcerZ"setAmoebaGlobalOutOfPlaneBendCubicZ$setAmoebaGlobalOutOfPlaneBendQuarticZ#setAmoebaGlobalOutOfPlaneBendPenticZ#setAmoebaGlobalOutOfPlaneBendSexticrOUT_OF_PLANE_BEND_FORCE_GROUPrZaddOutOfPlaneBendrrrrrr)rrrrrrFrFrGr s&      z%Structure.omm_out_of_plane_bend_forcec Cs|js dStjtj}t}||jxV|jD]L}|j dkrJt d| |j j |jj |jj |jj |jj |jj |j j|q4W|S)z Creates the Amoeba pi-torsion force Returns ------- AmoebaPiTorsionForce The pi-torsion force NzMissing pi-torsion parameters)rpr:rrrrZAmoebaPiTorsionForcerPI_TORSION_FORCE_GROUPrr Z addPiTorsionrrrrrrrr)rrrrrFrFrGr s    zStructure.omm_pi_torsion_forcec Cs|js dStjdd}tjtj}t}| |j xr|jD]h}|j dkrXt d| |jj|jj|jj|j j||j j||j jtjd|j j||j j|qBW|S)z Create the OpenMM Amoeba stretch-bend force for this system Returns ------- AmoebaStretchBendForce The stretch-bend force containing all terms in this system NgQD@zMissing stretch-bend parameters)rqrrr:r?rrrZAmoebaStretchBendForcerSTRETCH_BEND_FORCE_GROUPrr ZaddStretchBendrrrrZreq1Zreq2rZk1Zk2)rrrrZstrbndrFrFrGr s    z Structure.omm_stretch_bend_forcecCs|js dStddS)z Create the OpenMM Amoeba coupled-torsion (CMAP) force Returns ------- AmoebaTorsionTorsionForce The torsion-torsion (CMAP) force with all coupled-torsion parameters for this system Nz$Torsion-torsions not yet implemented)rrNotImplementedError)rrFrFrGr s z#Structure.omm_torsion_torsion_forcecCspxjttt|jD]T}|j|}|jdkr@|jdkr@|j|=q|jjdksX|jjdkr||j|=qWdS)z Gets rid of any empty bonds Nr)rr4rrgrrrr_)rrrirFrFrGr s  zStructure._prune_empty_bondscCs|xvttt|jD]`}|j|}|jdkrJ|jdkrJ|jdkrJ|j|=qd|jj|jj|jjfkr||j|=qWdS)z Gets rid of any empty angles Nr) rr4rrhrrrrr_)rrrrFrFrGr! s   zStructure._prune_empty_anglesricCsxtttt||D]|}t|||}|jdkr`|jdkr`|jdkr`|jdkr`t|||=qd|jj|jj|jj|jjfkr| t|||=qWdS)z! Gets rid of any empty dihedrals Nr) rr4rrrrrrrr_)rZdlistrZdihedrFrFrGr. s z Structure._prune_empty_dihedralscCs|ddS)z$ Gets rid of any empty R-B torsions rjN)r)rrFrFrGr; sz"Structure._prune_empty_rb_torsionscCspxjttt|jD]T}|j|}|jdkr@|jdkr@|j|=q|jjdksX|jjdkr||j|=qWdS)z* Gets rid of any empty Urey-Bradley terms Nr)rr4rrkrrrr_)rrrrFrFrGrA s  zStructure._prune_empty_ureyscCsxttt|jD]p}|j|}|jdkrT|jdkrT|jdkrT|jdkrT|j|=qd|jj|jj|jj|jjfkr| |j|=qWdS)z) Gets rid of any empty improper torsions Nr) rr4rrlrrrrrr_)rrrrFrFrGrM s (  z Structure._prune_empty_improperscCsxttt|jD]}|j|}|jdkr^|jdkr^|jdkr^|jdkr^|jdkr^|j|=qd|jj |jj |jj |jj |jj fkr| |j|=qWdS)z" Gets rid of any empty CMAP terms Nr) rr4rrmrrrrrrr_)rrrrFrFrGrY s  zStructure._prune_empty_cmapscCs|ddS)z' Gets rid of any empty trigonal angles rnN)r)rrFrFrGrg sz&Structure._prune_empty_trigonal_anglescCs|ddS)z* Gets rid of any empty out-of-plane bends roN)r)rrFrFrGrm sz)Structure._prune_empty_out_of_plane_bendscCsxttt|jD]}|j|}|jdkrh|jdkrh|jdkrh|jdkrh|jdkrh|j dkrh|j|=qd|jj |jj |jj |jj |jj |j j fkr|j|=qWdS)z# Gets rid of any empty pi-torsions Nr) rr4rrprrrrrrr)rrZpitrFrFrGrs s  z"Structure._prune_empty_pi_torsionscCs~xxttt|jD]b}|j|}|jdkrJ|jdkrJ|jdkrJ|j|=q|jjdksn|jjdksn|jjdkr|j|=qWdS)z* Gets rid of any empty stretch-bend terms Nr)rr4rrqrrrr)rrZsbrFrFrGr s   $z$Structure._prune_empty_stretch_bendscCsxttt|jD]}|j|}|jdkr^|jdkr^|jdkr^|jdkr^|jdkr^|j|=qd|jj |jj |jj |jj |jj fkr| |j|=qWdS)z- Gets rid of any empty torsion-torsion terms Nr) rr4rrrrrrrrrr_)rrrrFrFrGr s  &z'Structure._prune_empty_torsion_torsionscCshxbttt|jD]L}|j|}|jdks6|jdkr@|j|=q|jjdksX|jjdkr|j|=qWdS)z* Gets rid of any empty chiral frame terms Nr)rr4rrsrrr)rrZcfrFrFrGr s   z$Structure._prune_empty_chiral_framescCsHxBttt|jD],}|j|}|jdks8|jjdkr|j|=qWdS)z- Gets rid of any empty multipole frame terms Nr)rr4rrtrNr)rrZmfrFrFrGr s z'Structure._prune_empty_multipole_framescCshxbttt|jD]L}|j|}|jdks6|jdkr@|j|=q|jjdksX|jjdkr|j|=qWdS)z7 Gets rid of any empty nonbonded exception adjustments Nr)rr4rrurrr)rrZadjrFrFrGr s   zStructure._prune_empty_adjustsc s|tjkrdd|jD}dd|jD}xt|jD]\}}|jdkrRd||<nN|jdkrfd||<n:|jdkrzd ||<n&|jd krd ||<n|jd krd ||<||dkr6t|||<q6Wn|tjkr6dd|jD}dd|jD}dd|jD}dd|jD}dd|jD}xzt|jD]\}}|jdkrbd||<d||<d||<d||<n|jdkrd||<d||<d||<d||<n|jdkrd||<d||<d||<d||<nZ|jd krd ||<d!||<d"||<d#||<n,|jd krd$||<d!||<d"||<d#||<||s*t|||<q*Wndd%d|jD}d&d|jD}xBt|jD]4\}}||r~||sbt||\||<||<qbWt j t j fd'd|D}|tjkrt t|||||St t||S)(a} Gets the GB parameters for the requested GB model used by OpenMM Parameters ---------- implicitSolvent : app.HCT, app.OBC1, app.OBC2, app.GBn, or app.GBn2 The object specifying a particular GB model in OpenMM Returns ------- parameters : list of float List of parameters for the requested GB model cSsg|]}dqS)g?rF)rrNrFrFrGr sz0Structure._get_gb_parameters..cSsg|] }|jqSrF)solvent_radius)rrNrFrFrGr srHg,-?rgw#t?rIg1\^Yg?rPgdU ?rLgfE?rcSsg|]}dqS)g?rF)rrrFrFrGr scSsg|]}dqS)g?rF)rrrFrFrGr scSsg|]}dqS)gffffff@rF)rrrFrFrGr scSsg|]}dqS)g?rF)rrrFrFrGr scSsg|] }|jqSrF)r%)rrNrFrFrGr sg̰Q?gz?g>?4?g# Y?gZ ?gvۅ:?g"4?gr۾G?g,y?gzю?g@F?g4Op?g}?gY!?g>d ?g@-?gkтcSsg|] }|jqSrF)r%)rrNrFrFrGr scSsg|] }|jqSrF)screen)rrNrFrFrGr scsg|] }|qSrFrF)rr)rrFrGr s)r7rrerrerOrrVrZr:rrrrr5) rrr&ZradiirrNZalphaZbetaZgammarF)rrGr sx                       zStructure._get_gb_parameterscCst|dr|jr|jSt|S)NrT)rrTrepr)rrFrFrG__str__ szStructure.__str__cCst|}||7}|S)N)r)rothercprFrFrG__add__ szStructure.__add__c s>t|tstSd}tjyjdjd}Wntk rNd}YnXx:|jD]0}|j } t ||j |j ||j|j|jqXWfdd}||j|jjjddg||j|jjjddd g||j|jjjddd d d d g||j|jjjddd d d d g||j|jjjddg||j|jjjddd d g||j|jjjddd d d g||j|j jj ddd d g||j!|j"j!j"ddd d g||j#|j$j#j$ddd d d dg||j%|j&j%j&ddd g||j'|j(j'j(ddd d d g||j)gj)gdddg||j*gj*gdddddg||j+|j,j+j,ddg||j-gj-gddg||j.gj.gddg||j/gj/gdddg|dks|j0dkrd_0n`|d}t1|j2d|j2d}t3j4|d|ddddf|d|ddddffdd_0S)Nr!rrrc sdd|D}x|D]fdd|D}x2t|D]&\}}t|tr4j|j||<q4Wt} tdrjtkrt| d<n|rjdk r|jj| d<| t|| tdrj |d_ qW| |t|dr| dS) z Copies the valence terms from one list to another; oval=Other VALence; otyp=Other TYPe; sval=Self VALence; styp=Self TYPe; attrlist=ATTRibute LIST (atom1, atom2, ...) cSsg|] }t|qSrF)r)rrrFrFrGr6 szBStructure.__iadd__..copy_valence_terms..csg|]}t|qSrF)r)rr)rrFrGr8 srNrrr) rr@rrerr"rrrr=rextendr) r#r$r%r&r'r(r)rrr+)aoffsetr)rrGr-1 s"      z.Structure.__iadd__..copy_valence_termsrrrrrrrrrrNrrrrrr)Zaxis)5r@r[NotImplementedrxrrerfrr9rrrrTrrrrrgryrhrzrir{rjr~rkr|rlr}rmrrnrrorrprrqrrrrrsrtrurrwrvrxrminrvrZ concatenate) rr)ZmycrdroffsetrNrr-ZocrdZnframesrF)r-rrG__iadd__ st            @zStructure.__iadd__cCst|}|||S)z2 Replicates the current Structure `ncopies` times )r__imul__)rncopiesr*rFrFrG__mul__x szStructure.__mul__c st|tstSfdd}|dkr*t}xtt|dD]b}tj}jdjd}x:|jD]0}|j } t||j |j ||j |j|jqbW||j|jjddg||j|jjdddg||j|jjddddd d g||j|jjddddd d g||j|jjddg||j|jjddddg||j|jjddddd g||j|jjddddg||j |j j!ddddg||j"|j"j#ddddd d g||j$|j$j%dddg||j&|j&j'ddddd g||j(|j(gddd g||j)|j)gdddddg||j*|j*j+ddg||j,|j,gddg||j-|j-gddg||j.|j.gdddgq:Wj/dk rt01fddt|D_/S)Nc sx|D]fdd|D}x2t|D]&\}}t|tr&j|j|||<q&Wt}tdrtjtkrtt|d<n|rjdk r|jj|d<| t||tdrj |d_ qWdS)z Copies the valence terms from one list to another; oval=Other VALence; otyp=Other TYPe; sval=Self VALence; styp=Self TYPe; attrlist=ATTRibute LIST (atom1, atom2, ...) csg|]}t|qSrF)r)rr)rrFrGr szBStructure.__imul__..copy_valence_terms..rNrr) rr@rrerr"rrrr=r) r#r-r%r&r'r)rrr+)r)rrGr- s    z.Structure.__imul__..copy_valence_termsrrrrrrrrrrrrNrrrrrrcsg|] }jqSrF)r)rr)rrFrGr sz&Structure.__imul__..)2r@r1r.rr4rrerfrrrrTrrrrrgryrhrzrir{rjr~rkr|rlr}rmrrnrrorrprrqrrrrrsrtrurrwrvrxrrZhstack) rr3r)r-rr-r0rNrrF)rrGr2 sf               zStructure.__imul__cCst|jp|jp|jp|jp|jp|jp|jp|jp|j p|j p|j p|j p|j p|j p|jp|jp|jp|jp|jp|jp|jp|jp|jp|jp|jp|jp|jp|jp|jp|jp|jp|jp|j S)N)!r5rerfrgrhrirlrjrmrrrqrornrprkrsrtrurvrwrxryrzr{r|r}r~rrrrrr)rrFrFrG__bool__ s  zStructure.__bool__cCsJ|dkr dSt|ts t|tr.idxcs$g|]}|jj|jj|jfqSrF)rrrr)rr)rrFrGrsz*Structure.__getstate__..rgcs*g|]"}|jj|jj|jj|jfqSrF)rrrrr)rr)rrFrGrsrhc s8g|]0}|jj|jj|jj|jj|j|j|jfqSrF)rrrrrrrr)rr)rrFrGr srics0g|](}|jj|jj|jj|jj|jfqSrF)rrrrrr)rr)rrFrGr srlcs0g|](}|jj|jj|jj|jj|jfqSrF)rrrrrr)rr)rrFrGr srjcs$g|]}|jj|jj|jfqSrF)rrrr)rr:)rrFrGrsrkc s6g|].}|jj|jj|jj|jj|jj|jfqSrF)rrrrrrr)rr)rrFrGrsrmcs0g|](}|jj|jj|jj|jj|jfqSrF)rrrrrr)rr)rrFrGrsrncs0g|](}|jj|jj|jj|jj|jfqSrF)rrrrrr)rr)rrFrGrsroc s<g|]4}|jj|jj|jj|jj|jj|jj|jfqSrF)rrrrrrrr)rr)rrFrGrsrpcs*g|]"}|jj|jj|jj|jfqSrF)rrrrr)rr)rrFrGrsrqc s6g|].}|jj|jj|jj|jj|jj|jfqSrF)rrrrrrr)rr)rrFrGrsrrcSs g|]}|jj|jj|jfqSrF)rrrr)rrrFrFrGrsrscSs&g|]}|jj|j|j|j|jfqSrF)rNrrrrr)rr6rFrFrGrsrtcs$g|]}|jj|jj|jfqSrF)rrrr)rr\)rrFrGr!srucSsg|]}|jj|jjfqSrF)rrr)rrrFrFrGr"srvcSsg|]}|jj|jjfqSrF)rrr)rrrFrFrGr#srwcSs g|]}tdd|jDqS)css|] }|jVqdS)N)r)rr\rFrFrGr>$sz4Structure.__getstate__...)r7Z_exclusion_partners)rrrFrFrGr$s exclusions) experimentaljournalauthorskeywordsdoipmidjournal_authorsvolumeryearrrelated_entries)*r"rfryrzr{r|r}r~rrrrrrrrxrrrrrrrgrhrirlrjrkrmrnrorprqrrrsrtrurvrwrerrd)rZretdictr`rF)rrG __getstate__ sx              zStructure.__getstate__csx*dD]"}t|||t|qWx$dD]}||kr2t|||q2Wt_xjD]}j|jqbWddtfdd|dD_tfdd|d D_ tfd d|d D_ tfd d|d D_ tfdd|dD_ tfdd|dD_ tfdd|dD_tfdd|dD_tfdd|dD_tfdd|dD_tfdd|dD_tfdd|dD_tfdd|dD_tfd d|d!D_tfd"d|d#D_tfd$d|d%D_tfd&d|d'D_x>tj|d(D]*\}}x|D]}|j|qWqWdS))N)rfryrzr{r|r}r~rrrrrrxrr)r8r9r:r;r<r=r>r?rr@rrArrrrrrcSs|dkr dS||S)NrF)ZtypelistrrFrFrG assign_typeGsz+Structure.__setstate__..assign_typec3s<|]4}tj|dj|dj|ddVqdS)rrr\)rN)rrery)rr)rCrrFrGr>Nsz)Structure.__setstate__..rgc 3sH|]@}tj|dj|dj|dj|ddVqdS)rrr\r])rN)rrerz)rr)rCrrFrGr>Rsrhc 3s`|]X}tj|dj|dj|dj|d|d|dj|ddVqdS) rrr\r]r^r_rH)rrrN)rrer{)rr)rCrrFrGr>Wsric 3sT|]L}tj|dj|dj|dj|dj|ddVqdS)rrr\r]r^)rN)rrer})rr)rCrrFrGr>\srlc3s<|]4}tj|dj|dj|ddVqdS)rrr\)rN)r+rer|)rr)rCrrFrGr>`srkc 3sT|]L}tj|dj|dj|dj|dj|ddVqdS)rrr\r]r^)rN)rrer~)rr)rCrrFrGr>esrjc 3s`|]X}tj|dj|dj|dj|dj|dj|ddVqdS)rrr\r]r^r_)rN)rrer)rr)rCrrFrGr>jsrmc 3sT|]L}tj|dj|dj|dj|dj|ddVqdS)rrr\r]r^)rN)r(rer)rr)rCrrFrGr>osrnc 3sT|]L}tj|dj|dj|dj|dj|ddVqdS)rrr\r]r^)rN)r rer)rr)rCrrFrGr>tsroc 3sl|]d}tj|dj|dj|dj|dj|dj|dj|ddVqdS) rrr\r]r^r_rH)rN)r"rer)rr)rCrrFrGr>ysrpc 3sH|]@}tj|dj|dj|dj|ddVqdS)rrr\r])rN)r$rer)rr)rCrrFrGr>srqc 3s`|]X}tj|dj|dj|dj|dj|dj|ddVqdS)rrr\r]r^r_)rN)r&rer)rr)rCrrFrGr>srrc3s2|]*}tj|dj|d|dVqdS)rrr\N)rre)rr)rrFrGr>srsc3s.|]&}tj|df|ddVqdS)rrN)rre)rr)rrFrGr>srtc3s:|]2}tj|dj|dj|dVqdS)rrr\N)rrer)rr)rCrrFrGr>sruc3s,|]$}tj|dj|dVqdS)rrN)rre)rr)rrFrGr>srvc3s,|]$}tj|dj|dVqdS)rrN)rre)rr)rrFrGr>srwr7)setattrrrrrerfr,r'rgrhrirlrkrjrmrnrorprqrrrsrtrurvrwr5Zexclude)rrrr`rrNZexclrrF)rCrrG __setstate__1sx                 zStructure.__setstate__)rcrcrc)F)NF)r!)r!)NTT)NT)F)ri)N)hr __module__ __qualname____doc__rrrrrrrr r!r#ZTORSION_TORSION_FORCE_GROUPrrrrrrrrrrrrrrrrrr0r.propertyr<rGr]rarsetterr0rjrnrrxrryrXr~rr:r?rZlitersrrrrrrrrrrrrrrrrrrrrrrrrrrrrrrrrrrrrrr(r+r1r4__rmul__r2r5 __nonzero__ staticmethodrrBrErFrFrFrGr[sS5   i&.0  B & ;   3  "   V ) ? 1 -    &  ;mDr              R _ O  Fr[c@s eZdZdZddZddZdS)r;z= Class responsible for creating a StructureView when indexed cCs ||_dS)N)r,)rr,rFrFrGrsz _StructureViewerCreator.__init__c s|j}t|tr|j|St}||}|dkr6|St|trD|St|}|dkrX|S|dgx0tdt |D]} |d||qrWddt |Dt }xRt |jD]D\}}||sq|j ||j|krq|j |j||jqWfdd}||j|jddg||j|jddd g||j|jddd d g||j|jddd d g||j|jddg||j|jddd d g||j|jddd d d g||j|jddd d g||j|jddd d g||j|jddd d d d g||j|jddd g||j|jddd d d g||j|jddg||j|jd g||j|jddg||j |j ddg||j!|j!ddg|S)NrrcSsg|]\}}||qSrFrF)rrrrFrFrGrsz7_StructureViewerCreator.__getitem__..csJxD|D]<fdd|D}fdd|D}t|s8q|qWdS)z Adds the valence terms from one list to another; oval=Other VALence; otyp=Other TYPe; sval=Self VALence; styp=Self TYPe; attrlist=ATTRibute LIST (atom1, atom2, ...) csg|]}t|qSrF)r)rr)rrFrGrszR_StructureViewerCreator.__getitem__..add_valence_terms..cs g|]}t|tr|jqSrF)r@rr)rr)r rFrGrsN)r!r=)r#r%r'r)r*)r )rrGadd_valence_termss  z>_StructureViewerCreator.__getitem__..add_valence_termsrrrrrrrN)"r,r@r1re StructureViewr.rrr4rr=r5rrrrfrrgrhrirjrkrlrmrnrorprqrrrsrtrurwrv) rrr,r<r/rZsel_resrNrNrF)r rGr0sb             z#_StructureViewerCreator.__getitem__N)rrFrGrHrr0rFrFrFrGr;sr;c@s\eZdZdZddZddZddZedd Zed d Z d d Z ddZ e Z ddZ dS)rOaS A view of a Structure. In many cases, this can serve as a duck-typed Structure and it has many of the same attributes. However, none of its lists *own* their objects, and the lists of atoms, residues, and parameters/topologies are regular lists, rather than TrackedList instances. Therefore, the indexes correspond to the indexes from the *original* Structure from which this view was taken. Furthermore, there are no "type" lists, as they would be exactly equivalent to the type lists of the parent Structure instance. Attributes ---------- atoms : :class:`AtomList` List of all atoms in the structure residues : :class:`ResidueList` List of all residues in the structure bonds : :class:`TrackedList` (:class:`Bond`) List of all bonds in the structure angles : :class:`TrackedList` (:class:`Angle`) List of all angles in the structure dihedrals : :class:`TrackedList` (:class:`Dihedral`) List of all dihedrals in the structure rb_torsions : :class:`TrackedList` (:class:`RBTorsion`) List of all Ryckaert-Bellemans torsions in the structure urey_bradleys : :class:`TrackedList` (:class:`UreyBradley`) List of all Urey-Bradley angle bends in the structure impropers : :class:`TrackedList` (:class:`Improper`) List of all CHARMM-style improper torsions in the structure cmaps : :class:`TrackedList` (:class:`Cmap`) List of all CMAP objects in the structure trigonal_angles : :class:`TrackedList` (:class:`TrigonalAngle`) List of all AMOEBA-style trigonal angles in the structure out_of_plane_bends : :class:`TrackedList` (:class:`OutOfPlaneBends`) List of all AMOEBA-style out-of-plane bending angles pi_torsions : :class:`TrackedList` (:class:`PiTorsion`) List of all AMOEBA-style pi-torsion angles stretch_bends : :class:`TrackedList` (:class:`StretchBend`) List of all AMOEBA-style stretch-bend compound bond/angle terms torsion_torsions : :class:`TrackedList` (:class:`TorsionTorsion`) List of all AMOEBA-style coupled torsion-torsion terms chiral_frames : :class:`TrackedList` (:class:`ChiralFrame`) List of all AMOEBA-style chiral frames defined in the structure multipole_frames : :class:`TrackedList` (:class:`MultipoleFrame`) List of all AMOEBA-style multipole frames defined in the structure adjusts : :class:`TrackedList` (:class:`NonbondedException`) List of all nonbonded pair-exception rules acceptors : :class:`TrackedList` (:class:`AcceptorDonor`) List of all H-bond acceptors, if that information is present donors : :class:`TrackedList` (:class:`AcceptorDonor`) List of all H-bond donors, if that information is present positions : u.Quantity(list(Vec3), u.angstroms) Unit-bearing atomic coordinates. If not all atoms have coordinates, this property is None coordinates : np.ndarray of shape (nframes, natom, 3) If no coordinates are set, this is set to None. The first frame will match the coordinates present on the atoms. cCsvg|_g|_g|_g|_g|_g|_g|_g|_g|_g|_ g|_ g|_ g|_ g|_ g|_g|_g|_g|_g|_dS)N)rerfrgrhrirjrkrlrmrnrorprqrrrsrtrurvrw)rrFrFrGr*s&zStructureView.__init__cCsddlm}||S)aM Generates a DataFrame from the current Structure's atomic properties Returns ------- df : DataFrame DataFrame with all atomic properties See Also -------- :func:`parmed.utils.pandautils.create_dataframe` for full documentation of the generated DataFrame r)r)rr)rrrFrFrGr?s zStructureView.to_dataframecCsddlm}|||S)av Loads atomic properties from an input DataFrame Parameters ---------- df : pandas.DataFrame A pandas DataFrame with atomic properties that will be used to set the properties on the current list of atoms See Also -------- :func:`parmed.utils.pandautils.load_dataframe` for full documentation r)r)rr)rrrrFrFrGrOs zStructureView.load_dataframecCs0ytdd|jDStk r*dSXdS)NcSsg|]}|j|j|jgqSrF)rkrlrm)rrrFrFrGrbsz-StructureView.coordinates..)rrrerd)rrFrFrGr_szStructureView.coordinatescCs0ydd|jDtjStk r*dSXdS)NcSsg|]}t|j|j|jqSrF)r6rkrlrm)rrrFrFrGrisz+StructureView.positions..)rer:r?rd)rrFrFrGrnfszStructureView.positionscCs|t|jpx|jpx|jpx|jpx|jpx|jpx|jpx|jpx|j px|j px|j px|j px|j px|j px|jpx|jpx|jpx|jpx|jpx|jS)N)r5rerfrgrhrirlrjrmrrrqrornrprkrsrtrurvrw)rrFrFrGr5ms   zStructureView.__bool__cCst|j}t|j}tdd|jD}dt|j|fg}|dkrR|d||d|t|j}|d|d|S) NcSsg|]}t|tqSrF)r@r)rrrFrFrGr{sz*StructureView.__repr__..z <%s %d atomsrz [%d EPs]z ; %d residuesz ; %d bonds>rc) rrerfrrrr=rgr)rrrrrrrFrFrGrxs   zStructureView.__repr__cCs t|jS)N)rBre)rrFrFrGrszStructureView.__iter__N)rrFrGrHrrrrIrrnr5rrLrrFrFrFrGrOs9    rOcCsg}xv|jD]l}tdd|D}|dkr6|dq x@|D].}x&|jD]}|j|k rH|dPqHWqsz_settler..r]FTzIncorrect length)rfrr=rQrrrB)rrrZnarrFrFrFrGrs       rcCs|j|kr||jSt|jdkrDtj|jrD|tj|jjStj|jrf|tj|jjStj|jrtj|jjdkr|tj|jjSdS)z: Returns the residue template inside lib that matches res r]rN) rTrrrrrrZ DNAResidueZ RNAResidue)rrrFrFrGr s  "r )^rHZ __future__rrZloggingrrO collectionsrrZnumpyrrcrr:rZ constantsr r exceptionsr Zgeometryr r rrZtopologyobjectsrrrrrrrrrrrrrrrrr r!r"r#r$r%r&r'r(r)r*r+r,Zutilsr-r.r/Zutils.decoratorsr0Z utils.sixr1r2r3Zutils.six.movesr4r5Zvec3r6Z simtk.openmmr7Zsimtkr8rZ"simtk.openmm.app.internal.unitcellr9rZ getLoggerrr rCrOrRrSrVrZobjectr[r;rOrr rFrFrFrGs     |        $.N