# | Line 190 | Line 190 | namespace OpenMD { | |
---|---|---|
190 | // ObjectRestraints in that method: | |
191 | ||
192 | scanFile(); | |
193 | – | |
193 | ||
194 | readSet(); | |
195 | ||
# | Line 279 | Line 278 | namespace OpenMD { | |
278 | } | |
279 | ||
280 | int index = tokenizer.nextTokenAsInt(); | |
281 | + | |
282 | + | StuntDouble* integrableObject = info_->getIOIndexToIntegrableObject(index); |
283 | ||
284 | + | if (integrableObject == NULL) { |
285 | + | return; |
286 | + | } |
287 | ||
288 | < | if (index < 1116){ |
289 | < | std::string type = tokenizer.nextToken(); |
290 | < | |
291 | < | Vector3d pos; |
292 | < | |
293 | < | pos[0] = tokenizer.nextTokenAsDouble(); |
294 | < | pos[1] = tokenizer.nextTokenAsDouble(); |
295 | < | pos[2] = tokenizer.nextTokenAsDouble(); |
296 | < | |
297 | < | all_pos_[index] = pos; |
294 | < | }else{ |
295 | < | |
296 | < | bool exist = false; |
297 | < | int indexCount = 0; |
298 | < | |
299 | < | while ( (!exist) && (indexCount < stuntDoubleIndex_.size())){ |
300 | < | if (index == stuntDoubleIndex_[indexCount]) |
301 | < | exist = true; |
302 | < | indexCount++; |
303 | < | } |
304 | < | |
305 | < | StuntDouble* integrableObject; |
306 | < | |
307 | < | if (exist){ |
308 | < | |
309 | < | integrableObject = info_->getIOIndexToIntegrableObject(index); |
310 | < | |
311 | < | int compareindex = integrableObject->getGlobalIntegrableObjectIndex(); |
312 | < | |
313 | < | if (integrableObject == NULL) { |
314 | < | return; |
315 | < | } |
316 | < | |
317 | < | std::string type = tokenizer.nextToken(); |
318 | < | |
319 | < | Vector3d pos; |
320 | < | |
288 | > | std::string type = tokenizer.nextToken(); |
289 | > | int size = type.size(); |
290 | > | |
291 | > | Vector3d pos; |
292 | > | Quat4d q; |
293 | > | |
294 | > | for(int i = 0; i < size; ++i) { |
295 | > | switch(type[i]) { |
296 | > | |
297 | > | case 'p': { |
298 | pos[0] = tokenizer.nextTokenAsDouble(); | |
299 | pos[1] = tokenizer.nextTokenAsDouble(); | |
300 | pos[2] = tokenizer.nextTokenAsDouble(); | |
301 | < | |
301 | > | break; |
302 | > | } |
303 | > | case 'v' : { |
304 | Vector3d vel; | |
305 | vel[0] = tokenizer.nextTokenAsDouble(); | |
306 | vel[1] = tokenizer.nextTokenAsDouble(); | |
307 | vel[2] = tokenizer.nextTokenAsDouble(); | |
308 | < | |
309 | < | |
310 | < | Quat4d q; |
308 | > | break; |
309 | > | } |
310 | > | |
311 | > | case 'q' : { |
312 | if (integrableObject->isDirectional()) { | |
313 | ||
314 | q[0] = tokenizer.nextTokenAsDouble(); | |
315 | q[1] = tokenizer.nextTokenAsDouble(); | |
316 | q[2] = tokenizer.nextTokenAsDouble(); | |
317 | q[3] = tokenizer.nextTokenAsDouble(); | |
318 | < | } |
319 | < | // keep the position in case we need it for a molecular restraint: |
320 | < | |
321 | < | all_pos_[index] = pos; |
322 | < | |
323 | < | // is this io restrained? |
324 | < | GenericData* data = integrableObject->getPropertyByName("Restraint"); |
325 | < | ObjectRestraint* oRest; |
318 | > | |
319 | > | RealType qlen = q.length(); |
320 | > | if (qlen < OpenMD::epsilon) { //check quaternion is not equal to 0 |
321 | > | |
322 | > | sprintf(painCave.errMsg, |
323 | > | "RestReader Error: initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2) ~ 0\n"); |
324 | > | painCave.isFatal = 1; |
325 | > | simError(); |
326 | > | } |
327 | > | |
328 | > | q.normalize(); |
329 | > | } |
330 | > | break; |
331 | > | } |
332 | > | case 'j' : { |
333 | > | Vector3d ji; |
334 | > | if (integrableObject->isDirectional()) { |
335 | > | ji[0] = tokenizer.nextTokenAsDouble(); |
336 | > | ji[1] = tokenizer.nextTokenAsDouble(); |
337 | > | ji[2] = tokenizer.nextTokenAsDouble(); |
338 | > | } |
339 | > | break; |
340 | > | } |
341 | > | case 'f': { |
342 | > | Vector3d force; |
343 | > | force[0] = tokenizer.nextTokenAsDouble(); |
344 | > | force[1] = tokenizer.nextTokenAsDouble(); |
345 | > | force[2] = tokenizer.nextTokenAsDouble(); |
346 | > | break; |
347 | > | } |
348 | > | case 't' : { |
349 | > | Vector3d torque; |
350 | > | torque[0] = tokenizer.nextTokenAsDouble(); |
351 | > | torque[1] = tokenizer.nextTokenAsDouble(); |
352 | > | torque[2] = tokenizer.nextTokenAsDouble(); |
353 | > | break; |
354 | > | } |
355 | > | default: { |
356 | > | sprintf(painCave.errMsg, |
357 | > | "RestReader Error: %s is an unrecognized type\n", type.c_str()); |
358 | > | painCave.isFatal = 1; |
359 | > | simError(); |
360 | > | break; |
361 | > | } |
362 | > | } |
363 | > | // keep the position in case we need it for a molecular restraint: |
364 | > | |
365 | > | all_pos_[index] = pos; |
366 | ||
367 | < | if (data != NULL) { |
368 | < | // make sure we can reinterpret the generic data as restraint data: |
369 | < | RestraintData* restData= dynamic_cast<RestraintData*>(data); |
370 | < | if (restData != NULL) { |
371 | < | // make sure we can reinterpet the restraint data as a pointer to |
367 | > | // is this io restrained? |
368 | > | GenericData* data = integrableObject->getPropertyByName("Restraint"); |
369 | > | ObjectRestraint* oRest; |
370 | > | |
371 | > | if (data != NULL) { |
372 | > | // make sure we can reinterpret the generic data as restraint data: |
373 | > | RestraintData* restData= dynamic_cast<RestraintData*>(data); |
374 | > | if (restData != NULL) { |
375 | > | // make sure we can reinterpet the restraint data as a pointer to |
376 | // an ObjectRestraint: | |
377 | < | oRest = dynamic_cast<ObjectRestraint*>(restData->getData()); |
378 | < | if (oRest != NULL) { |
379 | < | if (integrableObject->isDirectional()) { |
380 | < | oRest->setReferenceStructure(pos, q.toRotationMatrix3()); |
381 | < | } else { |
382 | < | oRest->setReferenceStructure(pos); |
359 | < | } |
377 | > | oRest = dynamic_cast<ObjectRestraint*>(restData->getData()); |
378 | > | if (oRest != NULL) { |
379 | > | if (integrableObject->isDirectional()) { |
380 | > | oRest->setReferenceStructure(pos, q.toRotationMatrix3()); |
381 | > | } else { |
382 | > | oRest->setReferenceStructure(pos); |
383 | } | |
384 | } | |
385 | } |
– | Removed lines |
+ | Added lines |
< | Changed lines |
> | Changed lines |