OOPSE/libmdtools/DirectionalAtom.cpp

#include <cmath>

#include "Atom.hpp"
#include "simError.h"

void DirectionalAtom::zeroForces() {
  if( hasCoords ){
    frc[offsetX] = 0.0; 
    frc[offsetY] = 0.0; 
    frc[offsetZ] = 0.0;
    
    trq[offsetX] = 0.0; 
    trq[offsetY] = 0.0; 
    trq[offsetZ] = 0.0;
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to zero frc and trq for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }
}

void DirectionalAtom::setCoords(void){

  if( myConfig->isAllocated() ){

    myConfig->getAtomPointers( index,
                     &pos, 
                     &vel, 
                     &frc, 
                     &trq, 
                     &Amat,
                     &mu,  
                     &ul );
  }
  else{
    sprintf( painCave.errMsg,
             "Attempted to set Atom %d  coordinates with an unallocated "
             "SimState object.\n" );
    painCave.isFatal = 1;
    simError();
  }

  hasCoords = true;

  *mu = myMu;

}

double DirectionalAtom::getMu( void ) {

  if( hasCoords ){
    return *mu;
  }
  else{
    return myMu;
  }
  return 0;
}

void DirectionalAtom::setMu( double the_mu ) { 

  if( hasCoords ){
    *mu = the_mu;
    myMu = the_mu;
  }
  else{
    myMu = the_mu;
  }
}

void DirectionalAtom::setA( double the_A[3][3] ){

  if( hasCoords ){
    Amat[Axx] = the_A[0][0]; Amat[Axy] = the_A[0][1]; Amat[Axz] = the_A[0][2];
    Amat[Ayx] = the_A[1][0]; Amat[Ayy] = the_A[1][1]; Amat[Ayz] = the_A[1][2];
    Amat[Azx] = the_A[2][0]; Amat[Azy] = the_A[2][1]; Amat[Azz] = the_A[2][2];
    
    this->updateU();  
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to set Amat for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }
}

void DirectionalAtom::setI( double the_I[3][3] ){
  
  Ixx = the_I[0][0]; Ixy = the_I[0][1]; Ixz = the_I[0][2];
  Iyx = the_I[1][0]; Iyy = the_I[1][1]; Iyz = the_I[1][2];
  Izx = the_I[2][0]; Izy = the_I[2][1]; Izz = the_I[2][2];
}

void DirectionalAtom::setQ( double the_q[4] ){

  double q0Sqr, q1Sqr, q2Sqr, q3Sqr;

  if( hasCoords ){
    q0Sqr = the_q[0] * the_q[0];
    q1Sqr = the_q[1] * the_q[1];
    q2Sqr = the_q[2] * the_q[2];
    q3Sqr = the_q[3] * the_q[3];
    
    
    Amat[Axx] = q0Sqr + q1Sqr - q2Sqr - q3Sqr;
    Amat[Axy] = 2.0 * ( the_q[1] * the_q[2] + the_q[0] * the_q[3] );
    Amat[Axz] = 2.0 * ( the_q[1] * the_q[3] - the_q[0] * the_q[2] );
    
    Amat[Ayx] = 2.0 * ( the_q[1] * the_q[2] - the_q[0] * the_q[3] );
    Amat[Ayy] = q0Sqr - q1Sqr + q2Sqr - q3Sqr;
    Amat[Ayz] = 2.0 * ( the_q[2] * the_q[3] + the_q[0] * the_q[1] );
    
    Amat[Azx] = 2.0 * ( the_q[1] * the_q[3] + the_q[0] * the_q[2] );
    Amat[Azy] = 2.0 * ( the_q[2] * the_q[3] - the_q[0] * the_q[1] );
    Amat[Azz] = q0Sqr - q1Sqr -q2Sqr +q3Sqr;
    
    this->updateU();
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to set Q for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }

}

