ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE-4/src/integrators/SMIPDForceManager.cpp
(Generate patch)

Comparing trunk/OOPSE-4/src/integrators/SMIPDForceManager.cpp (file contents):
Revision 3461 by chuckv, Wed Oct 15 18:26:01 2008 UTC vs.
Revision 3463 by chuckv, Thu Oct 16 18:25:36 2008 UTC

# Line 218 | Line 218 | namespace oopse {
218                }
219              }
220            }
221 <          HydroProp* currHydroProp = currShape->getHydroProp(1.0,simParams->getTargetTemp());
221 >          HydroProp* currHydroProp = currShape->getHydroProp(simParams->getViscosity(),simParams->getTargetTemp());
222            std::map<std::string, HydroProp*>::iterator iter = hydroPropMap.find(integrableObject->getType());
223            if (iter != hydroPropMap.end())
224              hydroProps_.push_back(iter->second);
# Line 243 | Line 243 | namespace oopse {
243  
244      surfaceMesh_->computeHull(localSites_);
245      Area0_ = surfaceMesh_->getArea();
246 <    //variance_ = 2.0 * OOPSEConstant::kb*simParams->getTargetTemp()/simParams->getDt();
246 >    variance_ = 2.0 * OOPSEConstant::kb*simParams->getTargetTemp()/simParams->getDt();
247      
248    }  
249  
# Line 290 | Line 290 | namespace oopse {
290  
291       /* Compute variance for random forces */
292      
293 <     variance_ = sqrt(2.0*NumericConstant::PI)*((targetPressure_/OOPSEConstant::pressureConvert)*area/nSurfaceSDs)
293 >     RealType TD_variance = sqrt(2.0*NumericConstant::PI)*((targetPressure_/OOPSEConstant::pressureConvert)*area/nSurfaceSDs)
294         /OOPSEConstant::energyConvert;
295      
296      std::vector<Triangle*> sMesh = surfaceMesh_->getMesh();
297 <    std::vector<RealType>  randNums = genTriangleForces(sMesh.size(),variance_);
297 >    std::vector<RealType>  randNums = genTriangleForces(sMesh.size(),TD_variance);
298      
299      /* Loop over the mesh faces and apply random force to each of the faces*/
300      
# Line 340 | Line 340 | namespace oopse {
340          A = integrableObject->getA();
341          Atrans = A.transpose();
342          Vector3d rcrLab = Atrans * hydroProps_[index]->getCOR();  
343 +        //apply random force and torque at center of resistance
344 +
345 +        Vector3d randomForceBody;
346 +        Vector3d randomTorqueBody;
347 +        genRandomForceAndTorque(randomForceBody, randomTorqueBody, index, variance_);
348 +        Vector3d randomForceLab = Atrans * randomForceBody;
349 +        Vector3d randomTorqueLab = Atrans * randomTorqueBody;
350 +        integrableObject->addFrc(randomForceLab);            
351 +        integrableObject->addTrq(randomTorqueLab + cross(rcrLab, randomForceLab ));            
352 +        
353  
354          
355          Mat3x3d I = integrableObject->getI();
# Line 426 | Line 436 | namespace oopse {
436        } else {
437          //spherical atom
438  
439 <        
439 >        Vector3d randomForce;
440 >        Vector3d randomTorque;
441 >        genRandomForceAndTorque(randomForce, randomTorque, index, variance_);
442 >        integrableObject->addFrc(randomForce);  
443 >        
444          // What remains contains velocity explicitly, but the velocity required
445          // is at the full step: v(t + h), while we have initially the velocity
446          // at the half step: v(t + h/2).  We need to iterate to converge the
# Line 536 | Line 550 | namespace oopse {
550      return gaussRand;
551    }
552  
553 +
554 +
555 +
556 +
557   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines