#include <stdlib.h>
#include <malloc.h>
#include <errno.h>
#include "TclCommands.h"
#include "Vector.h"
#include "NamdTypes.h"
#include <tcl.h>
#include "Matrix4.C"
#include "TclVec.C"
Go to the source code of this file.
Functions | |
static int | get_3D_vector (Tcl_Interp *interp, Tcl_Obj *const list, Vector &result) |
static Tcl_Obj * | obj_3D_vector (const Vector &v) |
int | proc_getbond (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_getdihedral (ClientData, Tcl_Interp *interp, int argc, Tcl_Obj *const argv[]) |
int | proc_anglegrad (ClientData, Tcl_Interp *interp, int argc, Tcl_Obj *const argv[]) |
int | proc_dihedralgrad (ClientData, Tcl_Interp *interp, int argc, Tcl_Obj *const argv[]) |
int | tcl_vector_math_init (Tcl_Interp *interp) |
static int get_3D_vector | ( | Tcl_Interp * | interp, | |
Tcl_Obj *const | list, | |||
Vector & | result | |||
) | [static] |
Copyright (c) 1995, 1996, 1997, 1998, 1999, 2000 by The Board of Trustees of the University of Illinois. All rights reserved.
Definition at line 25 of file TclCommands.C.
References Vector::x, Vector::y, and Vector::z.
Referenced by proc_anglegrad(), proc_dihedralgrad(), proc_getangle(), proc_getbond(), and proc_getdihedral().
00026 { 00027 int num; 00028 Tcl_Obj **data; 00029 00030 if (Tcl_ListObjGetElements(interp, list, &num, &data) != TCL_OK) 00031 return 0; 00032 if (Tcl_GetDoubleFromObj(interp,data[0],&(result.x)) != TCL_OK) 00033 return 0; 00034 if (Tcl_GetDoubleFromObj(interp,data[1],&(result.y)) != TCL_OK) 00035 return 0; 00036 if (Tcl_GetDoubleFromObj(interp,data[2],&(result.z)) != TCL_OK) 00037 return 0; 00038 00039 return 1; 00040 }
static Tcl_Obj* obj_3D_vector | ( | const Vector & | v | ) | [static] |
Definition at line 44 of file TclCommands.C.
References Vector::x, Vector::y, and Vector::z.
Referenced by proc_anglegrad(), and proc_dihedralgrad().
int proc_anglegrad | ( | ClientData | , | |
Tcl_Interp * | interp, | |||
int | argc, | |||
Tcl_Obj *const | argv[] | |||
) |
Definition at line 130 of file TclCommands.C.
References get_3D_vector(), Vector::length(), and obj_3D_vector().
Referenced by tcl_vector_math_init().
00131 { 00132 Vector r1, r2, r3, r12, r32; 00133 00134 if (argc != 4) 00135 return TCL_ERROR; 00136 if (!get_3D_vector(interp,argv[1],r1)) 00137 return TCL_ERROR; 00138 if (!get_3D_vector(interp,argv[2],r2)) 00139 return TCL_ERROR; 00140 if (!get_3D_vector(interp,argv[3],r3)) 00141 return TCL_ERROR; 00142 00143 r12 = r1 - r2; 00144 BigReal d12 = r12.length(); 00145 r32 = r3 - r2; 00146 BigReal d32 = r32.length(); 00147 00148 BigReal cos_theta = (r12*r32)/(d12*d32); 00149 00150 // Normalize vector r12 and r32 00151 BigReal d12inv = 1. / d12; 00152 BigReal d32inv = 1. / d32; 00153 00154 BigReal sin_theta = sqrt(1.0 - cos_theta*cos_theta); 00155 BigReal diff = 1/sin_theta; 00156 BigReal c1 = diff * d12inv; 00157 BigReal c2 = diff * d32inv; 00158 00159 // Calculate the actual forces 00160 Force force1 = c1*(r12*(d12inv*cos_theta) - r32*d32inv); 00161 Force force2 = force1; 00162 Force force3 = c2*(r32*(d32inv*cos_theta) - r12*d12inv); 00163 force2 += force3; force2 *= -1; 00164 00165 Tcl_Obj* forcev[3]; 00166 00167 forcev[0] = obj_3D_vector(force1); 00168 forcev[1] = obj_3D_vector(force2); 00169 forcev[2] = obj_3D_vector(force3); 00170 00171 Tcl_SetObjResult(interp, Tcl_NewListObj(3, forcev)); 00172 00173 return TCL_OK; 00174 }
int proc_dihedralgrad | ( | ClientData | , | |
Tcl_Interp * | interp, | |||
int | argc, | |||
Tcl_Obj *const | argv[] | |||
) |
Definition at line 180 of file TclCommands.C.
References A, B, cross(), get_3D_vector(), Vector::length(), obj_3D_vector(), Vector::x, Vector::y, and Vector::z.
Referenced by tcl_vector_math_init().
00181 { 00182 BigReal K1; 00183 Vector r1, r2, r3, r4, r12, r23, r34; 00184 00185 if (argc != 5) 00186 return TCL_ERROR; 00187 if (!get_3D_vector(interp,argv[1],r1)) 00188 return TCL_ERROR; 00189 if (!get_3D_vector(interp,argv[2],r2)) 00190 return TCL_ERROR; 00191 if (!get_3D_vector(interp,argv[3],r3)) 00192 return TCL_ERROR; 00193 if (!get_3D_vector(interp,argv[4],r4)) 00194 return TCL_ERROR; 00195 00196 r12 = r1 - r2; 00197 r23 = r2 - r3; 00198 r34 = r3 - r4; 00199 00200 // Calculate the cross products and distances 00201 Vector A = cross(r12,r23); 00202 BigReal rA = A.length(); 00203 Vector B = cross(r23,r34); 00204 BigReal rB = B.length(); 00205 Vector C = cross(r23,A); 00206 BigReal rC = C.length(); 00207 00208 // Calculate the sin and cos 00209 BigReal cos_phi = (A*B)/(rA*rB); 00210 BigReal sin_phi = (C*B)/(rC*rB); 00211 00212 Force f1,f2,f3; 00213 00214 // Normalize B 00215 rB = 1.0/rB; 00216 B *= rB; 00217 00218 // We first need to figure out whether the 00219 // sin or cos form will be more stable. For this, 00220 // just look at the value of phi 00221 if (fabs(sin_phi) > 0.1) 00222 { 00223 // use the sin version to avoid 1/cos terms 00224 00225 // Normalize A 00226 rA = 1.0/rA; 00227 A *= rA; 00228 Vector dcosdA = rA*(cos_phi*A-B); 00229 Vector dcosdB = rB*(cos_phi*B-A); 00230 00231 K1 = -1/sin_phi; 00232 00233 f1.x = K1*(r23.y*dcosdA.z - r23.z*dcosdA.y); 00234 f1.y = K1*(r23.z*dcosdA.x - r23.x*dcosdA.z); 00235 f1.z = K1*(r23.x*dcosdA.y - r23.y*dcosdA.x); 00236 00237 f3.x = K1*(r23.z*dcosdB.y - r23.y*dcosdB.z); 00238 f3.y = K1*(r23.x*dcosdB.z - r23.z*dcosdB.x); 00239 f3.z = K1*(r23.y*dcosdB.x - r23.x*dcosdB.y); 00240 00241 f2.x = K1*(r12.z*dcosdA.y - r12.y*dcosdA.z 00242 + r34.y*dcosdB.z - r34.z*dcosdB.y); 00243 f2.y = K1*(r12.x*dcosdA.z - r12.z*dcosdA.x 00244 + r34.z*dcosdB.x - r34.x*dcosdB.z); 00245 f2.z = K1*(r12.y*dcosdA.x - r12.x*dcosdA.y 00246 + r34.x*dcosdB.y - r34.y*dcosdB.x); 00247 } 00248 else 00249 { 00250 // This angle is closer to 0 or 180 than it is to 00251 // 90, so use the cos version to avoid 1/sin terms 00252 00253 // Normalize C 00254 rC = 1.0/rC; 00255 C *= rC; 00256 Vector dsindC = rC*(sin_phi*C-B); 00257 Vector dsindB = rB*(sin_phi*B-C); 00258 00259 K1 = 1/cos_phi; 00260 00261 f1.x = K1*((r23.y*r23.y + r23.z*r23.z)*dsindC.x 00262 - r23.x*r23.y*dsindC.y 00263 - r23.x*r23.z*dsindC.z); 00264 f1.y = K1*((r23.z*r23.z + r23.x*r23.x)*dsindC.y 00265 - r23.y*r23.z*dsindC.z 00266 - r23.y*r23.x*dsindC.x); 00267 f1.z = K1*((r23.x*r23.x + r23.y*r23.y)*dsindC.z 00268 - r23.z*r23.x*dsindC.x 00269 - r23.z*r23.y*dsindC.y); 00270 00271 f3 = cross(K1,dsindB,r23); 00272 00273 f2.x = K1*(-(r23.y*r12.y + r23.z*r12.z)*dsindC.x 00274 +(2.0*r23.x*r12.y - r12.x*r23.y)*dsindC.y 00275 +(2.0*r23.x*r12.z - r12.x*r23.z)*dsindC.z 00276 +dsindB.z*r34.y - dsindB.y*r34.z); 00277 f2.y = K1*(-(r23.z*r12.z + r23.x*r12.x)*dsindC.y 00278 +(2.0*r23.y*r12.z - r12.y*r23.z)*dsindC.z 00279 +(2.0*r23.y*r12.x - r12.y*r23.x)*dsindC.x 00280 +dsindB.x*r34.z - dsindB.z*r34.x); 00281 f2.z = K1*(-(r23.x*r12.x + r23.y*r12.y)*dsindC.z 00282 +(2.0*r23.z*r12.x - r12.z*r23.x)*dsindC.x 00283 +(2.0*r23.z*r12.y - r12.z*r23.y)*dsindC.y 00284 +dsindB.y*r34.x - dsindB.x*r34.y); 00285 } 00286 00287 Tcl_Obj* forcev[4]; 00288 00289 forcev[0] = obj_3D_vector(f1); 00290 forcev[1] = obj_3D_vector(f2-f1); 00291 forcev[2] = obj_3D_vector(f3-f2); 00292 forcev[3] = obj_3D_vector(-f3); 00293 00294 Tcl_SetObjResult(interp, Tcl_NewListObj(4, forcev)); 00295 00296 return TCL_OK; 00297 }
int proc_getangle | ( | ClientData | , | |
Tcl_Interp * | interp, | |||
int | argc, | |||
Tcl_Obj *const | argv[] | |||
) |
Definition at line 75 of file TclCommands.C.
References get_3D_vector(), Vector::length(), and PI.
Referenced by tcl_vector_math_init().
00076 { 00077 Vector r1, r2, r3, r12, r32; 00078 00079 if (argc != 4) 00080 return TCL_ERROR; 00081 if (!get_3D_vector(interp,argv[1],r1)) 00082 return TCL_ERROR; 00083 if (!get_3D_vector(interp,argv[2],r2)) 00084 return TCL_ERROR; 00085 if (!get_3D_vector(interp,argv[3],r3)) 00086 return TCL_ERROR; 00087 r12 = r1 - r2; 00088 r32 = r3 - r2; 00089 Tcl_SetObjResult(interp, Tcl_NewDoubleObj( 00090 acos((r12*r32)/(r12.length()*r32.length()))*180/PI )); 00091 return TCL_OK; 00092 }
int proc_getbond | ( | ClientData | , | |
Tcl_Interp * | interp, | |||
int | argc, | |||
Tcl_Obj *const | argv[] | |||
) |
Definition at line 58 of file TclCommands.C.
References get_3D_vector().
Referenced by tcl_vector_math_init().
00059 { 00060 Vector r1, r2; 00061 00062 if (argc != 3) 00063 return TCL_ERROR; 00064 if (!get_3D_vector(interp,argv[1],r1)) 00065 return TCL_ERROR; 00066 if (!get_3D_vector(interp,argv[2],r2)) 00067 return TCL_ERROR; 00068 Tcl_SetObjResult(interp, Tcl_NewDoubleObj((r2-r1).length())); 00069 return TCL_OK; 00070 }
int proc_getdihedral | ( | ClientData | , | |
Tcl_Interp * | interp, | |||
int | argc, | |||
Tcl_Obj *const | argv[] | |||
) |
Definition at line 97 of file TclCommands.C.
References A, B, cross(), get_3D_vector(), Vector::length(), and PI.
Referenced by tcl_vector_math_init().
00098 { 00099 BigReal rA, rB, rC; 00100 Vector r1, r2, r3, r4, r12, r23, r34, A, B, C; 00101 00102 if (argc != 5) 00103 return TCL_ERROR; 00104 if (!get_3D_vector(interp,argv[1],r1)) 00105 return TCL_ERROR; 00106 if (!get_3D_vector(interp,argv[2],r2)) 00107 return TCL_ERROR; 00108 if (!get_3D_vector(interp,argv[3],r3)) 00109 return TCL_ERROR; 00110 if (!get_3D_vector(interp,argv[4],r4)) 00111 return TCL_ERROR; 00112 r12 = r1 - r2; 00113 r23 = r2 - r3; 00114 r34 = r3 - r4; 00115 A = cross(r12,r23); 00116 B = cross(r23,r34); 00117 C = cross(r23,A); 00118 rA = A.length(); 00119 rB = B.length(); 00120 rC = C.length(); 00121 Tcl_SetObjResult(interp, Tcl_NewDoubleObj( 00122 -atan2((C*B)/(rC*rB),(A*B)/(rA*rB))*180/PI )); 00123 return TCL_OK; 00124 }
int tcl_vector_math_init | ( | Tcl_Interp * | ) |
Copyright (c) 1995, 1996, 1997, 1998, 1999, 2000 by The Board of Trustees of the University of Illinois. All rights reserved.
Definition at line 299 of file TclCommands.C.
References proc_anglegrad(), proc_dihedralgrad(), proc_getangle(), proc_getbond(), proc_getdihedral(), and Vec_Init().
Referenced by ComputeTclBC::ComputeTclBC(), ScriptTcl::ScriptTcl(), and ScriptTcl::tclsh().
00299 { 00300 00301 // first import from TclVec.C stolen from VMD 00302 Vec_Init(interp); 00303 00304 Tcl_CreateObjCommand(interp, "getbond", proc_getbond, 00305 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL); 00306 Tcl_CreateObjCommand(interp, "getangle", proc_getangle, 00307 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL); 00308 Tcl_CreateObjCommand(interp, "getdihedral", proc_getdihedral, 00309 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL); 00310 Tcl_CreateObjCommand(interp, "anglegrad", proc_anglegrad, 00311 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL); 00312 Tcl_CreateObjCommand(interp, "dihedralgrad", proc_dihedralgrad, 00313 (ClientData) NULL, (Tcl_CmdDeleteProc *) NULL); 00314 00315 return TCL_OK; 00316 }