void DirectionalAtom::getA( double the_A[3][3] ){
  
  if( hasCoords ){
    the_A[0][0] = Amat[Axx];
    the_A[0][1] = Amat[Axy];
    the_A[0][2] = Amat[Axz];
    
    the_A[1][0] = Amat[Ayx];
    the_A[1][1] = Amat[Ayy];
    the_A[1][2] = Amat[Ayz];
    
    the_A[2][0] = Amat[Azx];
    the_A[2][1] = Amat[Azy];
    the_A[2][2] = Amat[Azz];
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to get Amat for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }

}

void DirectionalAtom::printAmatIndex( void ){

  if( hasCoords ){
    std::cerr << "Atom[" << index << "] index =>\n" 
              << "[ " << Axx << ", " << Axy << ", " << Axz << " ]\n"
              << "[ " << Ayx << ", " << Ayy << ", " << Ayz << " ]\n"
              << "[ " << Azx << ", " << Azy << ", " << Azz << " ]\n";
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to print Amat indices for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }
}


void DirectionalAtom::getU( double the_u[3] ){
  
  the_u[0] = sux;
  the_u[1] = suy;
  the_u[2] = suz;

  this->body2Lab( the_u );
}

void DirectionalAtom::getQ( double q[4] ){
  
  double t, s;
  double ad1, ad2, ad3;

  if( hasCoords ){
    
    t = Amat[Axx] + Amat[Ayy] + Amat[Azz] + 1.0;
    if( t > 0.0 ){
      
      s = 0.5 / sqrt( t );
      q[0] = 0.25 / s;
      q[1] = (Amat[Ayz] - Amat[Azy]) * s;
      q[2] = (Amat[Azx] - Amat[Axz]) * s;
      q[3] = (Amat[Axy] - Amat[Ayx]) * s;
    }
    else{
      
      ad1 = fabs( Amat[Axx] );
      ad2 = fabs( Amat[Ayy] );
      ad3 = fabs( Amat[Azz] );
      
      if( ad1 >= ad2 && ad1 >= ad3 ){
        
        s = 2.0 * sqrt( 1.0 + Amat[Axx] - Amat[Ayy] - Amat[Azz] );
        q[0] = (Amat[Ayz] + Amat[Azy]) / s;
        q[1] = 0.5 / s;
        q[2] = (Amat[Axy] + Amat[Ayx]) / s;
        q[3] = (Amat[Axz] + Amat[Azx]) / s;
      }
      else if( ad2 >= ad1 && ad2 >= ad3 ){
        
        s = sqrt( 1.0 + Amat[Ayy] - Amat[Axx] - Amat[Azz] ) * 2.0;
        q[0] = (Amat[Axz] + Amat[Azx]) / s;
        q[1] = (Amat[Axy] + Amat[Ayx]) / s;
        q[2] = 0.5 / s;
        q[3] = (Amat[Ayz] + Amat[Azy]) / s;
      }
      else{
        
        s = sqrt( 1.0 + Amat[Azz] - Amat[Axx] - Amat[Ayy] ) * 2.0;
        q[0] = (Amat[Axy] + Amat[Ayx]) / s;
        q[1] = (Amat[Axz] + Amat[Azx]) / s;
        q[2] = (Amat[Ayz] + Amat[Azy]) / s;
        q[3] = 0.5 / s;
      }
    }
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to get Q for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }
}


void DirectionalAtom::setEuler( double phi, double theta, double psi ){
  
  if( hasCoords ){
    Amat[Axx] = (cos(phi) * cos(psi)) - (sin(phi) * cos(theta) * sin(psi));
    Amat[Axy] = (sin(phi) * cos(psi)) + (cos(phi) * cos(theta) * sin(psi));
    Amat[Axz] = sin(theta) * sin(psi);
    
    Amat[Ayx] = -(cos(phi) * sin(psi)) - (sin(phi) * cos(theta) * cos(psi));
    Amat[Ayy] = -(sin(phi) * sin(psi)) + (cos(phi) * cos(theta) * cos(psi));
    Amat[Ayz] = sin(theta) * cos(psi);
    
    Amat[Azx] = sin(phi) * sin(theta);
    Amat[Azy] = -cos(phi) * sin(theta);
    Amat[Azz] = cos(theta);
    
    this->updateU();
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to set Euler angles for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }
}


