ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/branches/development/src/visitors/AtomVisitor.cpp
(Generate patch)

Comparing trunk/src/visitors/AtomVisitor.cpp (file contents):
Revision 507 by gezelter, Fri Apr 15 22:04:00 2005 UTC vs.
Revision 1008 by chrisfen, Wed Jul 19 12:35:31 2006 UTC

# Line 82 | Line 82 | namespace oopse {
82    void SSDAtomVisitor::visit(DirectionalAtom *datom) {
83      std::vector<AtomInfo*>atoms;
84  
85 <    //we need to convert SSD into 4 differnet atoms
86 <    //one oxygen atom, two hydrogen atoms and one pseudo atom which is the center of the mass
87 <    //of the water with a dipole moment
85 >    //we need to convert SSD into 4 different atoms
86 >    //one oxygen atom, two hydrogen atoms and one pseudo atom which is the center of
87 >    //the mass of the water with a dipole moment
88      Vector3d h1(0.0, -0.75695, 0.5206);
89      Vector3d h2(0.0, 0.75695, 0.5206);
90      Vector3d ox(0.0, 0.0, -0.0654);
# Line 214 | Line 214 | namespace oopse {
214      return result;
215    }
216  
217 +
218 +  bool TREDAtomVisitor::isTREDAtom(const std::string&atomType) {
219 +    std::set<std::string>::iterator strIter;
220 +    strIter = tredAtomType.find(atomType);
221 +    return strIter != tredAtomType.end() ? true : false;
222 +  }
223 +
224 +  void TREDAtomVisitor::visit(DirectionalAtom *datom) {
225 +    std::vector<AtomInfo*>atoms;
226 +
227 +    // we need to convert a TRED into 4 different atoms:
228 +    // one oxygen atom, two hydrogen atoms, and one atom which is the center of
229 +    // the mass of the water with a dipole moment
230 +    Vector3d h1(0.0, -0.75695, 0.5206);
231 +    Vector3d h2(0.0, 0.75695, 0.5206);
232 +    Vector3d ox(0.0, 0.0, -0.0654);
233 +    Vector3d u(0, 0, 1);
234 +    RotMat3x3d   rotMatrix;
235 +    RotMat3x3d   rotTrans;
236 +    AtomInfo *   atomInfo;
237 +    Vector3d     pos;
238 +    Vector3d     newVec;
239 +    Quat4d       q;
240 +    AtomData *   atomData;
241 +    GenericData *data;
242 +    bool         haveAtomData;
243 +
244 +    // if the atom is not a TRED atom, skip it
245 +    if (!isTREDAtom(datom->getType()))
246 +      return;
247 +
248 +    data = datom->getPropertyByName("ATOMDATA");
249 +
250 +    if (data != NULL) {
251 +      atomData = dynamic_cast<AtomData *>(data);
252 +
253 +      if (atomData == NULL) {
254 +        std::cerr << "can not get Atom Data from " << datom->getType() << std::endl;
255 +        atomData = new AtomData;
256 +        haveAtomData = false;
257 +      } else
258 +        haveAtomData = true;
259 +    } else {
260 +      atomData = new AtomData;
261 +      haveAtomData = false;
262 +    }
263 +
264 +    pos = datom->getPos();
265 +    q = datom->getQ();
266 +    rotMatrix = datom->getA();
267 +
268 +    // We need A^T to convert from body-fixed to space-fixed:
269 +    // transposeMat3(rotMatrix, rotTrans);
270 +    rotTrans = rotMatrix.transpose();
271 +
272 +    // center of mass of the water molecule
273 +    // matVecMul3(rotTrans, u, newVec);
274 +    newVec = rotTrans * u;
275 +
276 +    atomInfo = new AtomInfo;
277 +    atomInfo->atomTypeName = "TRED";
278 +    atomInfo->pos[0] = pos[0];
279 +    atomInfo->pos[1] = pos[1];
280 +    atomInfo->pos[2] = pos[2];
281 +    atomInfo->dipole[0] = newVec[0];
282 +    atomInfo->dipole[1] = newVec[1];
283 +    atomInfo->dipole[2] = newVec[2];
284 +
285 +    atomData->addAtomInfo(atomInfo);
286 +
287 +    // oxygen
288 +    // matVecMul3(rotTrans, ox, newVec);
289 +    newVec = rotTrans * ox;
290 +
291 +    atomInfo = new AtomInfo;
292 +    atomInfo->atomTypeName = "O";
293 +    atomInfo->pos[0] = pos[0] + newVec[0];
294 +    atomInfo->pos[1] = pos[1] + newVec[1];
295 +    atomInfo->pos[2] = pos[2] + newVec[2];
296 +    atomInfo->dipole[0] = 0.0;
297 +    atomInfo->dipole[1] = 0.0;
298 +    atomInfo->dipole[2] = 0.0;
299 +    atomData->addAtomInfo(atomInfo);
300 +
301 +    // hydrogen1
302 +    // matVecMul3(rotTrans, h1, newVec);
303 +    newVec = rotTrans * h1;
304 +    atomInfo = new AtomInfo;
305 +    atomInfo->atomTypeName = "H";
306 +    atomInfo->pos[0] = pos[0] + newVec[0];
307 +    atomInfo->pos[1] = pos[1] + newVec[1];
308 +    atomInfo->pos[2] = pos[2] + newVec[2];
309 +    atomInfo->dipole[0] = 0.0;
310 +    atomInfo->dipole[1] = 0.0;
311 +    atomInfo->dipole[2] = 0.0;
312 +    atomData->addAtomInfo(atomInfo);
313 +
314 +    // hydrogen2
315 +    // matVecMul3(rotTrans, h2, newVec);
316 +    newVec = rotTrans * h2;
317 +    atomInfo = new AtomInfo;
318 +    atomInfo->atomTypeName = "H";
319 +    atomInfo->pos[0] = pos[0] + newVec[0];
320 +    atomInfo->pos[1] = pos[1] + newVec[1];
321 +    atomInfo->pos[2] = pos[2] + newVec[2];
322 +    atomInfo->dipole[0] = 0.0;
323 +    atomInfo->dipole[1] = 0.0;
324 +    atomInfo->dipole[2] = 0.0;
325 +    atomData->addAtomInfo(atomInfo);
326 +
327 +    // add atom data into atom's property
328 +
329 +    if (!haveAtomData) {
330 +      atomData->setID("ATOMDATA");
331 +      datom->addProperty(atomData);
332 +    }
333 +
334 +    setVisited(datom);
335 +  }
336 +
337 +  const std::string TREDAtomVisitor::toString() {
338 +    char   buffer[65535];
339 +    std::string result;
340 +
341 +    sprintf(buffer,
342 +            "------------------------------------------------------------------\n");
343 +    result += buffer;
344 +
345 +    sprintf(buffer, "Visitor name: %s\n", visitorName.c_str());
346 +    result += buffer;
347 +
348 +    sprintf(buffer,
349 +            "Visitor Description: Convert the TRED atom into 4 different atoms\n");
350 +    result += buffer;
351 +
352 +    sprintf(buffer,
353 +            "------------------------------------------------------------------\n");
354 +    result += buffer;
355 +
356 +    return result;
357 +  }
358 +
359 +
360    bool LinearAtomVisitor::isLinearAtom(const std::string& atomType){
361      std::set<std::string>::iterator strIter;
362      strIter = linearAtomType.find(atomType);
# Line 221 | Line 364 | namespace oopse {
364      return strIter != linearAtomType.end() ? true : false;
365    }
366  
367 +  void LinearAtomVisitor::addGayBerneAtomType(const std::string& atomType){
368 +   linearAtomType.insert(atomType);
369 +  }
370 +
371    void LinearAtomVisitor::visit(DirectionalAtom* datom){
372      std::vector<AtomInfo*> atoms;
373      //we need to convert linear into 4 different atoms
# Line 237 | Line 384 | namespace oopse {
384      AtomData* atomData;
385      GenericData* data;
386      bool haveAtomData;
387 <
388 <    //if atom is not SSD atom, just skip it
389 <    if(!isLinearAtom(datom->getType()))
387 >    AtomType* atomType;
388 >    //if atom is not linear atom, just skip it
389 >    if(!isLinearAtom(datom->getType()) || !datom->getAtomType()->isGayBerne())
390        return;
391  
392 +    //setup GayBerne type in fortran side
393 +    data = datom->getAtomType()->getPropertyByName("GayBerne");
394 +    if (data != NULL) {
395 +       GayBerneParamGenericData* gayBerneData = dynamic_cast<GayBerneParamGenericData*>(data);
396 +
397 +       if (gayBerneData != NULL) {
398 +           GayBerneParam gayBerneParam = gayBerneData->getData();
399 +
400 +                          // double halfLen = gayBerneParam.GB_sigma * gayBerneParam.GB_l2b_ratio/2.0;
401 +                          double halfLen = gayBerneParam.GB_l/2.0;
402 +                          c1[2] = -halfLen;
403 +              c2[2] = -halfLen /2;
404 +              c3[2] = halfLen/2;
405 +              c4[2] = halfLen;
406 +                
407 +            }
408 +            
409 +              else {
410 +                    sprintf( painCave.errMsg,
411 +                           "Can not cast GenericData to GayBerneParam\n");
412 +                    painCave.severity = OOPSE_ERROR;
413 +                    painCave.isFatal = 1;
414 +                    simError();          
415 +        }            
416 +    }
417 +
418 +
419      data = datom->getPropertyByName("ATOMDATA");
420      if(data != NULL){
421        atomData = dynamic_cast<AtomData*>(data);  
# Line 339 | Line 513 | namespace oopse {
513      return result;
514    }
515  
516 +  bool GBLipidAtomVisitor::isGBLipidAtom(const std::string& atomType){
517 +    std::set<std::string>::iterator strIter;
518 +    strIter = GBLipidAtomType.find(atomType);
519 +
520 +    return strIter != GBLipidAtomType.end() ? true : false;
521 +  }
522 +
523 +  void GBLipidAtomVisitor::visit(DirectionalAtom* datom){
524 +    std::vector<AtomInfo*> atoms;
525 +    //we need to convert linear into 4 different atoms
526 +    Vector3d c1(0.0, 0.0, -6.25);
527 +    Vector3d c2(0.0, 0.0, -2.1);
528 +    Vector3d c3(0.0, 0.0,  2.1);
529 +    Vector3d c4(0.0, 0.0,  6.25);
530 +    RotMat3x3d rotMatrix;
531 +    RotMat3x3d rotTrans;
532 +    AtomInfo* atomInfo;
533 +    Vector3d pos;
534 +    Vector3d newVec;
535 +    Quat4d q;
536 +    AtomData* atomData;
537 +    GenericData* data;
538 +    bool haveAtomData;
539 +
540 +    //if atom is not GBlipid atom, just skip it
541 +    if(!isGBLipidAtom(datom->getType()))
542 +      return;
543 +
544 +    data = datom->getPropertyByName("ATOMDATA");
545 +    if(data != NULL){
546 +      atomData = dynamic_cast<AtomData*>(data);  
547 +      if(atomData == NULL){
548 +        std::cerr << "can not get Atom Data from " << datom->getType() << std::endl;
549 +        atomData = new AtomData;
550 +        haveAtomData = false;      
551 +      } else {
552 +        haveAtomData = true;
553 +      }
554 +    } else {
555 +      atomData = new AtomData;
556 +      haveAtomData = false;
557 +    }
558 +  
559 +  
560 +    pos = datom->getPos();
561 +    q = datom->getQ();
562 +    rotMatrix = datom->getA();
563 +
564 +    // We need A^T to convert from body-fixed to space-fixed:  
565 +    rotTrans = rotMatrix.transpose();
566 +
567 +    newVec = rotTrans * c1;
568 +    atomInfo = new AtomInfo;
569 +    atomInfo->atomTypeName = "K";
570 +    atomInfo->pos[0] = pos[0] + newVec[0];
571 +    atomInfo->pos[1] = pos[1] + newVec[1];
572 +    atomInfo->pos[2] = pos[2] + newVec[2];
573 +    atomInfo->dipole[0] = 0.0;
574 +    atomInfo->dipole[1] = 0.0;
575 +    atomInfo->dipole[2] = 0.0;
576 +    atomData->addAtomInfo(atomInfo);
577 +
578 +    newVec = rotTrans * c2;
579 +    atomInfo = new AtomInfo;
580 +    atomInfo->atomTypeName = "K";
581 +    atomInfo->pos[0] = pos[0] + newVec[0];
582 +    atomInfo->pos[1] = pos[1] + newVec[1];
583 +    atomInfo->pos[2] = pos[2] + newVec[2];
584 +    atomInfo->dipole[0] = 0.0;
585 +    atomInfo->dipole[1] = 0.0;
586 +    atomInfo->dipole[2] = 0.0;
587 +    atomData->addAtomInfo(atomInfo);
588 +
589 +    newVec = rotTrans * c3;
590 +    atomInfo = new AtomInfo;
591 +    atomInfo->atomTypeName = "K";
592 +    atomInfo->pos[0] = pos[0] + newVec[0];
593 +    atomInfo->pos[1] = pos[1] + newVec[1];
594 +    atomInfo->pos[2] = pos[2] + newVec[2];
595 +    atomInfo->dipole[0] = 0.0;
596 +    atomInfo->dipole[1] = 0.0;
597 +    atomInfo->dipole[2] = 0.0;
598 +    atomData->addAtomInfo(atomInfo);
599 +
600 +    newVec = rotTrans * c4;
601 +    atomInfo = new AtomInfo;
602 +    atomInfo->atomTypeName = "K";
603 +    atomInfo->pos[0] = pos[0] + newVec[0];
604 +    atomInfo->pos[1] = pos[1] + newVec[1];
605 +    atomInfo->pos[2] = pos[2] + newVec[2];
606 +    atomInfo->dipole[0] = 0.0;
607 +    atomInfo->dipole[1] = 0.0;
608 +    atomInfo->dipole[2] = 0.0;
609 +    atomData->addAtomInfo(atomInfo);
610 +
611 +    //add atom data into atom's property
612 +
613 +    if(!haveAtomData){
614 +      atomData->setID("ATOMDATA");
615 +      datom->addProperty(atomData);
616 +    }
617 +
618 +    setVisited(datom);
619 +
620 +  }
621 +
622 +  const std::string GBLipidAtomVisitor::toString(){
623 +    char buffer[65535];
624 +    std::string result;
625 +  
626 +    sprintf(buffer ,"------------------------------------------------------------------\n");
627 +    result += buffer;
628 +
629 +    sprintf(buffer ,"Visitor name: %s\n", visitorName.c_str());
630 +    result += buffer;
631 +
632 +    sprintf(buffer , "Visitor Description: Convert GBlipid into 4 different K atoms\n");
633 +    result += buffer;
634 +
635 +    sprintf(buffer ,"------------------------------------------------------------------\n");
636 +    result += buffer;
637 +
638 +    return result;
639 +  }
640 +
641    //----------------------------------------------------------------------------//
642  
643    void DefaultAtomVisitor::visit(Atom *atom) {
# Line 380 | Line 679 | namespace oopse {
679        return;
680  
681      pos = datom->getPos();
682 <    u = datom->getElectroFrame().getColumn(2);
683 <
682 >    if (datom->getAtomType()->isGayBerne()) {
683 >        u = datom->getA().transpose()*V3Z;        
684 >    } else if (datom->getAtomType()->isMultipole()) {
685 >        u = datom->getElectroFrame().getColumn(2);
686 >    }
687      atomData = new AtomData;
688      atomData->setID("ATOMDATA");
689      atomInfo = new AtomInfo;

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines