Commit b77ab693 authored by Mark Friedrichs's avatar Mark Friedrichs
Browse files

OutOfPlaneBend force fix

parent 25ce1664
...@@ -1163,16 +1163,26 @@ void kCalculateAmoebaLocalForces_kernel() ...@@ -1163,16 +1163,26 @@ void kCalculateAmoebaLocalForces_kernel()
float bkk3 = (rdb2 != 0.0f ) ? bkk2/rdb2 : 0.0f; float bkk3 = (rdb2 != 0.0f ) ? bkk2/rdb2 : 0.0f;
bkk2 = rdb2 - bkk2; bkk2 = rdb2 - bkk2;
float adXcd_0 = yad*zcd - zad*ycd;
float adXcd_1 = zad*xcd - xad*zcd;
float adXcd_2 = xad*ycd - yad*xcd;
float adXcd_nrm2 = adXcd_0*adXcd_0 + adXcd_1*adXcd_1 + adXcd_2*adXcd_2;
float adXcd_dot_db = xdb*adXcd_0 + ydb*adXcd_1 + zdb*adXcd_2;
adXcd_dot_db /= sqrt(rdb2*adXcd_nrm2);
float angle; float angle;
/*
if( bkk3 < 1.0e-6f ){ if( bkk3 < 1.0e-6f ){
angle = sqrtf( bkk3 ); angle = sqrtf( bkk3 );
} else if( bkk3 < 0.98f ){ } else if( bkk3 < 0.98f ){
float cosine = sqrtf(1.0f - bkk3); float cosine = sqrtf(1.0f - bkk3);
angle = acos(cosine); angle = acos(cosine);
} else { } else {
float sin = sqrtf(bkk3); angle = abs( asin(adXcd_dot_db) );
angle = asin(sin);
} }
*/
angle = abs( asin(adXcd_dot_db) );
// find the out-of-plane energy and master chain rule terms // find the out-of-plane energy and master chain rule terms
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment