ComputeCylindricalBC Class Reference

#include <ComputeCylindricalBC.h>

Inheritance diagram for ComputeCylindricalBC:

ComputeHomePatch Compute List of all members.

Public Member Functions

 ComputeCylindricalBC (ComputeID c, PatchID pid)
virtual ~ComputeCylindricalBC ()
virtual void doForce (FullAtom *p, Results *r)

Public Attributes

SubmitReductionreduction

Detailed Description

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 ComputeCylindricalBC.h.


Constructor & Destructor Documentation

ComputeCylindricalBC::ComputeCylindricalBC ( ComputeID  c,
PatchID  pid 
)

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 ComputeCylindricalBC.C.

References FALSE, Node::Object(), ReductionMgr::Object(), reduction, REDUCTIONS_BASIC, Node::simParameters, simParams, 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         axis = simParams->cylindricalBCAxis;
00034         r1 = simParams->cylindricalBCr1;
00035         r2 = simParams->cylindricalBCr2;
00036         r1_2 = r1*r1;
00037         r2_2 = r2*r2;
00038         k1 = simParams->cylindricalBCk1;
00039         k2 = simParams->cylindricalBCk2;
00040         exp1 = simParams->cylindricalBCexp1;
00041         exp2 = simParams->cylindricalBCexp2;
00042         //Additions for ends
00043         l1 = simParams->cylindricalBCl1;
00044         l2 = simParams->cylindricalBCl2;
00045         l1_2 = l1*l1;
00046         l2_2 = l2*l2;
00047 
00048         //  Check to see if this is one set of parameters or two
00049         if (r2 > -1.0)
00050         {
00051                 twoForces = TRUE;
00052         }
00053         else
00054         {
00055                 twoForces = FALSE;
00056         }
00057 
00058         center = simParams->cylindricalCenter;
00059 
00060 }

ComputeCylindricalBC::~ComputeCylindricalBC (  )  [virtual]

Definition at line 72 of file ComputeCylindricalBC.C.

References reduction.

00074 {
00075         delete reduction;
00076 }


Member Function Documentation

void ComputeCylindricalBC::doForce ( FullAtom p,
Results r 
) [virtual]

Implements ComputeHomePatch.

Definition at line 93 of file ComputeCylindricalBC.C.

References Results::f, f, forces, if(), SubmitReduction::item(), j, Results::normal, ComputeHomePatch::numAtoms, reduction, REDUCTION_BC_ENERGY, SubmitReduction::submit(), Vector::x, x, Vector::y, and Vector::z.

