#include <ComputeSphericalBC.h>
Public Member Functions | |
ComputeSphericalBC (ComputeID c, PatchID pid) | |
virtual | ~ComputeSphericalBC () |
virtual void | doForce (FullAtom *p, Results *r) |
Public Attributes | |
SubmitReduction * | reduction |
Copyright (c) 1995, 1996, 1997, 1998, 1999, 2000 by The Board of Trustees of the University of Illinois. All rights reserved.
Definition at line 13 of file ComputeSphericalBC.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 25 of file ComputeSphericalBC.C.
References FALSE, Node::Object(), ReductionMgr::Object(), reduction, REDUCTIONS_BASIC, Node::simParameters, simParams, SimParameters::sphericalBCexp1, SimParameters::sphericalBCexp2, SimParameters::sphericalBCk1, SimParameters::sphericalBCk2, SimParameters::sphericalBCr1, SimParameters::sphericalBCr2, SimParameters::sphericalCenter, TRUE, and ReductionMgr::willSubmit().
00026 : ComputeHomePatch(c,pid) 00027 { 00028 reduction = ReductionMgr::Object()->willSubmit(REDUCTIONS_BASIC); 00029 00030 SimParameters *simParams = Node::Object()->simParameters; 00031 00032 // Get parameters from the SimParameters object 00033 r1 = simParams->sphericalBCr1; 00034 r2 = simParams->sphericalBCr2; 00035 r1_2 = r1*r1; 00036 r2_2 = r2*r2; 00037 k1 = simParams->sphericalBCk1; 00038 k2 = simParams->sphericalBCk2; 00039 exp1 = simParams->sphericalBCexp1; 00040 exp2 = simParams->sphericalBCexp2; 00041 00042 // Check to see if this is one set of parameters or two 00043 if (r2 > -1.0) 00044 { 00045 twoForces = TRUE; 00046 } 00047 else 00048 { 00049 twoForces = FALSE; 00050 } 00051 00052 center = simParams->sphericalCenter; 00053 00054 }
ComputeSphericalBC::~ComputeSphericalBC | ( | ) | [virtual] |
Definition at line 66 of file ComputeSphericalBC.C.
References reduction.
00068 { 00069 delete reduction; 00070 }
Implements ComputeHomePatch.
Definition at line 86 of file ComputeSphericalBC.C.
References Results::f, forces, SubmitReduction::item(), j, Results::normal, ComputeHomePatch::numAtoms, CompAtom::position, reduction, REDUCTION_BC_ENERGY, SubmitReduction::submit(), Vector::x, Vector::y, and Vector::z.
00087 { 00088 Vector diff; // Distance from atom to center of sphere 00089 Vector f; // Calculated force vector 00090 int i, j; // Loop counters 00091 BigReal dist, dist_2; // Distance from atom to center, and this 00092 // distance squared 00093 BigReal rval; // Difference between distance from atom 00094 // to center and radius of sphere 00095 BigReal eval; // Energy value for this atom 00096 BigReal fval; // Force magnitude for this atom 00097 00098 // aliases to work with old code 00099 FullAtom *x = p; 00100 Force *forces = r->f[Results::normal]; 00101 BigReal energy = 0; 00102 00103 // Loop through and check each atom 00104 for (i=0; i<numAtoms; i++) 00105 { 00106 // Calculate the vector from the atom to the center of the 00107 // sphere 00108 diff.x = x[i].position.x - center.x; 00109 diff.y = x[i].position.y - center.y; 00110 diff.z = x[i].position.z - center.z; 00111 00112 // Calculate the distance squared 00113 dist_2 = diff.x*diff.x + diff.y*diff.y + diff.z*diff.z; 00114 00115 // Look to see if we are outside either radius 00116 if ( (dist_2 > r1_2) || (twoForces && (dist_2 > r2_2)) ) 00117 { 00118 // Calculate the distance to the center 00119 dist = sqrt(dist_2); 00120 00121 // Normalize the direction vector 00122 diff /= dist; 00123 00124 // Check to see if we are outside radius 1 00125 if (dist > r1) 00126 { 00127 // Assign the force vector to the 00128 // unit direction vector 00129 f.x = diff.x; 00130 f.y = diff.y; 00131 f.z = diff.z; 00132 00133 // Calculate the energy which is 00134 // e = k1*(r_i-r_center)^exp1 00135 eval = k1; 00136 rval = fabs(dist - r1); 00137 00138 for (j=0; j<exp1; j++) 00139 { 00140 eval *= rval; 00141 } 00142 00143 energy += eval; 00144 00145 // Now calculate the force which is 00146 // e = -k1*exp1*(r_i-r_center1)^(exp1-1) 00147 fval = -exp1*k1; 00148 00149 for (j=0; j<exp1-1; j++) 00150 { 00151 fval *= rval; 00152 } 00153 00154 // Multiply the force magnitude to the 00155 // unit direction vector to get the 00156 // resulting force 00157 f *= fval; 00158 00159 // Add the force to the force vectors 00160 forces[i].x += f.x; 00161 forces[i].y += f.y; 00162 forces[i].z += f.z; 00163 } 00164 00165 // Check to see if there is a second radius 00166 // and if we are outside of it 00167 if (twoForces && (dist > r2) ) 00168 { 00169 // Assign the force vector to the 00170 // unit direction vector 00171 f.x = diff.x; 00172 f.y = diff.y; 00173 f.z = diff.z; 00174 00175 // Calculate the energy which is 00176 // e = k2*(r_i-r_center2)^exp2 00177 eval = k2; 00178 rval = fabs(dist - r2); 00179 00180 for (j=0; j<exp2; j++) 00181 { 00182 eval *= rval; 00183 } 00184 00185 energy += eval; 00186 00187 // Now calculate the force which is 00188 // e = -k2*exp2*(r_i-r_center2)^(exp2-1) 00189 fval = -exp2*k2; 00190 00191 for (j=0; j<exp2-1; j++) 00192 { 00193 fval *= rval; 00194 } 00195 00196 // Multiply the force magnitude to the 00197 // unit direction vector to get the 00198 // resulting force 00199 f *= fval; 00200 00201 // Add the force to the force vectors 00202 forces[i].x += f.x; 00203 forces[i].y += f.y; 00204 forces[i].z += f.z; 00205 } 00206 } 00207 } 00208 00209 reduction->item(REDUCTION_BC_ENERGY) += energy; 00210 reduction->submit(); 00211 00212 }
Definition at line 32 of file ComputeSphericalBC.h.
Referenced by ComputeSphericalBC(), doForce(), and ~ComputeSphericalBC().