void DirectionalAtom::lab2Body( double r[3] ){

  double rl[3]; // the lab frame vector 
  
  if( hasCoords ){
    rl[0] = r[0];
    rl[1] = r[1];
    rl[2] = r[2];
    
    r[0] = (Amat[Axx] * rl[0]) + (Amat[Axy] * rl[1]) + (Amat[Axz] * rl[2]);
    r[1] = (Amat[Ayx] * rl[0]) + (Amat[Ayy] * rl[1]) + (Amat[Ayz] * rl[2]);
    r[2] = (Amat[Azx] * rl[0]) + (Amat[Azy] * rl[1]) + (Amat[Azz] * rl[2]);
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to convert lab2body for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }

}

void DirectionalAtom::body2Lab( double r[3] ){

  double rb[3]; // the body frame vector 
  
  if( hasCoords ){
    rb[0] = r[0];
    rb[1] = r[1];
    rb[2] = r[2];
    
    r[0] = (Amat[Axx] * rb[0]) + (Amat[Ayx] * rb[1]) + (Amat[Azx] * rb[2]);
    r[1] = (Amat[Axy] * rb[0]) + (Amat[Ayy] * rb[1]) + (Amat[Azy] * rb[2]);
    r[2] = (Amat[Axz] * rb[0]) + (Amat[Ayz] * rb[1]) + (Amat[Azz] * rb[2]);
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to convert body2lab for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }
}

void DirectionalAtom::updateU( void ){

  if( hasCoords ){
    ul[offsetX] = (Amat[Axx] * sux) + (Amat[Ayx] * suy) + (Amat[Azx] * suz);
    ul[offsetY] = (Amat[Axy] * sux) + (Amat[Ayy] * suy) + (Amat[Azy] * suz);
    ul[offsetZ] = (Amat[Axz] * sux) + (Amat[Ayz] * suy) + (Amat[Azz] * suz);
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to updateU for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }
}

void DirectionalAtom::getJ( double theJ[3] ){
  
  theJ[0] = jx;
  theJ[1] = jy;
  theJ[2] = jz;
}

void DirectionalAtom::setJ( double theJ[3] ){
  
  jx = theJ[0];
  jy = theJ[1];
  jz = theJ[2];
}

void DirectionalAtom::getTrq( double theT[3] ){
  
  if( hasCoords ){
    theT[0] = trq[offsetX];
    theT[1] = trq[offsetY];
    theT[2] = trq[offsetZ];
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to get Trq for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }
}

void DirectionalAtom::addTrq( double theT[3] ){
  
  if( hasCoords ){
    trq[offsetX] += theT[0];
    trq[offsetY] += theT[1];
    trq[offsetZ] += theT[2];
  }
  else{
    
    sprintf( painCave.errMsg,
             "Attempt to add Trq for atom %d before coords set.\n",
             index );
    painCave.isFatal = 1;
    simError();
  }
}