00094 {
00095         Vector diff;            //  Distance from atom to center of cylinder
00096         Vector f;               //  Calculated force vector
00097         int i, j;               //  Loop counters
00098         BigReal dist, dist_2;   //  Distance from atom to center, and this
00099                                 //  distance squared
00100         BigReal rval;           //  Difference between distance from atom
00101                                 //  to center and radius of cylinder
00102         BigReal eval;           //  Energy value for this atom
00103         BigReal fval;           //  Force magnitude for this atom
00104 
00105         // aliases to work with old code
00106         FullAtom *x = p;
00107         Force *forces = r->f[Results::normal];
00108         BigReal energy = 0;
00109 
00110         //  Loop through and check each atom
00111         for (i=0; i<numAtoms; i++)
00112         {
00113                 //  Calculate the vector from the atom to the center of the
00114                 //  cylinder
00115                 diff.x = ( axis == 'x' ? 0.0 : x[i].position.x - center.x );
00116                 diff.y = ( axis == 'y' ? 0.0 : x[i].position.y - center.y );
00117                 diff.z = ( axis == 'z' ? 0.0 : x[i].position.z - center.z );
00118                 
00119                 //  Calculate the distance squared
00120                 dist_2 = diff.x*diff.x + diff.y*diff.y + diff.z*diff.z;
00121 
00122                 //  Look to see if we are outside either radius
00123                 if ( (dist_2 > r1_2) || (twoForces && (dist_2 > r2_2)) )
00124                 {
00125                         //  Calculate the distance to the center
00126                         dist = sqrt(dist_2);
00127 
00128                         //  Normalize the direction vector
00129                         diff /= dist;
00130 
00131                         //  Check to see if we are outside radius 1
00132                         if (dist > r1)
00133                         {
00134                                 //  Assign the force vector to the
00135                                 //  unit direction vector
00136                                 f.x = diff.x;
00137                                 f.y = diff.y;
00138                                 f.z = diff.z;
00139 
00140                                 //  Calculate the energy which is
00141                                 //  e = k1*(r_i-r_center)^exp1
00142                                 eval = k1;
00143                                 rval = fabs(dist - r1);
00144 
00145                                 for (j=0; j<exp1; j++)
00146                                 {
00147                                         eval *= rval;
00148                                 }
00149 
00150                                 energy += eval;
00151 
00152                                 //  Now calculate the force which is
00153                                 //  e = -k1*exp1*(r_i-r_center1)^(exp1-1)
00154                                 fval = -exp1*k1;
00155 
00156                                 for (j=0; j<exp1-1; j++)
00157                                 {
00158                                         fval *= rval;
00159                                 }
00160 
00161                                 //  Multiply the force magnitude to the
00162                                 //  unit direction vector to get the
00163                                 //  resulting force
00164                                 f *= fval;
00165 
00166                                 //  Add the force to the force vectors
00167                                 forces[i].x += f.x;
00168                                 forces[i].y += f.y;
00169                                 forces[i].z += f.z;
00170                         }
00171 
00172                         //  Check to see if there is a second radius
00173                         //  and if we are outside of it
00174                         if (twoForces && (dist > r2) )
00175                         {
00176                                 //  Assign the force vector to the
00177                                 //  unit direction vector
00178                                 f.x = diff.x;
00179                                 f.y = diff.y;
00180                                 f.z = diff.z;
00181 
00182                                 //  Calculate the energy which is
00183                                 //  e = k2*(r_i-r_center2)^exp2
00184                                 eval = k2;
00185                                 rval = fabs(dist - r2);
00186 
00187                                 for (j=0; j<exp2; j++)
00188                                 {
00189                                         eval *= rval;
00190                                 }
00191 
00192                                 energy += eval;
00193 
00194                                 //  Now calculate the force which is
00195                                 //  e = -k2*exp2*(r_i-r_center2)^(exp2-1)
00196                                 fval = -exp2*k2;
00197 
00198                                 for (j=0; j<exp2-1; j++)
00199                                 {
00200                                         fval *= rval;
00201                                 }
00202 
00203                                 //  Multiply the force magnitude to the
00204                                 //  unit direction vector to get the
00205                                 //  resulting force
00206                                 f *= fval;
00207 
00208                                 //  Add the force to the force vectors
00209                                 forces[i].x += f.x;
00210                                 forces[i].y += f.y;
00211                                 forces[i].z += f.z;
00212                         }
00213                 }
00214         }
00215        //  Loop through and check each atom
00216         for (i=0; i<numAtoms; i++)
00217         {
00218                 //  Calculate the vector from the atom to the center of the
00219                 //  cylinder
00220                 diff.x = ( axis != 'x' ? 0.0 : x[i].position.x - center.x );
00221                 diff.y = ( axis != 'y' ? 0.0 : x[i].position.y - center.y );
00222                 diff.z = ( axis != 'z' ? 0.0 : x[i].position.z - center.z );
00223 
00224                 //  Calculate the distance squared
00225                 dist_2 = diff.x*diff.x + diff.y*diff.y + diff.z*diff.z;
00226 
00227                 //  Look to see if we are outside either radius
00228                 if ( (dist_2 > l1_2)  || (twoForces && (dist_2 > l2_2)) )
00229                 {
00230                         //  Calculate the distance to the center
00231                         dist = sqrt(dist_2);
00232 
00233                         //  Normalize the direction vector
00234                         diff /= dist;
00235 
00236                         //  Check to see if we are outside radius 1
00237                         if (dist > l1)
00238                         {
00239 //printf ("Do faces one force z=%f\n", dist);
00240                                 //  Assign the force vector to the
00241                                 //  unit direction vector
00242                                 f.x = diff.x;
00243                                 f.y = diff.y;
00244                                 f.z = diff.z;
00245 
00246                                 //  Calculate the energy which is
00247                                 //  e = k1*(r_i-r_center)^exp1
00248                                 eval = k1;
00249                                 rval = fabs(dist - l1);
00250 
00251                                 for (j=0; j<exp1; j++)
00252                                 {
00253                                         eval *= rval;
00254                                 }
00255 
00256                                 energy += eval;
00257 
00258                                 //  Now calculate the force which is
00259                                 //  e = -k1*exp1*(r_i-r_center1)^(exp1-1)
00260                                 fval = -exp1*k1;
00261 
00262                                 for (j=0; j<exp1-1; j++)
00263                                 {
00264                                         fval *= rval;
00265                                 }
00266 
00267                                 //  Multiply the force magnitude to the
00268                                 //  unit direction vector to get the
00269                                 //  resulting force
00270                                 f *= fval;
00271 
00272                                 //  Add the force to the force vectors
00273                                 forces[i].x += f.x;
00274                                 forces[i].y += f.y;
00275                                 forces[i].z += f.z;
00276                         }
00277 
00278                         //  Check to see if there is a second radius
00279                         //  and if we are outside of it
00280                         if (twoForces && (dist > l2) )
00281                         {
00282 //printf ("Do faces two force z=%f\n", dist);
00283                                 //  Assign the force vector to the
00284                                 //  unit direction vector
00285                                 f.x = diff.x;
00286                                 f.y = diff.y;
00287                                 f.z = diff.z;
00288 
00289                                 //  Calculate the energy which is
00290                                 //  e = k2*(r_i-r_center2)^exp2
00291                                 eval = k2;
00292                                 rval = fabs(dist - l2);
00293 
00294                                 for (j=0; j<exp2; j++)
00295                                 {
00296                                         eval *= rval;
00297                                 }
00298 
00299                                 energy += eval;
00300 
00301                                 //  Now calculate the force which is
00302                                 //  e = -k2*exp2*(r_i-r_center2)^(exp2-1)
00303                                 fval = -exp2*k2;
00304 
00305                                 for (j=0; j<exp2-1; j++)
00306                                 {
00307                                         fval *= rval;
00308                                 }
00309 
00310                                 //  Multiply the force magnitude to the
00311                                 //  unit direction vector to get the
00312                                 //  resulting force
00313                                 f *= fval;
00314 
00315                                 //  Add the force to the force vectors
00316                                 forces[i].x += f.x;
00317                                 forces[i].y += f.y;
00318                                 forces[i].z += f.z;
00319                         }
00320                 }
00321         }
00322 
00323     reduction->item(REDUCTION_BC_ENERGY) += energy;
00324     reduction->submit();
00325 
00326 }


Member Data Documentation

SubmitReduction* ComputeCylindricalBC::reduction

Definition at line 37 of file ComputeCylindricalBC.h.

Referenced by ComputeCylindricalBC(), doForce(), and ~ComputeCylindricalBC().


The documentation for this class was generated from the following files:
Generated on Sun Sep 24 01:17:16 2017 for NAMD by  doxygen 1.4.7