ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/tcProps/directorWhole.c
(Generate patch)

Comparing trunk/tcProps/directorWhole.c (file contents):
Revision 1062 by mmeineke, Fri Feb 20 21:20:37 2004 UTC vs.
Revision 1073 by mmeineke, Sat Feb 28 16:45:57 2004 UTC

# Line 74 | Line 74 | void calcDirWholeCorr(double startTime, struct atomCoo
74      percentComplete =
75        (int)( 100.0 * (double)framesFinished / (double) corrFrames );
76      
77 < //     fprintf( stdout,
78 < //           "\rDirector bilayer corr %3d%% complete.",
79 < //           percentComplete );
80 < //     fflush( stdout );
77 >    fprintf( stdout,
78 >             "\rDirector bilayer corr %3d%% complete.",
79 >             percentComplete );
80 >    fflush( stdout );
81  
82      readFrame( i, atoms, Hmat );
83      
# Line 88 | Line 88 | void calcDirWholeCorr(double startTime, struct atomCoo
88      framesFinished++;
89    }
90  
91 < //   sprintf( outName, "%s.dirHead", outPrefix );
92 < //   outFile = fopen( outName, "w" );
91 >  sprintf( outName, "%s.dirWhole", outPrefix );
92 >  outFile = fopen( outName, "w" );
93  
94 +  fprintf( outFile,
95 +           "#time\torderParam\tx\ty\tz\n");
96 +  
97 +  for(i=0;i<corrFrames;i++){
98 +    fprintf( outFile,
99 +             "%6G\t%6G\t%6G\t%6G\t%6G\n",
100 +             directorWhole[i].time,
101 +             directorWhole[i].order,
102 +             directorWhole[i].u[0],
103 +             directorWhole[i].u[1],
104 +             directorWhole[i].u[2]);
105 +  }
106  
107    
108  
109 < //   fflush(outFile);
110 < //   fclose(outFile);
109 >  fflush(outFile);
110 >  fclose(outFile);
111  
112    percentComplete =
113      (int)( 100.0 * (double)framesFinished / (double) corrFrames );
# Line 105 | Line 117 | void calcDirWholeCorr(double startTime, struct atomCoo
117             "done.\n",
118             percentComplete );
119    fflush( stdout );
120 +
121 +  free(directorWhole);
122 +  directorWhole = NULL;
123   }
124  
125  
# Line 127 | Line 142 | void accumDWFrame( int index, struct atomCoord *atoms
142  
143    double oMat[3][3];  
144    double iMat[3][3];
145 +  double matWrap[9];
146    double com[3];
147    double totMass;
148    double lenPrinc, distl;
# Line 172 | Line 188 | void accumDWFrame( int index, struct atomCoord *atoms
188      for(j=0;j<nLipAtoms;j++){
189        l=k+j;
190        
191 <      for(m=0;m<3;m++){
191 >      for(m=0;m<3;m++)
192          com[m] += atoms[l].mass * atoms[l].pos[m];
193 <        totMass += atoms[l].mass;
194 <      }
193 >      
194 >      totMass += atoms[l].mass;
195      }
196      for(m=0;m<3;m++) com[m] /= totMass;
197  
# Line 197 | Line 213 | void accumDWFrame( int index, struct atomCoord *atoms
213        
214        
215      }
216 +    
217 + //     fprintf(stderr, "iMat (before) = %6G, %6G, %6G\n"
218 + //                  "                %6G, %6G, %6G\n"
219 + //                  "                %6G, %6G, %6G\n\n",
220 + //          iMat[0][0], iMat[0][1], iMat[0][2],
221 + //          iMat[1][0], iMat[1][1], iMat[1][2],
222 + //          iMat[2][0], iMat[2][1], iMat[2][2]);
223  
224 +    // matWrap to fortran convention
225 +
226 +    for(j=0;j<3;j++)
227 +      for(l=0;l<3;l++)
228 +        matWrap[l+3*j] = iMat[l][j];
229 +
230      ifail = 0;
231 <    dsyev(&job, &uplo, &nfilled, iMat, &ndiag, evals, work, &lWork, &ifail);
231 >    dsyev(&job, &uplo, &nfilled, matWrap, &ndiag, evals, work, &lWork, &ifail);
232      
233      if (ifail) {
234        fprintf(stderr, "dsyev screwed something up!\n");
235        exit(0);
236      }
237  
238 <    fprintf(stderr, "evals = %6G, %6G, %6G\n\n",
210 <            evals[0], evals[1], evals[2]);
238 >    // matWrap from fortran convention
239  
240 <    fprintf(stderr, "iMat = %6G, %6G, %6G\n"
241 <                    "       %6G, %6G, %6G\n"
242 <                    "       %6G, %6G, %6G\n\n",
215 <            iMat[0][0], iMat[0][1], iMat[0][2],
216 <            iMat[1][0], iMat[1][1], iMat[1][2],
217 <            iMat[2][0], iMat[2][1], iMat[2][2]);
240 >    for(j=0;j<3;j++)
241 >      for(l=0;l<3;l++)
242 >        iMat[l][j] = matWrap[l+3*j];
243  
244 + //     fprintf(stderr, "evals = %6G, %6G, %6G\n\n",
245 + //          evals[0], evals[1], evals[2]);
246  
247 + //     fprintf(stderr, "iMat (after) = %6G, %6G, %6G\n"
248 + //                  "               %6G, %6G, %6G\n"
249 + //                  "               %6G, %6G, %6G\n\n",
250 + //          iMat[0][0], iMat[0][1], iMat[0][2],
251 + //          iMat[1][0], iMat[1][1], iMat[1][2],
252 + //          iMat[2][0], iMat[2][1], iMat[2][2]);
253 +
254 +
255        
256      smallest = fabs(evals[0]);
257      which = 0;
# Line 230 | Line 265 | void accumDWFrame( int index, struct atomCoord *atoms
265      lenPrinc = 0.0;
266      for(j=0;j<3;j++)
267        lenPrinc +=  iMat[j][which] * iMat[j][which];
268 +    lenPrinc = sqrt( lenPrinc );
269      
270      for(j=0;j<3;j++)
271        principal[i].u[j] = iMat[j][which]/lenPrinc;
# Line 259 | Line 295 | void accumDWFrame( int index, struct atomCoord *atoms
295      for(j=0;j<3;j++)
296        oMat[i][j] /= (double)nLipids;
297  
298 +  // matWrap to fortran convention
299 +  
300 +  for(j=0;j<3;j++)
301 +    for(l=0;l<3;l++)
302 +      matWrap[l+3*j] = oMat[l][j];
303 +  
304    ifail = 0;
305 +  dsyev(&job, &uplo, &nfilled, matWrap, &ndiag, evals, work, &lWork, &ifail);
306  
264 –  dsyev(&job, &uplo, &nfilled, oMat, &ndiag, evals, work, &lWork, &ifail);
265 –
307    if (ifail) {
308      fprintf(stderr, "dsyev screwed something up!\n");
309      exit(0);
310    }
311 +
312 +  // matWrap from fortran convention
313    
314 +  for(j=0;j<3;j++)
315 +    for(l=0;l<3;l++)
316 +      oMat[l][j] = matWrap[l+3*j];
317 +  
318    max = 0.0;
319    for (i=0; i<3;i++) {
320      if (fabs(evals[i]) > max) {
# Line 282 | Line 329 | void accumDWFrame( int index, struct atomCoord *atoms
329  
330    directorWhole[index].order = 1.5 * max;
331  
332 <  fprintf(stderr,
333 <          "frame[%d] => order = %6G; < %6G, %6G, %6G >\n",
334 <          index,
335 <          directorWhole[index].order,
336 <          directorWhole[index].u[0],
337 <          directorWhole[index].u[1],
338 <          directorWhole[index].u[2] );
332 > //   fprintf(stderr,
333 > //        "frame[%d] => order = %6G; < %6G, %6G, %6G >\n",
334 > //        index,
335 > //        directorWhole[index].order,
336 > //        directorWhole[index].u[0],
337 > //        directorWhole[index].u[1],
338 > //        directorWhole[index].u[2] );
339  
340   //   fprintf(stdout,
341   //        "%6g\t%6G\n",

Diff Legend

– Removed lines
+ Added lines
< Changed lines
> Changed lines