| 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); |
| 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 |
|
|
| 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 |
|
|
| 311 |
|
Vector3d unitNormal = thisTriangle->getNormal(); |
| 312 |
|
unitNormal.normalize(); |
| 313 |
|
Vector3d randomForce = -randNums[thisNumber] * unitNormal; |
| 314 |
< |
|
| 314 |
> |
Vector3d centroid = thisTriangle->getCentroid(); |
| 315 |
> |
|
| 316 |
|
for (vertex = vertexSDs.begin(); vertex != vertexSDs.end(); ++vertex){ |
| 317 |
|
|
| 318 |
|
// mass = integrableObject->getMass(); |
| 319 |
< |
|
| 320 |
< |
(*vertex)->addFrc(randomForce/3.0); |
| 319 |
> |
Vector3d vertexForce = randomForce/3.0; |
| 320 |
> |
(*vertex)->addFrc(vertexForce); |
| 321 |
> |
if (integrableObject->isDirectional()){ |
| 322 |
> |
Vector3d vertexPos = (*vertex)->getPos(); |
| 323 |
> |
Vector3d vertexCentroidVector = vertexPos - centroid; |
| 324 |
> |
(*vertex)->addTrq(cross(vertexCentroidVector,vertexForce)); |
| 325 |
> |
} |
| 326 |
> |
|
| 327 |
|
} |
| 328 |
|
} |
| 329 |
|
|
| 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(); |
| 356 |
|
Vector3d omegaBody; |
| 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 |
| 550 |
|
return gaussRand; |
| 551 |
|
} |
| 552 |
|
|
| 553 |
+ |
|
| 554 |
+ |
|
| 555 |
+ |
|
| 556 |
+ |
|
| 557 |
|
} |
