ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/OpenMD/trunk/src/brains/Thermo.cpp
(Generate patch)

Comparing trunk/src/brains/Thermo.cpp (file contents):
Revision 2021 by gezelter, Wed Feb 26 14:14:50 2014 UTC vs.
Revision 2022 by gezelter, Fri Sep 26 22:22:28 2014 UTC

# Line 416 | Line 416 | namespace OpenMD {
416      }
417  
418      return snap->getSystemDipole();
419 +  }
420 +
421 +
422 +  Mat3x3d Thermo::getSystemQuadrupole() {
423 +    Snapshot* snap = info_->getSnapshotManager()->getCurrentSnapshot();
424 +
425 +    if (!snap->hasSystemQuadrupole) {
426 +      SimInfo::MoleculeIterator miter;
427 +      vector<Atom*>::iterator aiter;
428 +      Molecule* mol;
429 +      Atom* atom;
430 +      RealType charge;
431 +      Vector3d ri(0.0);
432 +      Vector3d dipole(0.0);
433 +      Mat3x3d qpole(0.0);
434 +      
435 +      RealType chargeToC = 1.60217733e-19;
436 +      RealType angstromToM = 1.0e-10;
437 +      RealType debyeToCm = 3.33564095198e-30;
438 +      
439 +      for (mol = info_->beginMolecule(miter); mol != NULL;
440 +           mol = info_->nextMolecule(miter)) {
441 +        
442 +        for (atom = mol->beginAtom(aiter); atom != NULL;
443 +             atom = mol->nextAtom(aiter)) {
444 +
445 +          ri = atom->getPos();
446 +          snap->wrapVector(ri);
447 +          ri *= angstromToM;
448 +          
449 +          charge = 0.0;
450 +          
451 +          FixedChargeAdapter fca = FixedChargeAdapter(atom->getAtomType());
452 +          if ( fca.isFixedCharge() ) {
453 +            charge = fca.getCharge();
454 +          }
455 +          
456 +          FluctuatingChargeAdapter fqa = FluctuatingChargeAdapter(atom->getAtomType());
457 +          if ( fqa.isFluctuatingCharge() ) {
458 +            charge += atom->getFlucQPos();
459 +          }
460 +          
461 +          charge *= chargeToC;
462 +          
463 +          qpole += 0.5 * charge * outProduct(ri, ri);
464 +
465 +          MultipoleAdapter ma = MultipoleAdapter(atom->getAtomType());
466 +          
467 +          if ( ma.isDipole() ) {
468 +            dipole = atom->getDipole() * debyeToCm;
469 +            qpole += 0.5 * outProduct( dipole, ri );
470 +          }
471 +
472 +          if ( ma.isQuadrupole() ) {
473 +            qpole += atom->getQuadrupole() * debyeToCm * angstromToM;          
474 +          }
475 +        }
476 +      }
477 +        
478 + #ifdef IS_MPI
479 +      MPI_Allreduce(MPI_IN_PLACE, qpole.getArrayPointer(),
480 +                    9, MPI_REALTYPE, MPI_SUM, MPI_COMM_WORLD);
481 + #endif
482 +      
483 +      snap->setSystemQuadrupole(qpole);
484 +    }
485 +    
486 +    return snap->getSystemQuadrupole();
487    }
488  
489    // Returns the Heat Flux Vector for the system

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines