| 30 | 
  | 
} | 
| 31 | 
  | 
 | 
| 32 | 
  | 
   | 
| 33 | 
< | 
Verlet::Verlet( SimInfo &info ){ | 
| 33 | 
> | 
Verlet::Verlet( SimInfo &info, ForceField* the_ff ){ | 
| 34 | 
  | 
   | 
| 35 | 
  | 
  // get what information we need from the SimInfo object | 
| 36 | 
  | 
 | 
| 37 | 
  | 
  entry_plug = &info; | 
| 38 | 
+ | 
  myFF = the_ff; | 
| 39 | 
  | 
 | 
| 40 | 
+ | 
 | 
| 41 | 
  | 
  c_natoms = info.n_atoms; | 
| 42 | 
  | 
  c_atoms = info.atoms; | 
| 43 | 
  | 
  c_sr_interactions = info.sr_interactions; | 
| 122 | 
  | 
   | 
| 123 | 
  | 
  delete[] c_mass; | 
| 124 | 
  | 
  c_mass = 0; | 
| 123 | 
– | 
} | 
| 124 | 
– | 
 | 
| 125 | 
– | 
 | 
| 126 | 
– | 
 | 
| 127 | 
– | 
void Verlet::integrate_b( double time_length, double dt,  | 
| 128 | 
– | 
                            int n_bond_0, int n_bond_f, | 
| 129 | 
– | 
                            int n_bend_0, int n_bend_f, | 
| 130 | 
– | 
                            int n_torsion_0, int n_torsion_f, | 
| 131 | 
– | 
                            bool do_bonds, bool do_bends, bool do_torsions, | 
| 132 | 
– | 
                            bool do_LRI ){ | 
| 133 | 
– | 
   | 
| 134 | 
– | 
//  double percent_tolerance = 0.001; | 
| 135 | 
– | 
//  int max_iterations = 10000; | 
| 136 | 
– | 
 | 
| 137 | 
– | 
  int i, j; /* loop counters */ | 
| 138 | 
– | 
  double n_loops = time_length / dt; | 
| 139 | 
– | 
   | 
| 140 | 
– | 
  // the first time integrate is called, the forces need to be initialized  | 
| 141 | 
– | 
 | 
| 142 | 
– | 
  if(is_first){  | 
| 143 | 
– | 
    is_first = 0; | 
| 144 | 
– | 
     | 
| 145 | 
– | 
    for(i = 0; i < c_natoms; i++){ | 
| 146 | 
– | 
      c_atoms[i]->zeroForces(); | 
| 147 | 
– | 
    } | 
| 148 | 
– | 
     | 
| 149 | 
– | 
    if( do_bonds ){ | 
| 150 | 
– | 
      for(i = n_bond_0; i <= n_bond_f; i++){ | 
| 151 | 
– | 
        c_sr_interactions[i]->calc_forces(); | 
| 152 | 
– | 
      } | 
| 153 | 
– | 
    } | 
| 154 | 
– | 
 | 
| 155 | 
– | 
    if( do_bends ){ | 
| 156 | 
– | 
      for(i = n_bend_0; i <= n_bend_f; i++){ | 
| 157 | 
– | 
        c_sr_interactions[i]->calc_forces(); | 
| 158 | 
– | 
      } | 
| 159 | 
– | 
    } | 
| 160 | 
– | 
     | 
| 161 | 
– | 
    if( do_torsions ){ | 
| 162 | 
– | 
      for(i = n_torsion_0; i <= n_torsion_f; i++){ | 
| 163 | 
– | 
        c_sr_interactions[i]->calc_forces(); | 
| 164 | 
– | 
      } | 
| 165 | 
– | 
    } | 
| 166 | 
– | 
  | 
| 167 | 
– | 
    if( do_LRI ) longRange->calc_forces(); | 
| 168 | 
– | 
  } | 
| 169 | 
– | 
 | 
| 170 | 
– | 
  for(i = 0; i < n_loops; i++){ | 
| 171 | 
– | 
     | 
| 172 | 
– | 
    move_a( dt ); | 
| 173 | 
– | 
     | 
| 174 | 
– | 
    // calculate the forces | 
| 175 | 
– | 
     | 
| 176 | 
– | 
    for(j = 0; j < c_natoms; j++){ | 
| 177 | 
– | 
      c_atoms[j]->zeroForces(); | 
| 178 | 
– | 
    } | 
| 179 | 
– | 
   | 
| 180 | 
– | 
 | 
| 181 | 
– | 
    if( do_bonds ){ | 
| 182 | 
– | 
      for(i = n_bond_0; i <= n_bond_f; i++){ | 
| 183 | 
– | 
        c_sr_interactions[i]->calc_forces(); | 
| 184 | 
– | 
      } | 
| 185 | 
– | 
    } | 
| 186 | 
– | 
 | 
| 187 | 
– | 
    if( do_bends ){ | 
| 188 | 
– | 
      for(i = n_bend_0; i <= n_bend_f; i++){ | 
| 189 | 
– | 
        c_sr_interactions[i]->calc_forces(); | 
| 190 | 
– | 
      } | 
| 191 | 
– | 
    } | 
| 192 | 
– | 
     | 
| 193 | 
– | 
    if( do_torsions ){ | 
| 194 | 
– | 
      for(i = n_torsion_0; i <= n_torsion_f; i++){ | 
| 195 | 
– | 
        c_sr_interactions[i]->calc_forces(); | 
| 196 | 
– | 
      } | 
| 197 | 
– | 
    } | 
| 198 | 
– | 
 | 
| 199 | 
– | 
    if( do_LRI ) longRange->calc_forces(); | 
| 200 | 
– | 
     | 
| 201 | 
– | 
 | 
| 202 | 
– | 
    // complete the verlet move | 
| 203 | 
– | 
 | 
| 204 | 
– | 
    move_b( dt ); | 
| 205 | 
– | 
  } | 
| 125 | 
  | 
} | 
| 126 | 
  | 
 | 
| 127 | 
  | 
 | 
| 164 | 
  | 
  // the first time integrate is called, the forces need to be initialized  | 
| 165 | 
  | 
 | 
| 166 | 
  | 
   | 
| 167 | 
< | 
  for(i = 0; i < c_natoms; i++){ | 
| 249 | 
< | 
    c_atoms[i]->zeroForces(); | 
| 250 | 
< | 
  } | 
| 167 | 
> | 
  myFF->doForces(); | 
| 168 | 
  | 
   | 
| 252 | 
– | 
  for(i = 0; i < c_n_SRI; i++){ | 
| 253 | 
– | 
    c_sr_interactions[i]->calc_forces(); | 
| 254 | 
– | 
  } | 
| 255 | 
– | 
   | 
| 256 | 
– | 
  longRange->calc_forces(); | 
| 257 | 
– | 
   | 
| 169 | 
  | 
  if( entry_plug->setTemp ){ | 
| 170 | 
  | 
    tStats->velocitize(); | 
| 171 | 
  | 
  } | 
| 213 | 
  | 
 | 
| 214 | 
  | 
      // calculate the forces | 
| 215 | 
  | 
       | 
| 216 | 
< | 
      for(j = 0; j < c_natoms; j++){ | 
| 306 | 
< | 
        c_atoms[j]->zeroForces(); | 
| 307 | 
< | 
      } | 
| 216 | 
> | 
      myFF->doForces(); | 
| 217 | 
  | 
       | 
| 309 | 
– | 
      for(j = 0; j < c_n_SRI; j++){ | 
| 310 | 
– | 
        c_sr_interactions[j]->calc_forces(); | 
| 311 | 
– | 
      } | 
| 312 | 
– | 
       | 
| 313 | 
– | 
      longRange->calc_forces(); | 
| 314 | 
– | 
       | 
| 218 | 
  | 
      // finish the constrain move ( same as above. ) | 
| 219 | 
  | 
 | 
| 220 | 
  | 
      for( j=0; j<c_natoms; j++ ){ | 
| 265 | 
  | 
       | 
| 266 | 
  | 
      // calculate the forces | 
| 267 | 
  | 
       | 
| 268 | 
< | 
      for(j = 0; j < c_natoms; j++){ | 
| 366 | 
< | 
        c_atoms[j]->zeroForces(); | 
| 367 | 
< | 
      } | 
| 368 | 
< | 
       | 
| 369 | 
< | 
      for(j = 0; j < c_n_SRI; j++){ | 
| 370 | 
< | 
        c_sr_interactions[j]->calc_forces(); | 
| 371 | 
< | 
      } | 
| 372 | 
< | 
       | 
| 373 | 
< | 
      longRange->calc_forces(); | 
| 268 | 
> | 
      myFF->doForces(); | 
| 269 | 
  | 
             | 
| 270 | 
  | 
      // complete the verlet move | 
| 271 | 
  | 
       | 
| 330 | 
  | 
  const double e_convert = 4.184e-4; // converts kcal/mol -> amu*A^2/fs^2 | 
| 331 | 
  | 
   | 
| 332 | 
  | 
  double vx, vy, vz; | 
| 438 | 
– | 
  double v_sqr; | 
| 333 | 
  | 
  int mb; | 
| 334 | 
  | 
  double h_dt = 0.5 * dt; | 
| 335 | 
  | 
   |