30 if (Tcl_ListObjGetElements(interp, list, &num, &data) != TCL_OK)
32 if (Tcl_GetDoubleFromObj(interp,data[0],&(result.
x)) != TCL_OK)
34 if (Tcl_GetDoubleFromObj(interp,data[1],&(result.
y)) != TCL_OK)
36 if (Tcl_GetDoubleFromObj(interp,data[2],&(result.
z)) != TCL_OK)
48 doublev[0] = Tcl_NewDoubleObj(v.
x);
49 doublev[1] = Tcl_NewDoubleObj(v.
y);
50 doublev[2] = Tcl_NewDoubleObj(v.
z);
52 return Tcl_NewListObj(3,doublev);
58 int proc_getbond(ClientData, Tcl_Interp *interp,
int argc, Tcl_Obj *
const argv[])
68 Tcl_SetObjResult(interp, Tcl_NewDoubleObj((r2-r1).length()));
75 int proc_getangle(ClientData, Tcl_Interp *interp,
int argc, Tcl_Obj *
const argv[])
77 Vector r1, r2, r3, r12, r32;
89 Tcl_SetObjResult(interp, Tcl_NewDoubleObj(
90 acos((r12*r32)/(r12.
length()*r32.length()))*180/
PI ));
100 Vector r1, r2, r3, r4, r12, r23, r34, A, B, C;
121 Tcl_SetObjResult(interp, Tcl_NewDoubleObj(
122 -atan2((C*B)/(rC*rB),(A*B)/(rA*rB))*180/
PI ));
130 int proc_anglegrad(ClientData, Tcl_Interp *interp,
int argc, Tcl_Obj *
const argv[])
132 Vector r1, r2, r3, r12, r32;
148 BigReal cos_theta = (r12*r32)/(d12*d32);
154 BigReal sin_theta = sqrt(1.0 - cos_theta*cos_theta);
160 Force force1 = c1*(r12*(d12inv*cos_theta) - r32*d32inv);
161 Force force2 = force1;
162 Force force3 = c2*(r32*(d32inv*cos_theta) - r12*d12inv);
163 force2 += force3; force2 *= -1;
171 Tcl_SetObjResult(interp, Tcl_NewListObj(3, forcev));
183 Vector r1, r2, r3, r4, r12, r23, r34;
201 Vector A = cross(r12,r23);
203 Vector B = cross(r23,r34);
209 BigReal cos_phi = (A*B)/(rA*rB);
210 BigReal sin_phi = (C*B)/(rC*rB);
221 if (fabs(sin_phi) > 0.1)
228 Vector dcosdA = rA*(cos_phi*A-B);
229 Vector dcosdB = rB*(cos_phi*B-A);
233 f1.
x = K1*(r23.
y*dcosdA.
z - r23.
z*dcosdA.
y);
234 f1.
y = K1*(r23.
z*dcosdA.
x - r23.
x*dcosdA.
z);
235 f1.
z = K1*(r23.
x*dcosdA.
y - r23.
y*dcosdA.
x);
237 f3.
x = K1*(r23.
z*dcosdB.y - r23.
y*dcosdB.z);
238 f3.
y = K1*(r23.
x*dcosdB.z - r23.
z*dcosdB.x);
239 f3.
z = K1*(r23.
y*dcosdB.x - r23.
x*dcosdB.y);
241 f2.
x = K1*(r12.
z*dcosdA.
y - r12.
y*dcosdA.
z 242 + r34.
y*dcosdB.z - r34.
z*dcosdB.y);
243 f2.
y = K1*(r12.
x*dcosdA.
z - r12.
z*dcosdA.
x 244 + r34.
z*dcosdB.x - r34.
x*dcosdB.z);
245 f2.
z = K1*(r12.
y*dcosdA.
x - r12.
x*dcosdA.
y 246 + r34.
x*dcosdB.y - r34.
y*dcosdB.x);
256 Vector dsindC = rC*(sin_phi*C-B);
257 Vector dsindB = rB*(sin_phi*B-C);
261 f1.
x = K1*((r23.
y*r23.
y + r23.
z*r23.
z)*dsindC.
x 262 - r23.
x*r23.
y*dsindC.
y 263 - r23.
x*r23.
z*dsindC.
z);
264 f1.
y = K1*((r23.
z*r23.
z + r23.
x*r23.
x)*dsindC.
y 265 - r23.
y*r23.
z*dsindC.
z 266 - r23.
y*r23.
x*dsindC.
x);
267 f1.
z = K1*((r23.
x*r23.
x + r23.
y*r23.
y)*dsindC.
z 268 - r23.
z*r23.
x*dsindC.
x 269 - r23.
z*r23.
y*dsindC.
y);
271 f3 = cross(K1,dsindB,r23);
273 f2.
x = K1*(-(r23.
y*r12.
y + r23.
z*r12.
z)*dsindC.
x 274 +(2.0*r23.
x*r12.
y - r12.
x*r23.
y)*dsindC.
y 275 +(2.0*r23.
x*r12.
z - r12.
x*r23.
z)*dsindC.
z 276 +dsindB.z*r34.
y - dsindB.y*r34.
z);
277 f2.
y = K1*(-(r23.
z*r12.
z + r23.
x*r12.
x)*dsindC.
y 278 +(2.0*r23.
y*r12.
z - r12.
y*r23.
z)*dsindC.
z 279 +(2.0*r23.
y*r12.
x - r12.
y*r23.
x)*dsindC.
x 280 +dsindB.x*r34.
z - dsindB.z*r34.
x);
281 f2.
z = K1*(-(r23.
x*r12.
x + r23.
y*r12.
y)*dsindC.
z 282 +(2.0*r23.
z*r12.
x - r12.
z*r23.
x)*dsindC.
x 283 +(2.0*r23.
z*r12.
y - r12.
z*r23.
y)*dsindC.
y 284 +dsindB.y*r34.
x - dsindB.x*r34.
y);
294 Tcl_SetObjResult(interp, Tcl_NewListObj(4, forcev));
305 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL);
307 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL);
309 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL);
311 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL);
313 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL);
static int get_3D_vector(Tcl_Interp *interp, Tcl_Obj *const list, Vector &result)
int Vec_Init(Tcl_Interp *interp)
int proc_dihedralgrad(ClientData, Tcl_Interp *interp, int argc, Tcl_Obj *const argv[])
NAMD_HOST_DEVICE BigReal length(void) const
int proc_getdihedral(ClientData, Tcl_Interp *interp, int argc, Tcl_Obj *const argv[])
static Tcl_Obj * obj_3D_vector(const Vector &v)
int tcl_vector_math_init(Tcl_Interp *interp)
int proc_anglegrad(ClientData, Tcl_Interp *interp, int argc, Tcl_Obj *const argv[])
int proc_getangle(ClientData, Tcl_Interp *interp, int argc, Tcl_Obj *const argv[])
int proc_getbond(ClientData, Tcl_Interp *interp, int argc, Tcl_Obj *const argv[])