| 307 | 
  | 
 | 
| 308 | 
  | 
 | 
| 309 | 
  | 
  if( simnfo->n_SRI ){ | 
| 310 | 
< | 
    the_sris = new SRI*[simnfo->n_SRI]; | 
| 311 | 
< | 
    the_excludes = new int[2 * simnfo->n_SRI]; | 
| 310 | 
> | 
    Exclude::createArray(simnfo->n_SRI); | 
| 311 | 
> | 
    the_excludes = new Exclude*[simnfo->n_SRI]; | 
| 312 | 
  | 
    simnfo->globalExcludes = new int; | 
| 313 | 
  | 
    simnfo->n_exclude = tot_SRI; | 
| 314 | 
  | 
  } | 
| 315 | 
  | 
  else{ | 
| 316 | 
  | 
     | 
| 317 | 
< | 
    the_excludes = new int[2]; | 
| 318 | 
< | 
    the_excludes[0] = 0; | 
| 319 | 
< | 
    the_excludes[1] = 0; | 
| 317 | 
> | 
    Exclude::createArray( 1 ); | 
| 318 | 
> | 
    the_excludes = new Exclude*; | 
| 319 | 
> | 
    the_excludes[0] = new Exclude(0); | 
| 320 | 
> | 
    the_excludes[0]->setPair( 0,0 ); | 
| 321 | 
  | 
    simnfo->globalExcludes = new int; | 
| 322 | 
  | 
    simnfo->globalExcludes[0] = 0; | 
| 323 | 
< | 
 | 
| 323 | 
< | 
    simnfo->n_exclude = 1; | 
| 323 | 
> | 
    simnfo->n_exclude = 0; | 
| 324 | 
  | 
  } | 
| 325 | 
  | 
 | 
| 326 | 
  | 
  // set the arrays into the SimInfo object | 
| 674 | 
  | 
 | 
| 675 | 
  | 
void SimSetup::makeMolecules( void ){ | 
| 676 | 
  | 
 | 
| 677 | 
< | 
  int i, j, exI, exJ, tempEx, stampID, atomOffset; | 
| 677 | 
> | 
  int i, j, exI, exJ, tempEx, stampID, atomOffset, excludeOffset; | 
| 678 | 
  | 
  molInit info; | 
| 679 | 
  | 
  DirectionalAtom* dAtom; | 
| 680 | 
+ | 
  LinkedAssign* extras; | 
| 681 | 
+ | 
  LinkedAssign* current_extra; | 
| 682 | 
  | 
  AtomStamp* currentAtom; | 
| 683 | 
  | 
  BondStamp* currentBond; | 
| 684 | 
  | 
  BendStamp* currentBend; | 
| 692 | 
  | 
  // init the molecules | 
| 693 | 
  | 
 | 
| 694 | 
  | 
  atomOffset = 0; | 
| 695 | 
+ | 
  excludeOffset = 0; | 
| 696 | 
  | 
  for(i=0; i<simnfo->n_mol; i++){ | 
| 697 | 
  | 
     | 
| 698 | 
  | 
    stampID = the_molecules[i].getStampID(); | 
| 701 | 
  | 
    info.nBonds    = comp_stamps[stampID]->getNBonds(); | 
| 702 | 
  | 
    info.nBends    = comp_stamps[stampID]->getNBends(); | 
| 703 | 
  | 
    info.nTorsions = comp_stamps[stampID]->getNTorsions(); | 
| 704 | 
< | 
     | 
| 704 | 
> | 
    info.nExcludes = info.nBonds + info.nBends + info.nTorsions; | 
| 705 | 
> | 
 | 
| 706 | 
  | 
    info.myAtoms = &the_atoms[atomOffset]; | 
| 707 | 
+ | 
    info.myExcludes = &the_excludes[excludeOffset]; | 
| 708 | 
  | 
    info.myBonds = new Bond*[info.nBonds]; | 
| 709 | 
  | 
    info.myBends = new Bend*[info.nBends]; | 
| 710 | 
  | 
    info.myTorsions = new Torsions*[info.nTorsions]; | 
| 752 | 
  | 
    }  | 
| 753 | 
  | 
     | 
| 754 | 
  | 
    // make the bonds | 
| 755 | 
< | 
    for(j=0; j<nBonds; j++){ | 
| 755 | 
> | 
    for(j=0; j<info.nBonds; j++){ | 
| 756 | 
  | 
       | 
| 757 | 
  | 
      currentBond = comp_stamps[stampID]->getBond( j ); | 
| 758 | 
  | 
      theBonds[j].a = currentBond->getA() + atomOffset; | 
| 768 | 
  | 
        exJ = tempEx; | 
| 769 | 
  | 
      } | 
| 770 | 
  | 
#ifdef IS_MPI | 
| 771 | 
< | 
       | 
| 772 | 
< | 
      the_excludes[index*2] =      | 
| 773 | 
< | 
        the_atoms[exI]->getGlobalIndex() + 1; | 
| 774 | 
< | 
      the_excludes[index*2 + 1] =  | 
| 770 | 
< | 
        the_atoms[exJ]->getGlobalIndex() + 1; | 
| 771 | 
> | 
      tempEx = exI; | 
| 772 | 
> | 
      exI = the_atoms[tempEx]->getGlobalIndex() + 1; | 
| 773 | 
> | 
      tempEx = exJ; | 
| 774 | 
> | 
      exJ = the_atoms[tempEx]->getGlobalIndex() + 1; | 
| 775 | 
  | 
       | 
| 776 | 
+ | 
      the_excludes[j+excludeOffset]->setPair( exI, exJ ); | 
| 777 | 
  | 
#else  // isn't MPI | 
| 778 | 
< | 
       | 
| 774 | 
< | 
      the_excludes[index*2] =     exI + 1;  | 
| 775 | 
< | 
      the_excludes[index*2 + 1] = exJ + 1; | 
| 776 | 
< | 
      // fortran index from 1 (hence the +1 in the indexing) | 
| 777 | 
< | 
 | 
| 778 | 
> | 
      the_excludes[j+excludeOffset]->setPair( (exI+1), (exJ+1) ); | 
| 779 | 
  | 
#endif  //is_mpi | 
| 780 | 
< | 
     | 
| 781 | 
< | 
  } | 
| 781 | 
< | 
 | 
| 782 | 
< | 
 | 
| 783 | 
< | 
 | 
| 784 | 
< | 
 | 
| 785 | 
< | 
 | 
| 786 | 
< | 
 | 
| 787 | 
< | 
 | 
| 788 | 
< | 
 | 
| 789 | 
< | 
 | 
| 790 | 
< | 
 | 
| 791 | 
< | 
 | 
| 792 | 
< | 
 | 
| 793 | 
< | 
 | 
| 794 | 
< | 
 | 
| 795 | 
< | 
 | 
| 796 | 
< | 
void SimSetup::makeAtoms( void ){ | 
| 797 | 
< | 
 | 
| 798 | 
< | 
  int i, j, k, index; | 
| 799 | 
< | 
  double ux, uy, uz, uSqr, u; | 
| 800 | 
< | 
  AtomStamp* current_atom; | 
| 801 | 
< | 
 | 
| 802 | 
< | 
  DirectionalAtom* dAtom; | 
| 803 | 
< | 
  int molIndex, molStart, molEnd, nMemb, lMolIndex; | 
| 780 | 
> | 
    } | 
| 781 | 
> | 
    excludeOffset += info.nBonds; | 
| 782 | 
  | 
 | 
| 783 | 
< | 
  lMolIndex = 0; | 
| 784 | 
< | 
  molIndex = 0; | 
| 785 | 
< | 
  index = 0; | 
| 786 | 
< | 
  for( i=0; i<n_components; i++ ){ | 
| 787 | 
< | 
 | 
| 788 | 
< | 
    for( j=0; j<components_nmol[i]; j++ ){ | 
| 789 | 
< | 
 | 
| 812 | 
< | 
#ifdef IS_MPI | 
| 813 | 
< | 
      if( mpiSim->getMyMolStart() <= molIndex && | 
| 814 | 
< | 
          molIndex <= mpiSim->getMyMolEnd() ){ | 
| 815 | 
< | 
#endif // is_mpi         | 
| 816 | 
< | 
 | 
| 817 | 
< | 
        molStart = index; | 
| 818 | 
< | 
        nMemb = comp_stamps[i]->getNAtoms(); | 
| 819 | 
< | 
        for( k=0; k<comp_stamps[i]->getNAtoms(); k++ ){ | 
| 783 | 
> | 
    //make the bends | 
| 784 | 
> | 
    for(j=0; j<info.nBends; j++){ | 
| 785 | 
> | 
       | 
| 786 | 
> | 
      currentBend = comp_stamps[stampID]->getBend( j ); | 
| 787 | 
> | 
      theBends[j].a = currentBend->getA() + atomOffset; | 
| 788 | 
> | 
      theBends[j].b = currentBend->getB() + atomOffset; | 
| 789 | 
> | 
      theBends[j].c = currentBend->getC() + atomOffset; | 
| 790 | 
  | 
           | 
| 791 | 
< | 
          current_atom = comp_stamps[i]->getAtom( k ); | 
| 822 | 
< | 
          if( current_atom->haveOrientation() ){ | 
| 791 | 
> | 
      if( currentBend->haveExtras() ){ | 
| 792 | 
  | 
             | 
| 793 | 
< | 
            dAtom = new DirectionalAtom(index); | 
| 794 | 
< | 
            simnfo->n_oriented++; | 
| 826 | 
< | 
            the_atoms[index] = dAtom; | 
| 793 | 
> | 
        extras = current_bend->getExtras(); | 
| 794 | 
> | 
        current_extra = extras; | 
| 795 | 
  | 
             | 
| 796 | 
< | 
            ux = current_atom->getOrntX(); | 
| 797 | 
< | 
            uy = current_atom->getOrntY(); | 
| 830 | 
< | 
            uz = current_atom->getOrntZ(); | 
| 831 | 
< | 
             | 
| 832 | 
< | 
            uSqr = (ux * ux) + (uy * uy) + (uz * uz); | 
| 833 | 
< | 
             | 
| 834 | 
< | 
            u = sqrt( uSqr ); | 
| 835 | 
< | 
            ux = ux / u; | 
| 836 | 
< | 
            uy = uy / u; | 
| 837 | 
< | 
            uz = uz / u; | 
| 838 | 
< | 
             | 
| 839 | 
< | 
            dAtom->setSUx( ux ); | 
| 840 | 
< | 
            dAtom->setSUy( uy ); | 
| 841 | 
< | 
            dAtom->setSUz( uz ); | 
| 842 | 
< | 
          } | 
| 843 | 
< | 
          else{ | 
| 844 | 
< | 
            the_atoms[index] = new GeneralAtom(index); | 
| 845 | 
< | 
          } | 
| 846 | 
< | 
          the_atoms[index]->setType( current_atom->getType() ); | 
| 847 | 
< | 
          the_atoms[index]->setIndex( index ); | 
| 848 | 
< | 
           | 
| 849 | 
< | 
          // increment the index and repeat; | 
| 850 | 
< | 
          index++; | 
| 851 | 
< | 
        } | 
| 852 | 
< | 
         | 
| 853 | 
< | 
        molEnd = index -1; | 
| 854 | 
< | 
        the_molecules[lMolIndex].setNMembers( nMemb ); | 
| 855 | 
< | 
        the_molecules[lMolIndex].setStartAtom( molStart ); | 
| 856 | 
< | 
        the_molecules[lMolIndex].setEndAtom( molEnd ); | 
| 857 | 
< | 
        the_molecules[lMolIndex].setStampID( i ); | 
| 858 | 
< | 
        lMolIndex++; | 
| 859 | 
< | 
 | 
| 860 | 
< | 
#ifdef IS_MPI | 
| 861 | 
< | 
      } | 
| 862 | 
< | 
#endif //is_mpi | 
| 863 | 
< | 
       | 
| 864 | 
< | 
      molIndex++; | 
| 865 | 
< | 
    } | 
| 866 | 
< | 
  } | 
| 867 | 
< | 
 | 
| 868 | 
< | 
#ifdef IS_MPI | 
| 869 | 
< | 
    for( i=0; i<mpiSim->getMyNlocal(); i++ ) the_atoms[i]->setGlobalIndex( globalIndex[i] ); | 
| 870 | 
< | 
     | 
| 871 | 
< | 
    delete[] globalIndex; | 
| 872 | 
< | 
 | 
| 873 | 
< | 
    mpiSim->mpiRefresh(); | 
| 874 | 
< | 
#endif //IS_MPI | 
| 875 | 
< | 
           | 
| 876 | 
< | 
  the_ff->initializeAtoms(); | 
| 877 | 
< | 
} | 
| 878 | 
< | 
 | 
| 879 | 
< | 
void SimSetup::makeBonds( void ){ | 
| 880 | 
< | 
 | 
| 881 | 
< | 
  int i, j, k, index, offset, molIndex, exI, exJ, tempEx; | 
| 882 | 
< | 
  bond_pair* the_bonds; | 
| 883 | 
< | 
  BondStamp* current_bond; | 
| 884 | 
< | 
 | 
| 885 | 
< | 
  the_bonds = new bond_pair[tot_bonds]; | 
| 886 | 
< | 
  index = 0; | 
| 887 | 
< | 
  offset = 0; | 
| 888 | 
< | 
  molIndex = 0; | 
| 889 | 
< | 
 | 
| 890 | 
< | 
  for( i=0; i<n_components; i++ ){ | 
| 891 | 
< | 
 | 
| 892 | 
< | 
    for( j=0; j<components_nmol[i]; j++ ){ | 
| 893 | 
< | 
 | 
| 894 | 
< | 
#ifdef IS_MPI | 
| 895 | 
< | 
      if( mpiSim->getMyMolStart() <= molIndex && | 
| 896 | 
< | 
          molIndex <= mpiSim->getMyMolEnd() ){ | 
| 897 | 
< | 
#endif // is_mpi         | 
| 898 | 
< | 
         | 
| 899 | 
< | 
        for( k=0; k<comp_stamps[i]->getNBonds(); k++ ){ | 
| 900 | 
< | 
           | 
| 901 | 
< | 
          current_bond = comp_stamps[i]->getBond( k ); | 
| 902 | 
< | 
          the_bonds[index].a = current_bond->getA() + offset; | 
| 903 | 
< | 
          the_bonds[index].b = current_bond->getB() + offset; | 
| 904 | 
< | 
 | 
| 905 | 
< | 
          exI = the_bonds[index].a; | 
| 906 | 
< | 
          exJ = the_bonds[index].b; | 
| 907 | 
< | 
 | 
| 908 | 
< | 
          // exclude_I must always be the smaller of the pair | 
| 909 | 
< | 
          if( exI > exJ ){ | 
| 910 | 
< | 
            tempEx = exI; | 
| 911 | 
< | 
            exI = exJ; | 
| 912 | 
< | 
            exJ = tempEx; | 
| 913 | 
< | 
          } | 
| 914 | 
< | 
 | 
| 915 | 
< | 
           | 
| 916 | 
< | 
#ifdef IS_MPI | 
| 917 | 
< | 
 | 
| 918 | 
< | 
          the_excludes[index*2] =      | 
| 919 | 
< | 
            the_atoms[exI]->getGlobalIndex() + 1; | 
| 920 | 
< | 
          the_excludes[index*2 + 1] =  | 
| 921 | 
< | 
            the_atoms[exJ]->getGlobalIndex() + 1; | 
| 922 | 
< | 
 | 
| 923 | 
< | 
#else  // isn't MPI | 
| 924 | 
< | 
           | 
| 925 | 
< | 
          the_excludes[index*2] =     exI + 1;  | 
| 926 | 
< | 
          the_excludes[index*2 + 1] = exJ + 1; | 
| 927 | 
< | 
          // fortran index from 1 (hence the +1 in the indexing) | 
| 928 | 
< | 
#endif  //is_mpi | 
| 929 | 
< | 
           | 
| 930 | 
< | 
          // increment the index and repeat; | 
| 931 | 
< | 
          index++; | 
| 932 | 
< | 
        } | 
| 933 | 
< | 
        offset += comp_stamps[i]->getNAtoms(); | 
| 934 | 
< | 
         | 
| 935 | 
< | 
#ifdef IS_MPI | 
| 936 | 
< | 
      } | 
| 937 | 
< | 
#endif //is_mpi | 
| 938 | 
< | 
       | 
| 939 | 
< | 
      molIndex++; | 
| 940 | 
< | 
    }       | 
| 941 | 
< | 
  } | 
| 942 | 
< | 
 | 
| 943 | 
< | 
  the_ff->initializeBonds( the_bonds ); | 
| 944 | 
< | 
} | 
| 945 | 
< | 
 | 
| 946 | 
< | 
void SimSetup::makeBends( void ){ | 
| 947 | 
< | 
 | 
| 948 | 
< | 
  int i, j, k, index, offset, molIndex, exI, exJ, tempEx; | 
| 949 | 
< | 
  bend_set* the_bends; | 
| 950 | 
< | 
  BendStamp* current_bend; | 
| 951 | 
< | 
  LinkedAssign* extras; | 
| 952 | 
< | 
  LinkedAssign* current_extra; | 
| 953 | 
< | 
   | 
| 954 | 
< | 
 | 
| 955 | 
< | 
  the_bends = new bend_set[tot_bends]; | 
| 956 | 
< | 
  index = 0; | 
| 957 | 
< | 
  offset = 0; | 
| 958 | 
< | 
  molIndex = 0; | 
| 959 | 
< | 
  for( i=0; i<n_components; i++ ){ | 
| 960 | 
< | 
 | 
| 961 | 
< | 
    for( j=0; j<components_nmol[i]; j++ ){ | 
| 962 | 
< | 
 | 
| 963 | 
< | 
#ifdef IS_MPI | 
| 964 | 
< | 
      if( mpiSim->getMyMolStart() <= molIndex && | 
| 965 | 
< | 
          molIndex <= mpiSim->getMyMolEnd() ){ | 
| 966 | 
< | 
#endif // is_mpi         | 
| 967 | 
< | 
 | 
| 968 | 
< | 
        for( k=0; k<comp_stamps[i]->getNBends(); k++ ){ | 
| 969 | 
< | 
           | 
| 970 | 
< | 
          current_bend = comp_stamps[i]->getBend( k ); | 
| 971 | 
< | 
          the_bends[index].a = current_bend->getA() + offset; | 
| 972 | 
< | 
          the_bends[index].b = current_bend->getB() + offset; | 
| 973 | 
< | 
          the_bends[index].c = current_bend->getC() + offset; | 
| 974 | 
< | 
           | 
| 975 | 
< | 
          if( current_bend->haveExtras() ){ | 
| 976 | 
< | 
             | 
| 977 | 
< | 
            extras = current_bend->getExtras(); | 
| 978 | 
< | 
            current_extra = extras; | 
| 979 | 
< | 
             | 
| 980 | 
< | 
            while( current_extra != NULL ){ | 
| 981 | 
< | 
              if( !strcmp( current_extra->getlhs(), "ghostVectorSource" )){ | 
| 796 | 
> | 
        while( current_extra != NULL ){ | 
| 797 | 
> | 
          if( !strcmp( current_extra->getlhs(), "ghostVectorSource" )){ | 
| 798 | 
  | 
                 | 
| 799 | 
< | 
                switch( current_extra->getType() ){ | 
| 799 | 
> | 
            switch( current_extra->getType() ){ | 
| 800 | 
> | 
               | 
| 801 | 
> | 
            case 0: | 
| 802 | 
> | 
              theBends[j].ghost = | 
| 803 | 
> | 
                current_extra->getInt() + atomOffset; | 
| 804 | 
> | 
              theBends[j].isGhost = 1; | 
| 805 | 
> | 
              break; | 
| 806 | 
  | 
                   | 
| 807 | 
< | 
                case 0: | 
| 808 | 
< | 
                  the_bends[index].ghost = | 
| 809 | 
< | 
                    current_extra->getInt() + offset; | 
| 810 | 
< | 
                  the_bends[index].isGhost = 1; | 
| 811 | 
< | 
                  break; | 
| 990 | 
< | 
                   | 
| 991 | 
< | 
                case 1: | 
| 992 | 
< | 
                  the_bends[index].ghost =  | 
| 993 | 
< | 
                    (int)current_extra->getDouble() + offset; | 
| 994 | 
< | 
                  the_bends[index].isGhost = 1; | 
| 995 | 
< | 
                  break; | 
| 996 | 
< | 
                   | 
| 997 | 
< | 
                default: | 
| 998 | 
< | 
                  sprintf( painCave.errMsg, | 
| 999 | 
< | 
                           "SimSetup Error: ghostVectorSource was neiter a " | 
| 1000 | 
< | 
                           "double nor an int.\n" | 
| 1001 | 
< | 
                           "-->Bend[%d] in %s\n", | 
| 1002 | 
< | 
                           k, comp_stamps[i]->getID() ); | 
| 1003 | 
< | 
                  painCave.isFatal = 1; | 
| 1004 | 
< | 
                  simError(); | 
| 1005 | 
< | 
                } | 
| 1006 | 
< | 
              } | 
| 807 | 
> | 
            case 1: | 
| 808 | 
> | 
              theBends[j].ghost =  | 
| 809 | 
> | 
                (int)current_extra->getDouble() + atomOffset; | 
| 810 | 
> | 
              theBends[j].isGhost = 1; | 
| 811 | 
> | 
              break; | 
| 812 | 
  | 
               | 
| 813 | 
< | 
              else{ | 
| 814 | 
< | 
                 | 
| 815 | 
< | 
                sprintf( painCave.errMsg, | 
| 816 | 
< | 
                         "SimSetup Error: unhandled bend assignment:\n" | 
| 817 | 
< | 
                         "    -->%s in Bend[%d] in %s\n", | 
| 818 | 
< | 
                         current_extra->getlhs(), | 
| 819 | 
< | 
                         k, comp_stamps[i]->getID() ); | 
| 820 | 
< | 
                painCave.isFatal = 1; | 
| 1016 | 
< | 
                simError(); | 
| 1017 | 
< | 
              } | 
| 1018 | 
< | 
               | 
| 1019 | 
< | 
              current_extra = current_extra->getNext(); | 
| 813 | 
> | 
            default: | 
| 814 | 
> | 
              sprintf( painCave.errMsg, | 
| 815 | 
> | 
                       "SimSetup Error: ghostVectorSource was neiter a " | 
| 816 | 
> | 
                       "double nor an int.\n" | 
| 817 | 
> | 
                       "-->Bend[%d] in %s\n", | 
| 818 | 
> | 
                       j, comp_stamps[stampID]->getID() ); | 
| 819 | 
> | 
              painCave.isFatal = 1; | 
| 820 | 
> | 
              simError(); | 
| 821 | 
  | 
            } | 
| 822 | 
  | 
          } | 
| 823 | 
  | 
           | 
| 1023 | 
– | 
          if( !the_bends[index].isGhost ){ | 
| 1024 | 
– | 
             | 
| 1025 | 
– | 
            exI = the_bends[index].a; | 
| 1026 | 
– | 
            exJ = the_bends[index].c; | 
| 1027 | 
– | 
          } | 
| 824 | 
  | 
          else{ | 
| 825 | 
  | 
             | 
| 826 | 
< | 
            exI = the_bends[index].a; | 
| 827 | 
< | 
            exJ = the_bends[index].b; | 
| 826 | 
> | 
            sprintf( painCave.errMsg, | 
| 827 | 
> | 
                     "SimSetup Error: unhandled bend assignment:\n" | 
| 828 | 
> | 
                     "    -->%s in Bend[%d] in %s\n", | 
| 829 | 
> | 
                     current_extra->getlhs(), | 
| 830 | 
> | 
                     j, comp_stamps[stampID]->getID() ); | 
| 831 | 
> | 
            painCave.isFatal = 1; | 
| 832 | 
> | 
            simError(); | 
| 833 | 
  | 
          } | 
| 834 | 
  | 
           | 
| 835 | 
< | 
          // exclude_I must always be the smaller of the pair | 
| 1035 | 
< | 
          if( exI > exJ ){ | 
| 1036 | 
< | 
            tempEx = exI; | 
| 1037 | 
< | 
            exI = exJ; | 
| 1038 | 
< | 
            exJ = tempEx; | 
| 1039 | 
< | 
          } | 
| 1040 | 
< | 
 | 
| 1041 | 
< | 
 | 
| 1042 | 
< | 
#ifdef IS_MPI | 
| 1043 | 
< | 
 | 
| 1044 | 
< | 
          the_excludes[(index + tot_bonds)*2] =      | 
| 1045 | 
< | 
            the_atoms[exI]->getGlobalIndex() + 1; | 
| 1046 | 
< | 
          the_excludes[(index + tot_bonds)*2 + 1] =  | 
| 1047 | 
< | 
            the_atoms[exJ]->getGlobalIndex() + 1; | 
| 1048 | 
< | 
           | 
| 1049 | 
< | 
#else  // isn't MPI | 
| 1050 | 
< | 
           | 
| 1051 | 
< | 
          the_excludes[(index + tot_bonds)*2] =     exI + 1;  | 
| 1052 | 
< | 
          the_excludes[(index + tot_bonds)*2 + 1] = exJ + 1; | 
| 1053 | 
< | 
          // fortran index from 1 (hence the +1 in the indexing) | 
| 1054 | 
< | 
#endif  //is_mpi | 
| 1055 | 
< | 
           | 
| 1056 | 
< | 
           | 
| 1057 | 
< | 
          // increment the index and repeat; | 
| 1058 | 
< | 
          index++; | 
| 835 | 
> | 
          current_extra = current_extra->getNext(); | 
| 836 | 
  | 
        } | 
| 837 | 
< | 
        offset += comp_stamps[i]->getNAtoms(); | 
| 837 | 
> | 
      } | 
| 838 | 
> | 
           | 
| 839 | 
> | 
      if( !theBends[j].isGhost ){ | 
| 840 | 
> | 
             | 
| 841 | 
> | 
        exI = theBends[j].a; | 
| 842 | 
> | 
        exJ = theBends[j].c; | 
| 843 | 
> | 
      } | 
| 844 | 
> | 
      else{ | 
| 845 | 
  | 
         | 
| 846 | 
< | 
#ifdef IS_MPI | 
| 846 | 
> | 
        exI = theBends[j].a; | 
| 847 | 
> | 
        exJ = theBends[j].b; | 
| 848 | 
  | 
      } | 
| 849 | 
< | 
#endif //is_mpi | 
| 850 | 
< | 
 | 
| 851 | 
< | 
      molIndex++; | 
| 852 | 
< | 
    } | 
| 853 | 
< | 
  } | 
| 849 | 
> | 
       | 
| 850 | 
> | 
      // exclude_I must always be the smaller of the pair | 
| 851 | 
> | 
      if( exI > exJ ){ | 
| 852 | 
> | 
        tempEx = exI; | 
| 853 | 
> | 
        exI = exJ; | 
| 854 | 
> | 
        exJ = tempEx; | 
| 855 | 
> | 
      } | 
| 856 | 
> | 
#ifdef IS_MPI | 
| 857 | 
> | 
      tempEx = exI; | 
| 858 | 
> | 
      exI = the_atoms[tempEx]->getGlobalIndex() + 1; | 
| 859 | 
> | 
      tempEx = exJ; | 
| 860 | 
> | 
      exJ = the_atoms[tempEx]->getGlobalIndex() + 1; | 
| 861 | 
> | 
       | 
| 862 | 
> | 
      the_excludes[j+excludeOffset]->setPair( exI, exJ ); | 
| 863 | 
> | 
#else  // isn't MPI | 
| 864 | 
> | 
      the_excludes[j+excludeOffset]->setPair( (exI+1), (exJ+1) ); | 
| 865 | 
> | 
#endif  //is_mpi | 
| 866 | 
> | 
    } | 
| 867 | 
> | 
    excludeOffset += info.nBends; | 
| 868 | 
  | 
 | 
| 869 | 
< | 
#ifdef IS_MPI | 
| 870 | 
< | 
  sprintf( checkPointMsg, | 
| 871 | 
< | 
           "Successfully created the bends list.\n" ); | 
| 872 | 
< | 
  MPIcheckPoint(); | 
| 873 | 
< | 
#endif // is_mpi | 
| 874 | 
< | 
   | 
| 869 | 
> | 
    for(j=0; j<info.nTorsions; j++){ | 
| 870 | 
> | 
       | 
| 871 | 
> | 
      currentTorsion = comp_stamps[stampID]->getTorsion( j ); | 
| 872 | 
> | 
      theTorsions[j].a = currentTorsion->getA() + atomOffset; | 
| 873 | 
> | 
      theTorsions[j].b = currentTorsion->getB() + atomOffset; | 
| 874 | 
> | 
      theTorsions[j].c = currentTorsion->getC() + atomOffset; | 
| 875 | 
> | 
      theTorsions[j].d = currentTorsion->getD() + atomOffset; | 
| 876 | 
> | 
       | 
| 877 | 
> | 
      exI = theTorsions[j].a; | 
| 878 | 
> | 
      exJ = theTorsions[j].d; | 
| 879 | 
  | 
 | 
| 880 | 
< | 
  the_ff->initializeBends( the_bends ); | 
| 881 | 
< | 
} | 
| 882 | 
< | 
 | 
| 883 | 
< | 
void SimSetup::makeTorsions( void ){ | 
| 884 | 
< | 
 | 
| 885 | 
< | 
  int i, j, k, index, offset, molIndex, exI, exJ, tempEx; | 
| 1083 | 
< | 
  torsion_set* the_torsions; | 
| 1084 | 
< | 
  TorsionStamp* current_torsion; | 
| 1085 | 
< | 
 | 
| 1086 | 
< | 
  the_torsions = new torsion_set[tot_torsions]; | 
| 1087 | 
< | 
  index = 0; | 
| 1088 | 
< | 
  offset = 0; | 
| 1089 | 
< | 
  molIndex = 0; | 
| 1090 | 
< | 
  for( i=0; i<n_components; i++ ){ | 
| 1091 | 
< | 
 | 
| 1092 | 
< | 
    for( j=0; j<components_nmol[i]; j++ ){ | 
| 1093 | 
< | 
 | 
| 880 | 
> | 
      // exclude_I must always be the smaller of the pair | 
| 881 | 
> | 
      if( exI > exJ ){ | 
| 882 | 
> | 
        tempEx = exI; | 
| 883 | 
> | 
        exI = exJ; | 
| 884 | 
> | 
        exJ = tempEx; | 
| 885 | 
> | 
      } | 
| 886 | 
  | 
#ifdef IS_MPI | 
| 887 | 
< | 
      if( mpiSim->getMyMolStart() <= molIndex && | 
| 888 | 
< | 
          molIndex <= mpiSim->getMyMolEnd() ){ | 
| 889 | 
< | 
#endif // is_mpi         | 
| 890 | 
< | 
 | 
| 891 | 
< | 
      for( k=0; k<comp_stamps[i]->getNTorsions(); k++ ){ | 
| 892 | 
< | 
 | 
| 1101 | 
< | 
        current_torsion = comp_stamps[i]->getTorsion( k ); | 
| 1102 | 
< | 
        the_torsions[index].a = current_torsion->getA() + offset; | 
| 1103 | 
< | 
        the_torsions[index].b = current_torsion->getB() + offset; | 
| 1104 | 
< | 
        the_torsions[index].c = current_torsion->getC() + offset; | 
| 1105 | 
< | 
        the_torsions[index].d = current_torsion->getD() + offset; | 
| 1106 | 
< | 
 | 
| 1107 | 
< | 
        exI = the_torsions[index].a; | 
| 1108 | 
< | 
        exJ = the_torsions[index].d; | 
| 1109 | 
< | 
 | 
| 1110 | 
< | 
         | 
| 1111 | 
< | 
        // exclude_I must always be the smaller of the pair | 
| 1112 | 
< | 
        if( exI > exJ ){ | 
| 1113 | 
< | 
          tempEx = exI; | 
| 1114 | 
< | 
          exI = exJ; | 
| 1115 | 
< | 
          exJ = tempEx; | 
| 1116 | 
< | 
        } | 
| 1117 | 
< | 
 | 
| 1118 | 
< | 
 | 
| 1119 | 
< | 
#ifdef IS_MPI | 
| 1120 | 
< | 
         | 
| 1121 | 
< | 
        the_excludes[(index + tot_bonds + tot_bends)*2] =      | 
| 1122 | 
< | 
          the_atoms[exI]->getGlobalIndex() + 1; | 
| 1123 | 
< | 
        the_excludes[(index + tot_bonds + tot_bends)*2 + 1] =  | 
| 1124 | 
< | 
          the_atoms[exJ]->getGlobalIndex() + 1; | 
| 1125 | 
< | 
         | 
| 887 | 
> | 
      tempEx = exI; | 
| 888 | 
> | 
      exI = the_atoms[tempEx]->getGlobalIndex() + 1; | 
| 889 | 
> | 
      tempEx = exJ; | 
| 890 | 
> | 
      exJ = the_atoms[tempEx]->getGlobalIndex() + 1; | 
| 891 | 
> | 
       | 
| 892 | 
> | 
      the_excludes[j+excludeOffset]->setPair( exI, exJ ); | 
| 893 | 
  | 
#else  // isn't MPI | 
| 894 | 
< | 
         | 
| 1128 | 
< | 
        the_excludes[(index + tot_bonds + tot_bends)*2] =     exI + 1;  | 
| 1129 | 
< | 
        the_excludes[(index + tot_bonds + tot_bends)*2 + 1] = exJ + 1; | 
| 1130 | 
< | 
        // fortran indexes from 1 (hence the +1 in the indexing) | 
| 894 | 
> | 
      the_excludes[j+excludeOffset]->setPair( (exI+1), (exJ+1) ); | 
| 895 | 
  | 
#endif  //is_mpi | 
| 896 | 
< | 
         | 
| 896 | 
> | 
    } | 
| 897 | 
> | 
    excludeOffset += info.nTorsions; | 
| 898 | 
  | 
 | 
| 899 | 
< | 
        // increment the index and repeat; | 
| 900 | 
< | 
        index++; | 
| 1136 | 
< | 
      } | 
| 1137 | 
< | 
      offset += comp_stamps[i]->getNAtoms(); | 
| 899 | 
> | 
     | 
| 900 | 
> | 
    // send the arrays off to the forceField for init. | 
| 901 | 
  | 
 | 
| 902 | 
< | 
#ifdef IS_MPI | 
| 903 | 
< | 
      } | 
| 904 | 
< | 
#endif //is_mpi       | 
| 902 | 
> | 
    the_ff->initializeAtoms( info.nAtoms, info.myAtoms ); | 
| 903 | 
> | 
    the_ff->initializeBonds( info.nBonds, info.myBonds, theBonds ); | 
| 904 | 
> | 
    the_ff->initializeBends( info.nBends, info.myBends, theBends ); | 
| 905 | 
> | 
    the_ff->initializeTorsions( info.nTorsions, info.myTorsions, theTorsions ); | 
| 906 | 
  | 
 | 
| 907 | 
< | 
      molIndex++; | 
| 908 | 
< | 
    } | 
| 907 | 
> | 
 | 
| 908 | 
> | 
    the_molecules[i].initialize( info ); | 
| 909 | 
> | 
    atomOffset += info.nAtoms; | 
| 910 | 
  | 
  } | 
| 911 | 
  | 
 | 
| 912 | 
< | 
  the_ff->initializeTorsions( the_torsions ); | 
| 912 | 
> | 
  // clean up the forcefield | 
| 913 | 
> | 
  the_ff->calcRcut(); | 
| 914 | 
> | 
  the_ff->cleanMe(); | 
| 915 | 
  | 
} | 
| 916 | 
  | 
 | 
| 917 | 
  | 
void SimSetup::initFromBass( void ){ |