#include <ComputeRestraints.h>
Inheritance diagram for ComputeRestraints:
Public Member Functions | |
| ComputeRestraints (ComputeID c, PatchID pid) | |
| virtual | ~ComputeRestraints () |
| virtual void | doForce (CompAtom *p, CompAtomExt *pExt, Results *r) |
Public Attributes | |
| SubmitReduction * | reduction |
Definition at line 13 of file ComputeRestraints.h.
|
||||||||||||
|
Copyright (c) 1995, 1996, 1997, 1998, 1999, 2000 by The Board of Trustees of the University of Illinois. All rights reserved. Definition at line 20 of file ComputeRestraints.C. References SimParameters::constraintExp, SimParameters::constrXOn, SimParameters::constrYOn, SimParameters::constrZOn, SimParameters::movingConstraintsOn, SimParameters::movingConsVel, Node::Object(), ReductionMgr::Object(), reduction, REDUCTIONS_BASIC, SimParameters::rotConsAxis, SimParameters::rotConsPivot, SimParameters::rotConstraintsOn, SimParameters::rotConsVel, SimParameters::selectConstraintsOn, Node::simParameters, simParams, and ReductionMgr::willSubmit(). 00021 : ComputePatch(c,pid) 00022 { 00023 reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_BASIC); 00024 00025 SimParameters *simParams = Node::Object()->simParameters; 00026 00027 // Get parameters from the SimParameters object 00028 consExp = simParams->constraintExp; 00029 00030 //****** BEGIN selective restraints (X,Y,Z) changes 00031 consSelectOn = simParams->selectConstraintsOn; 00032 if (consSelectOn) { 00033 consSelectX = simParams->constrXOn; 00034 consSelectY = simParams->constrYOn; 00035 consSelectZ = simParams->constrZOn; 00036 } 00037 //****** END selective restraints (X,Y,Z) changes 00038 00039 //****** BEGIN moving constraints changes 00040 consMoveOn = simParams->movingConstraintsOn; 00041 if (consMoveOn) { 00042 moveVel = simParams->movingConsVel; 00043 } 00044 //****** END moving constraints changes 00045 //****** BEGIN rotating constraints changes 00046 consRotOn = simParams->rotConstraintsOn; 00047 if (consRotOn) { 00048 rotVel = simParams->rotConsVel; 00049 rotAxis = simParams->rotConsAxis; 00050 rotPivot = simParams->rotConsPivot; 00051 } 00052 //****** END rotating constraints changes 00053 00054 }
|
|
|
Definition at line 65 of file ComputeRestraints.C. 00067 {
00068 delete reduction;
00069 }
|
|
||||||||||||||||
|
Implements ComputePatch. Definition at line 77 of file ComputeRestraints.C. References ADD_TENSOR_OBJECT, ADD_VECTOR_OBJECT, BigReal, SimParameters::constraintScaling, Lattice::delta(), Results::f, Patch::flags, Force, Molecule::get_cons_params(), Molecule::is_atom_constrained(), SubmitReduction::item(), Patch::lattice, Vector::length2(), mat_multiply_vec(), Node::molecule, Node::Object(), outer(), Real, reduction, REDUCTION_BC_ENERGY, REDUCTION_EXT_FORCE_NORMAL, REDUCTION_VIRIAL_NORMAL, Node::simParameters, Flags::step, SubmitReduction::submit(), vec_rotation_matrix(), Vector::x, Vector::y, and Vector::z. 00078 {
00079 Molecule *molecule = Node::Object()->molecule;
00080 Real k; // Force constant
00081 BigReal scaling = Node::Object()->simParameters->constraintScaling;
00082 Vector refPos; // Reference position
00083 BigReal r, r2; // r=distance between atom position and the
00084 // reference position, r2 = r^2
00085 Vector Rij; // vector between current position and reference pos
00086 BigReal value; // Current calculated value
00087
00088 // aliases to work with old code
00089 Force *f = res->f[Results::normal];
00090 BigReal energy = 0;
00091 BigReal m[9];
00092 Tensor virial;
00093 Vector netForce = 0;
00094
00095 // BEGIN moving and rotating constraint changes ******
00096
00097 // This version only allows one atom to be moved
00098 // and only ALL ref positions to be rotated
00099
00100 int currentTime = patch->flags.step;
00101 if (consRotOn) {
00102 vec_rotation_matrix(rotVel * currentTime, rotAxis, m);
00103 }
00104
00105 // END moving and rotating constraint changes ******
00106
00107
00108 if (scaling != 0.) for (int localID=0; localID<numAtoms; ++localID)
00109 {
00110 if (molecule->is_atom_constrained(pExt[localID].id))
00111 {
00112 molecule->get_cons_params(k, refPos, pExt[localID].id);
00113
00114 k *= scaling;
00115
00116 // BEGIN moving and rotating constraint changes ******
00117
00118 if (consMoveOn) {
00119 refPos += currentTime * moveVel;
00120 }
00121 else if(consRotOn) {
00122 refPos = mat_multiply_vec(refPos - rotPivot, m) + rotPivot;
00123 }
00124
00125 // END moving and rotating constraint changes *******
00126
00127 Rij = patch->lattice.delta(refPos,p[localID].position);
00128 Vector vpos = refPos - Rij;
00129
00130 //****** BEGIN selective restraints (X,Y,Z) changes
00131 if (consSelectOn) { // check which components we want to restrain:
00132 if (!consSelectX) {Rij.x=0.0;} // by setting the appropriate
00133 if (!consSelectY) {Rij.y=0.0;} // Cartesian component to zero
00134 if (!consSelectZ) {Rij.z=0.0;} // we will avoid a restoring force
00135 } // along that component, and the
00136 // energy will also be correct
00137 //****** END selective restraints (X,Y,Z) changes
00138
00139
00140 // Calculate the distance and the distance squared
00141 r2 = Rij.length2();
00142 r = sqrt(r2);
00143
00144 // Only calculate the energy and the force if the distance is
00145 // non-zero. Otherwise, you could end up dividing by 0, which
00146 // is bad
00147 if (r>0.0)
00148 {
00149 value=k;
00150
00151 // Loop through and multiple k by r consExp times.
00152 // i.e. calculate kr^e
00153 // I know, I could use pow(), but I don't trust it.
00154 for (int k=0; k<consExp; ++k)
00155 {
00156 value *= r;
00157 }
00158
00159 // Add to the energy total
00160 energy += value;
00161
00162 // Now calculate the force, which is ekr^(e-1). Also, divide
00163 // by another factor of r to normalize the vector before we
00164 // multiple it by the magnitude
00165 value *= consExp;
00166 value /= r2;
00167
00168 Rij *= value;
00169 //iout << iINFO << "restraining force" << Rij << "\n";
00170
00171 f[localID] += Rij;
00172 netForce += Rij;
00173 virial += outer(Rij,vpos);
00174 }
00175 }
00176 }
00177
00178 reduction->item(REDUCTION_BC_ENERGY) += energy;
00179 ADD_TENSOR_OBJECT(reduction,REDUCTION_VIRIAL_NORMAL,virial);
00180 ADD_VECTOR_OBJECT(reduction,REDUCTION_EXT_FORCE_NORMAL,netForce);
00181 reduction->submit();
00182
00183 }
|
|
|
Definition at line 39 of file ComputeRestraints.h. Referenced by ComputeRestraints(), and doForce(). |
1.3.9.1