void DirectionalAtom::getI( double the_I[3][3] ){
  
  the_I[0][0] = Ixx;
  the_I[0][1] = Ixy;
  the_I[0][2] = Ixz;

  the_I[1][0] = Iyx;
  the_I[1][1] = Iyy;
  the_I[1][2] = Iyz;

  the_I[2][0] = Izx;
  the_I[2][1] = Izy;
  the_I[2][2] = Izz;
}
Revision:	690
Committed:	Tue Aug 12 21:44:06 2003 UTC (22 years, 2 months ago) by mmeineke
Original Path:	*trunk/OOPSE/libmdtools/DirectionalAtom.cpp*
File size:	8918 byte(s)
Log Message:	fixed a really annoying bug in Directional Atom, where mu was getting written to pseudorandom memory location.
#	User	Rev	Content
1	mmeineke	377	#include <cmath>
2
3			#include "Atom.hpp"
4	mmeineke	670	#include "simError.h"
5	mmeineke	377
6	mmeineke	670	void DirectionalAtom::zeroForces() {
7			if( hasCoords ){
8			frc[offsetX] = 0.0;
9			frc[offsetY] = 0.0;
10			frc[offsetZ] = 0.0;
11
12			trq[offsetX] = 0.0;
13			trq[offsetY] = 0.0;
14			trq[offsetZ] = 0.0;
15			}
16			else{
17
18			sprintf( painCave.errMsg,
19			"Attempt to zero frc and trq for atom %d before coords set.\n",
20			index );
21			painCave.isFatal = 1;
22			simError();
23			}
24			}
25	mmeineke	377
26	tim	689	void DirectionalAtom::setCoords(void){
27	mmeineke	414
28	tim	689	if( myConfig->isAllocated() ){
29
30			myConfig->getAtomPointers( index,
31			&pos,
32			&vel,
33			&frc,
34			&trq,
35			&Amat,
36			&mu,
37			&ul );
38	mmeineke	670	}
39			else{
40			sprintf( painCave.errMsg,
41	tim	689	"Attempted to set Atom %d coordinates with an unallocated "
42			"SimState object.\n" );
43	mmeineke	670	painCave.isFatal = 1;
44			simError();
45			}
46	tim	689
47			hasCoords = true;
48
49	mmeineke	690	*mu = myMu;
50	tim	689
51			}
52
53			double DirectionalAtom::getMu( void ) {
54
55			if( hasCoords ){
56	mmeineke	690	return *mu;
57	tim	689	}
58			else{
59			return myMu;
60			}
61	mmeineke	670	return 0;
62			}
63
64			void DirectionalAtom::setMu( double the_mu ) {
65
66			if( hasCoords ){
67	mmeineke	690	*mu = the_mu;
68	tim	689	myMu = the_mu;
69	mmeineke	670	}
70			else{
71	tim	689	myMu = the_mu;
72	mmeineke	670	}
73			}
74
75	mmeineke	377	void DirectionalAtom::setA( double the_A[3][3] ){
76
77	mmeineke	670	if( hasCoords ){
78			Amat[Axx] = the_A[0][0]; Amat[Axy] = the_A[0][1]; Amat[Axz] = the_A[0][2];
79			Amat[Ayx] = the_A[1][0]; Amat[Ayy] = the_A[1][1]; Amat[Ayz] = the_A[1][2];
80			Amat[Azx] = the_A[2][0]; Amat[Azy] = the_A[2][1]; Amat[Azz] = the_A[2][2];
81
82			this->updateU();
83			}
84			else{
85
86			sprintf( painCave.errMsg,
87			"Attempt to set Amat for atom %d before coords set.\n",
88			index );
89			painCave.isFatal = 1;
90			simError();
91			}
92	mmeineke	377	}
93
94			void DirectionalAtom::setI( double the_I[3][3] ){
95
96			Ixx = the_I[0][0]; Ixy = the_I[0][1]; Ixz = the_I[0][2];
97			Iyx = the_I[1][0]; Iyy = the_I[1][1]; Iyz = the_I[1][2];
98			Izx = the_I[2][0]; Izy = the_I[2][1]; Izz = the_I[2][2];
99			}
100
101			void DirectionalAtom::setQ( double the_q[4] ){
102
103			double q0Sqr, q1Sqr, q2Sqr, q3Sqr;
104
105	mmeineke	670	if( hasCoords ){
106			q0Sqr = the_q[0] * the_q[0];
107			q1Sqr = the_q[1] * the_q[1];
108			q2Sqr = the_q[2] * the_q[2];
109			q3Sqr = the_q[3] * the_q[3];
110
111
112			Amat[Axx] = q0Sqr + q1Sqr - q2Sqr - q3Sqr;
113			Amat[Axy] = 2.0 * ( the_q[1] * the_q[2] + the_q[0] * the_q[3] );
114			Amat[Axz] = 2.0 * ( the_q[1] * the_q[3] - the_q[0] * the_q[2] );
115
116			Amat[Ayx] = 2.0 * ( the_q[1] * the_q[2] - the_q[0] * the_q[3] );
117			Amat[Ayy] = q0Sqr - q1Sqr + q2Sqr - q3Sqr;
118			Amat[Ayz] = 2.0 * ( the_q[2] * the_q[3] + the_q[0] * the_q[1] );
119
120			Amat[Azx] = 2.0 * ( the_q[1] * the_q[3] + the_q[0] * the_q[2] );
121			Amat[Azy] = 2.0 * ( the_q[2] * the_q[3] - the_q[0] * the_q[1] );
122			Amat[Azz] = q0Sqr - q1Sqr -q2Sqr +q3Sqr;
123
124			this->updateU();
125			}
126			else{
127
128			sprintf( painCave.errMsg,
129			"Attempt to set Q for atom %d before coords set.\n",
130			index );
131			painCave.isFatal = 1;
132			simError();
133			}
134	mmeineke	377
135			}
136
137			void DirectionalAtom::getA( double the_A[3][3] ){
138
139	mmeineke	670	if( hasCoords ){
140			the_A[0][0] = Amat[Axx];
141			the_A[0][1] = Amat[Axy];
142			the_A[0][2] = Amat[Axz];
143
144			the_A[1][0] = Amat[Ayx];
145			the_A[1][1] = Amat[Ayy];
146			the_A[1][2] = Amat[Ayz];
147
148			the_A[2][0] = Amat[Azx];
149			the_A[2][1] = Amat[Azy];
150			the_A[2][2] = Amat[Azz];
151			}
152			else{
153
154			sprintf( painCave.errMsg,
155			"Attempt to get Amat for atom %d before coords set.\n",
156			index );
157			painCave.isFatal = 1;
158			simError();
159			}
160	mmeineke	377
161			}
162
163	mmeineke	597	void DirectionalAtom::printAmatIndex( void ){
164	mmeineke	377
165	mmeineke	670	if( hasCoords ){
166			std::cerr << "Atom[" << index << "] index =>\n"
167			<< "[ " << Axx << ", " << Axy << ", " << Axz << " ]\n"
168			<< "[ " << Ayx << ", " << Ayy << ", " << Ayz << " ]\n"
169			<< "[ " << Azx << ", " << Azy << ", " << Azz << " ]\n";
170			}
171			else{
172
173			sprintf( painCave.errMsg,
174			"Attempt to print Amat indices for atom %d before coords set.\n",
175			index );
176			painCave.isFatal = 1;
177			simError();
178			}
179	mmeineke	597	}
180
181
182	mmeineke	377	void DirectionalAtom::getU( double the_u[3] ){
183
184			the_u[0] = sux;
185			the_u[1] = suy;
186			the_u[2] = suz;
187
188			this->body2Lab( the_u );
189			}
190
191			void DirectionalAtom::getQ( double q[4] ){
192
193			double t, s;
194			double ad1, ad2, ad3;
195
196	mmeineke	670	if( hasCoords ){
197	mmeineke	377
198	mmeineke	670	t = Amat[Axx] + Amat[Ayy] + Amat[Azz] + 1.0;
199			if( t > 0.0 ){
200	mmeineke	377
201	mmeineke	670	s = 0.5 / sqrt( t );
202			q[0] = 0.25 / s;
203			q[1] = (Amat[Ayz] - Amat[Azy]) * s;
204			q[2] = (Amat[Azx] - Amat[Axz]) * s;
205			q[3] = (Amat[Axy] - Amat[Ayx]) * s;
206	mmeineke	377	}
207			else{
208
209	mmeineke	670	ad1 = fabs( Amat[Axx] );
210			ad2 = fabs( Amat[Ayy] );
211			ad3 = fabs( Amat[Azz] );
212
213			if( ad1 >= ad2 && ad1 >= ad3 ){
214
215			s = 2.0 * sqrt( 1.0 + Amat[Axx] - Amat[Ayy] - Amat[Azz] );
216			q[0] = (Amat[Ayz] + Amat[Azy]) / s;
217			q[1] = 0.5 / s;
218			q[2] = (Amat[Axy] + Amat[Ayx]) / s;
219			q[3] = (Amat[Axz] + Amat[Azx]) / s;
220			}
221			else if( ad2 >= ad1 && ad2 >= ad3 ){
222
223			s = sqrt( 1.0 + Amat[Ayy] - Amat[Axx] - Amat[Azz] ) * 2.0;
224			q[0] = (Amat[Axz] + Amat[Azx]) / s;
225			q[1] = (Amat[Axy] + Amat[Ayx]) / s;
226			q[2] = 0.5 / s;
227			q[3] = (Amat[Ayz] + Amat[Azy]) / s;
228			}
229			else{
230
231			s = sqrt( 1.0 + Amat[Azz] - Amat[Axx] - Amat[Ayy] ) * 2.0;
232			q[0] = (Amat[Axy] + Amat[Ayx]) / s;
233			q[1] = (Amat[Axz] + Amat[Azx]) / s;
234			q[2] = (Amat[Ayz] + Amat[Azy]) / s;
235			q[3] = 0.5 / s;
236			}
237	mmeineke	377	}
238			}
239	mmeineke	670	else{
240
241			sprintf( painCave.errMsg,
242			"Attempt to get Q for atom %d before coords set.\n",
243			index );
244			painCave.isFatal = 1;
245			simError();
246			}
247	mmeineke	377	}
248
249
250			void DirectionalAtom::setEuler( double phi, double theta, double psi ){
251
252	mmeineke	670	if( hasCoords ){
253			Amat[Axx] = (cos(phi) * cos(psi)) - (sin(phi) * cos(theta) * sin(psi));
254			Amat[Axy] = (sin(phi) * cos(psi)) + (cos(phi) * cos(theta) * sin(psi));
255			Amat[Axz] = sin(theta) * sin(psi);
256
257			Amat[Ayx] = -(cos(phi) * sin(psi)) - (sin(phi) * cos(theta) * cos(psi));
258			Amat[Ayy] = -(sin(phi) * sin(psi)) + (cos(phi) * cos(theta) * cos(psi));
259			Amat[Ayz] = sin(theta) * cos(psi);
260
261			Amat[Azx] = sin(phi) * sin(theta);
262			Amat[Azy] = -cos(phi) * sin(theta);
263			Amat[Azz] = cos(theta);
264
265			this->updateU();
266			}
267			else{
268
269			sprintf( painCave.errMsg,
270			"Attempt to set Euler angles for atom %d before coords set.\n",
271			index );
272			painCave.isFatal = 1;
273			simError();
274			}
275	mmeineke	377	}
276
277
278			void DirectionalAtom::lab2Body( double r[3] ){
279
280			double rl[3]; // the lab frame vector
281
282	mmeineke	670	if( hasCoords ){
283			rl[0] = r[0];
284			rl[1] = r[1];
285			rl[2] = r[2];
286
287			r[0] = (Amat[Axx] * rl[0]) + (Amat[Axy] * rl[1]) + (Amat[Axz] * rl[2]);
288			r[1] = (Amat[Ayx] * rl[0]) + (Amat[Ayy] * rl[1]) + (Amat[Ayz] * rl[2]);
289			r[2] = (Amat[Azx] * rl[0]) + (Amat[Azy] * rl[1]) + (Amat[Azz] * rl[2]);
290			}
291			else{
292
293			sprintf( painCave.errMsg,
294			"Attempt to convert lab2body for atom %d before coords set.\n",
295			index );
296			painCave.isFatal = 1;
297			simError();
298			}
299	mmeineke	377
300			}
301
302			void DirectionalAtom::body2Lab( double r[3] ){
303
304			double rb[3]; // the body frame vector
305
306	mmeineke	670	if( hasCoords ){
307			rb[0] = r[0];
308			rb[1] = r[1];
309			rb[2] = r[2];
310
311			r[0] = (Amat[Axx] * rb[0]) + (Amat[Ayx] * rb[1]) + (Amat[Azx] * rb[2]);
312			r[1] = (Amat[Axy] * rb[0]) + (Amat[Ayy] * rb[1]) + (Amat[Azy] * rb[2]);
313			r[2] = (Amat[Axz] * rb[0]) + (Amat[Ayz] * rb[1]) + (Amat[Azz] * rb[2]);
314			}
315			else{
316
317			sprintf( painCave.errMsg,
318			"Attempt to convert body2lab for atom %d before coords set.\n",
319			index );
320			painCave.isFatal = 1;
321			simError();
322			}
323	mmeineke	377	}
324
325			void DirectionalAtom::updateU( void ){
326
327	mmeineke	670	if( hasCoords ){
328			ul[offsetX] = (Amat[Axx] * sux) + (Amat[Ayx] * suy) + (Amat[Azx] * suz);
329			ul[offsetY] = (Amat[Axy] * sux) + (Amat[Ayy] * suy) + (Amat[Azy] * suz);
330			ul[offsetZ] = (Amat[Axz] * sux) + (Amat[Ayz] * suy) + (Amat[Azz] * suz);
331			}
332			else{
333
334			sprintf( painCave.errMsg,
335			"Attempt to updateU for atom %d before coords set.\n",
336			index );
337			painCave.isFatal = 1;
338			simError();
339			}
340	mmeineke	377	}
341
342	mmeineke	599	void DirectionalAtom::getJ( double theJ[3] ){
343
344			theJ[0] = jx;
345			theJ[1] = jy;
346			theJ[2] = jz;
347			}
348
349			void DirectionalAtom::setJ( double theJ[3] ){
350
351			jx = theJ[0];
352			jy = theJ[1];
353			jz = theJ[2];
354			}
355
356			void DirectionalAtom::getTrq( double theT[3] ){
357
358	mmeineke	670	if( hasCoords ){
359			theT[0] = trq[offsetX];
360			theT[1] = trq[offsetY];
361			theT[2] = trq[offsetZ];
362			}
363			else{
364
365			sprintf( painCave.errMsg,
366			"Attempt to get Trq for atom %d before coords set.\n",
367			index );
368			painCave.isFatal = 1;
369			simError();
370			}
371	mmeineke	599	}
372
373			void DirectionalAtom::addTrq( double theT[3] ){
374
375	mmeineke	670	if( hasCoords ){
376			trq[offsetX] += theT[0];
377			trq[offsetY] += theT[1];
378			trq[offsetZ] += theT[2];
379			}
380			else{
381
382			sprintf( painCave.errMsg,
383			"Attempt to add Trq for atom %d before coords set.\n",
384			index );
385			painCave.isFatal = 1;
386			simError();
387			}
388	mmeineke	599	}
389
390
391			void DirectionalAtom::getI( double the_I[3][3] ){
392
393			the_I[0][0] = Ixx;
394			the_I[0][1] = Ixy;
395			the_I[0][2] = Ixz;
396
397			the_I[1][0] = Iyx;
398			the_I[1][1] = Iyy;
399			the_I[1][2] = Iyz;
400
401			the_I[2][0] = Izx;
402			the_I[2][1] = Izy;
403			the_I[2][2] = Izz;
404			}