ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libBASS/MoleculeStamp.cpp
(Generate patch)

Comparing:
branches/mmeineke/OOPSE/libBASS/MoleculeStamp.cpp (file contents), Revision 377 by mmeineke, Fri Mar 21 17:42:12 2003 UTC vs.
trunk/OOPSE/libBASS/MoleculeStamp.cpp (file contents), Revision 1268 by tim, Fri Jun 11 17:16:21 2004 UTC

# Line 1 | Line 1
1 < #include <cstdlib>
2 < #include <cstdio>
3 < #include <cstring>
1 > #include <stdlib.h>
2 > #include <stdio.h>
3 > #include <string.h>
4   #include <iostream>
5  
6   #include "MoleculeStamp.hpp"
# Line 13 | Line 13 | MoleculeStamp::MoleculeStamp(){
13    n_bonds = 0;
14    n_bends = 0;
15    n_torsions = 0;
16 +  n_rigidbodies = 0;
17 +  n_cutoffgroups = 0;
18 +  n_integrable = 0;
19  
20    unhandled = NULL;
21    atoms = NULL;
22    bonds = NULL;
23    bends = NULL;
24    torsions = NULL;
25 +  rigidBodies = NULL;
26 +  cutoffGroups = NULL;
27  
28    have_name = 0;
29    have_atoms = 0;
30    have_bonds = 0;
31    have_bends = 0;
32    have_torsions = 0;
33 +  have_rigidbodies = 0;
34 +  have_cutoffgroups = 0;
35  
36   }
37  
# Line 32 | Line 39 | MoleculeStamp::~MoleculeStamp(){
39    int i;
40    
41    if( unhandled != NULL) delete unhandled;
42 +  
43 +  if( rigidBodies != NULL ) {
44 +    for( i=0; i<n_rigidbodies; i++ ) delete rigidBodies[i];
45 +  }
46 +  delete[] rigidBodies;
47  
48 +  if( cutoffGroups != NULL ) {
49 +    for( i=0; i<n_cutoffgroups; i++ ) delete cutoffGroups[i];
50 +  }
51 +  delete[] cutoffGroups;
52 +  
53    if( atoms != NULL ){
54      for( i=0; i<n_atoms; i++ ) delete atoms[i];
55    }
56 +  delete[] atoms;
57    
58    if( bonds != NULL ){
59      for( i=0; i<n_bonds; i++ ) delete bonds[i];
60    }
61 +  delete[] bonds;
62    
63    if( bends != NULL ){
64      for( i=0; i<n_bends; i++ ) delete bends[i];
65    }
66 +  delete[] bends;
67    
68    if( torsions != NULL ){
69      for( i=0; i<n_torsions; i++ ) delete torsions[i];
70    }
71 +  delete[] torsions;
72    
73 +
74 +
75 +
76   }
77  
78   char* MoleculeStamp::assignString( char* lhs, char* rhs ){
# Line 74 | Line 98 | char* MoleculeStamp::assignDouble( char* lhs, double r
98      if( have_atoms ){
99        sprintf( errMsg,
100                 "MoleculeStamp error, n_atoms already declared"
101 <               "for molecule: %s\n",
101 >               " for molecule: %s\n",
102                 name);
103        return strdup( errMsg );
104      }
# Line 127 | Line 151 | char* MoleculeStamp::assignDouble( char* lhs, double r
151      torsions = new TorsionStamp*[n_torsions];
152      for( i=0; i<n_torsions; i++ ) torsions[i] = NULL;
153    }
154 +
155 +  else if( !strcmp( lhs, "nRigidBodies" ) ){
156 +    n_rigidbodies = (int)rhs;
157 +
158 +    if( have_rigidbodies ){
159 +      sprintf( errMsg,
160 +               "MoleculeStamp error, n_rigidbodies already declared for"
161 +               " molecule: %s\n",
162 +               name );
163 +      return strdup( errMsg );
164 +    }
165 +    have_rigidbodies = 1;
166 +    rigidBodies = new RigidBodyStamp*[n_rigidbodies];
167 +    for( i=0; i<n_rigidbodies; i++ ) rigidBodies[i] = NULL;
168 +  }
169 +
170 +  else if( !strcmp( lhs, "nCutoffGroups" ) ){
171 +    n_cutoffgroups = (int)rhs;
172 +
173 +    if( have_cutoffgroups ){
174 +      sprintf( errMsg,
175 +               "MoleculeStamp error, n_cutoffgroups already declared for"
176 +               " molecule: %s\n",
177 +               name );
178 +      return strdup( errMsg );
179 +    }
180 +    have_cutoffgroups = 1;
181 +    cutoffGroups = new CutoffGroupStamp*[n_cutoffgroups];
182 +    for( i=0; i<n_cutoffgroups; i++ ) cutoffGroups[i] = NULL;
183 +  }
184 +  
185    else{
186      if( unhandled == NULL ) unhandled = new LinkedAssign( lhs, rhs );
187      else unhandled->add( lhs, rhs );
# Line 197 | Line 252 | char*  MoleculeStamp::assignInt( char* lhs, int rhs ){
252      torsions = new TorsionStamp*[n_torsions];
253      for( i=0; i<n_torsions; i++ ) torsions[i] = NULL;
254    }
255 +
256 +  else if( !strcmp( lhs, "nRigidBodies" ) ){
257 +    n_rigidbodies = rhs;
258 +
259 +    if( have_rigidbodies ){
260 +      sprintf( errMsg,
261 +               "MoleculeStamp error, n_rigidbodies already declared for"
262 +               " molecule: %s\n",
263 +               name);
264 +      return strdup( errMsg );
265 +    }
266 +    have_rigidbodies = 1;
267 +    rigidBodies = new RigidBodyStamp*[n_rigidbodies];
268 +    for( i=0; i<n_rigidbodies; i++ ) rigidBodies[i] = NULL;
269 +  }
270 +  else if( !strcmp( lhs, "nCutoffGroups" ) ){
271 +    n_cutoffgroups = rhs;
272 +
273 +    if( have_cutoffgroups ){
274 +      sprintf( errMsg,
275 +               "MoleculeStamp error, n_cutoffgroups already declared for"
276 +               " molecule: %s\n",
277 +               name);
278 +      return strdup( errMsg );
279 +    }
280 +    have_cutoffgroups = 1;
281 +    cutoffGroups = new CutoffGroupStamp*[n_cutoffgroups];
282 +    for( i=0; i<n_cutoffgroups; i++ ) cutoffGroups[i] = NULL;
283 +  }
284    else{
285      if( unhandled == NULL ) unhandled = new LinkedAssign( lhs, rhs );
286      else unhandled->add( lhs, rhs );
# Line 208 | Line 292 | char* MoleculeStamp::addAtom( AtomStamp* the_atom, int
292   char* MoleculeStamp::addAtom( AtomStamp* the_atom, int atomIndex ){
293    
294    if( have_atoms && atomIndex < n_atoms ) atoms[atomIndex] = the_atom;
295 <  else{
295 >  else {
296      if( have_atoms ){
297        sprintf( errMsg, "MoleculeStamp error, %d out of nAtoms range",
298                 atomIndex );
299        return strdup( errMsg );
300      }
301      else return strdup("MoleculeStamp error, nAtoms not given before"
302 <                       "first atom declaration." );
302 >                       " first atom declaration." );
303    }
304  
305    return NULL;
306   }
307  
308 + char* MoleculeStamp::addRigidBody( RigidBodyStamp* the_rigidbody,
309 +                                   int rigidBodyIndex ){
310 +  
311 +  if( have_rigidbodies && rigidBodyIndex < n_rigidbodies )
312 +    rigidBodies[rigidBodyIndex] = the_rigidbody;
313 +  else {
314 +    if( have_rigidbodies ){
315 +      sprintf( errMsg, "MoleculeStamp error, %d out of nRigidBodies range",
316 +               rigidBodyIndex );
317 +      return strdup( errMsg );
318 +    }
319 +    else return strdup("MoleculeStamp error, nRigidBodies not given before"
320 +                       " first rigidBody declaration." );
321 +  }
322 +  
323 +  return NULL;
324 + }
325 +
326 + char* MoleculeStamp::addCutoffGroup( CutoffGroupStamp* the_cutoffgroup,
327 +                                     int cutoffGroupIndex ){
328 +  
329 +  if( have_cutoffgroups && cutoffGroupIndex < n_cutoffgroups )
330 +    cutoffGroups[cutoffGroupIndex] = the_cutoffgroup;
331 +  else {
332 +    if( have_cutoffgroups ){
333 +      sprintf( errMsg, "MoleculeStamp error, %d out of nCutoffGroups range",
334 +               cutoffGroupIndex );
335 +      return strdup( errMsg );
336 +    }
337 +    else return strdup("MoleculeStamp error, nCutoffGroups not given before"
338 +                       " first CutoffGroup declaration." );
339 +  }
340 +  
341 +  return NULL;
342 + }
343 +
344   char* MoleculeStamp::addBond( BondStamp* the_bond, int bondIndex ){
345    
346    
# Line 276 | Line 396 | char* MoleculeStamp::checkMe( void ){
396   char* MoleculeStamp::checkMe( void ){
397    
398    int i;
399 <  short int no_atom;
399 >  short int no_atom, no_rigidbody, no_cutoffgroup;
400 >
401 >  if( !have_name ) return strdup( "MoleculeStamp error. Molecule's name"
402 >                                  " was not given.\n" );
403    
404 <  if( !have_name || !have_atoms ){
405 <    if( !have_name ) return strdup( "MoleculeStamp error. Molecule's name"
283 <                                    " was not given.\n" );
284 <    else return strdup( "MoleculeStamp error. Molecule contains no atoms." );
404 >  if( !have_atoms ){
405 >    return strdup( "MoleculeStamp error. Molecule contains no atoms." );
406    }
407    
408 +  no_rigidbody = 0;
409 +  for( i=0; i<n_rigidbodies; i++ ){
410 +    if( rigidBodies[i] == NULL ) no_rigidbody = 1;
411 +  }
412 +
413 +  if( no_rigidbody ){
414 +    sprintf( errMsg,
415 +             "MoleculeStamp error. Not all of the RigidBodies were"
416 +             " declared in molecule \"%s\".\n", name );
417 +    return strdup( errMsg );
418 +  }
419 +
420 +  no_cutoffgroup = 0;
421 +  for( i=0; i<n_cutoffgroups; i++ ){
422 +    if( cutoffGroups[i] == NULL ) no_cutoffgroup = 1;
423 +  }
424 +
425 +  if( no_cutoffgroup ){
426 +    sprintf( errMsg,
427 +             "MoleculeStamp error. Not all of the CutoffGroups were"
428 +             " declared in molecule \"%s\".\n", name );
429 +    return strdup( errMsg );
430 +  }
431 +  
432    no_atom = 0;
433    for( i=0; i<n_atoms; i++ ){
434      if( atoms[i] == NULL ) no_atom = 1;
# Line 296 | Line 441 | char* MoleculeStamp::checkMe( void ){
441      return strdup( errMsg );
442    }
443  
444 +  n_integrable = n_atoms;
445 +  for (i = 0; i < n_rigidbodies; i++)
446 +    n_integrable = n_integrable - rigidBodies[i]->getNMembers() + 1; //rigidbody is an integrable object
447 +  
448 +  if (n_integrable <= 0 || n_integrable > n_atoms) {
449 +    sprintf( errMsg,
450 +             "MoleculeStamp error. n_integrable is either <= 0 or"
451 +             " greater than n_atoms in molecule \"%s\".\n", name );
452 +    return strdup( errMsg );
453 +  }
454 +
455    return NULL;
456 + }  
457 +
458 +
459 + //Function Name: isBondInSameRigidBody
460 + //Return true is both atoms of the bond belong to the same rigid body, otherwise return false
461 + bool MoleculeStamp::isBondInSameRigidBody(BondStamp* bond){
462 +  int rbA;
463 +  int rbB;
464 +  int consAtomA;
465 +  int consAtomB;
466 +
467 +  if (!isAtomInRigidBody(bond->getA(),rbA, consAtomA))
468 +    return false;
469 +
470 +  if(!isAtomInRigidBody(bond->getB(),rbB, consAtomB) )
471 +    return false;
472 +
473 +  if(rbB == rbA)
474 +    return true;
475 +  else
476 +    return false;
477   }
478 +
479 +
480 + // Function Name isAtomInRigidBody
481 + //return false if atom does not belong to a rigid body otherwise return true and set whichRigidBody
482 + //and consAtomIndex
483 + //atomIndex : the index of atom in component
484 + //whichRigidBody: the index of rigidbody in component
485 + //consAtomIndex:  the position of joint atom apears in  rigidbody's definition
486 + bool MoleculeStamp::isAtomInRigidBody(int atomIndex, int& whichRigidBody, int& consAtomIndex){
487 +  RigidBodyStamp* rbStamp;
488 +  int numRb;
489 +  int numAtom;
490 +
491 +  numRb = this->getNRigidBodies();
492 +  
493 +  for(int i = 0 ; i < numRb; i++){
494 +    rbStamp = this->getRigidBody(i);
495 +    numAtom = rbStamp->getNMembers();
496 +    for(int j = 0; j < numAtom; j++)
497 +      if (rbStamp->getMember(j) == atomIndex){
498 +        whichRigidBody = i;
499 +        consAtomIndex = j;
500 +        return true;
501 +      }
502 +  }
503 +
504 +  return false;
505 +  
506 + }
507 +
508 + //return the position of joint atom apears in  rigidbody's definition
509 + //for the time being, we will use the most inefficient algorithm, the complexity is O(N2)
510 + //actually we could improve the complexity to O(NlgN) by sorting the atom index in rigid body first
511 + vector<pair<int, int> > MoleculeStamp::getJointAtoms(int rb1, int rb2){
512 +  RigidBodyStamp* rbStamp1;
513 +  RigidBodyStamp* rbStamp2;
514 +  int natomInRb1;
515 +  int natomInRb2;
516 +  int atomIndex1;
517 +  int atomIndex2;
518 +  vector<pair<int, int> > jointAtomIndexPair;
519 +  
520 +  rbStamp1 = this->getRigidBody(rb1);
521 +  natomInRb1 =rbStamp1->getNMembers();
522 +
523 +  rbStamp2 = this->getRigidBody(rb2);
524 +  natomInRb2 =rbStamp2->getNMembers();
525 +
526 +  for(int i = 0; i < natomInRb1; i++){
527 +    atomIndex1 = rbStamp1->getMember(i);
528 +      
529 +    for(int j= 0; j < natomInRb1; j++){
530 +      atomIndex2 = rbStamp2->getMember(j);
531 +
532 +      if(atomIndex1 == atomIndex2){
533 +        jointAtomIndexPair.push_back(make_pair(i, j));
534 +        break;
535 +      }
536 +      
537 +    }//end for(j =0)
538 +
539 +  }//end for (i = 0)
540 +
541 +  return jointAtomIndexPair;
542 + }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines