psi::extrema::internals Class Reference

Abstract internal coordinate derived class. More...

#include <internals.h>

Inheritance diagram for psi::extrema::internals:

psi::extrema::coord_base psi::extrema::coord_base_carts psi::extrema::math_tools psi::extrema::deloc psi::extrema::zmat

Protected Member Functions

 internals ()
 Dummy constructor.
void mem_alloc ()
 Allocates memory.
 ~internals ()
 Destructor. Frees memory.
virtual void compute_B ()=0
 Derived classes must provide a function to compute the B matrix though calls to B_row...() functions.
virtual void cart_to_internal (double **)=0
void back_transform (double *, double *)
 Performs iterative back transformation from internals to cartesians.
double * B_row_bond (double *, int, int)
 Computes a row of B corrresponding to a simple bonding coordinate.
double * B_row_angle (double *, int, int, int)
 Computes a row of B corresponding to a simple bending coordinate.
double * B_row_tors (double *, int, int, int, int)
 Computes row of B matrix corresponding to a simple torsion coordinate.
void print_B ()
 Prints the B matrix.
void compute_G ()
 Computes the G=B.B^t matrix.
void compute_A ()
 Computes the A=u.B^t.G matrix. ---------------------------------------------------------------------------.
virtual void grad_trans ()
 Performs gradient transformation from cartesians to internals.

Protected Attributes

int fnum_coords
int B_dim
double ** full_geom
double * fcoords
double * fcoords_old
double * fgrads
double ** B
double ** B_red
double ** G
double ** A
double ** u

Detailed Description

Abstract internal coordinate derived class.

Member functions and data common to internal coordinate types. Primarily, this class contains the B matrix and related matrices, functions for computing these matrices, and functions which apply these matrices in gradient and coordinate transformations

"full" variables refer to the complete set of coordinates, including any equivalent coordinates. In the iterative back transformation to cartesians, the complete set of coordinates is utilized.

Definition at line 28 of file extrema/Internals/internals.h.


Constructor & Destructor Documentation

psi::extrema::internals::internals (  )  [inline, protected]

Dummy constructor.

This is a dummy constructor which exists to make compilers happy. It does nothing. The actual constructor may depend on data which is not know when a derived class is initialized and may be called later in that class's constructor.

Definition at line 44 of file extrema/Internals/internals.h.

00044 { return; };


Member Function Documentation

void internals::mem_alloc (  )  [protected]

Allocates memory.

Must be called by top-level class once dimensions determined.

Reimplemented from psi::extrema::coord_base.

Definition at line 25 of file extrema/Internals/internals.cc.

References A, B, B_red, fcoords, fcoords_old, fgrads, fnum_coords, full_geom, G, init_array(), init_matrix(), psi::extrema::coord_base::mem_alloc(), psi::extrema::coord_base::num_coords, psi::extrema::coord_base_carts::num_entries, and u.

Referenced by psi::extrema::deloc::deloc(), and psi::extrema::zmat::zmat().

00025                           {
00026 
00027     coord_base::mem_alloc();
00028 
00029     full_geom = init_matrix(num_entries,3);
00030     fgrads = init_array(fnum_coords);
00031     fcoords = init_array(fnum_coords);
00032     fcoords_old = init_array(fnum_coords);
00033     B = (double**) malloc(fnum_coords*sizeof(double*));
00034     B_red = init_matrix(num_coords,3*num_entries);
00035     G = init_matrix(fnum_coords, fnum_coords);
00036     A = init_matrix(3*num_entries, fnum_coords);
00037     u = init_matrix(3*num_entries, 3*num_entries);
00038     
00039     return;
00040 }

void internals::back_transform ( double *  c_new,
double *  c_old 
) [protected]

Performs iterative back transformation from internals to cartesians.

This function requires the full array of coordinate values.

Parameters:
*c_new new internal coordinate array
*c_old old internal coordinate array

Definition at line 84 of file extrema/Internals/internals.cc.

References A, ANOTHER_ZERO, BT_CONV, BT_LOOP, psi::extrema::coord_base_carts::carts, compute_A(), compute_B(), compute_G(), fnum_coords, init_array(), NORMAL_PRINT, psi::extrema::coord_base_carts::num_entries, psi::extrema::coord_base_carts::print_carts(), psi::extrema::coord_base::print_lvl, and RIDICULOUS_PRINT.

Referenced by psi::extrema::zmat::back_transform().

00084                                                               {
00085 
00086   int i, j, pos;
00087   double conv=1.0, *dq, *dx;
00088 
00089   dq = init_array(fnum_coords);
00090   dx = init_array(3*num_entries);         
00091 
00092   int loop=0;
00093   int converged=0;
00094   double criteria = 1.0e-14;
00095   double dx_sum; 
00096 
00097   fprintf(outfile,"\n\n\n");
00098   fprintf(outfile,"  --------------------------------------");
00099   fprintf(outfile,"--------------------------------------\n");
00100   fprintf(outfile,"  Performing iterative transformation to find");
00101   fprintf(outfile," new cartesian coordinates\n");
00102   fprintf(outfile,"  --------------------------------------");
00103   fprintf(outfile,"--------------------------------------");
00104   if(print_lvl > NORMAL_PRINT) {
00105       fprintf(outfile,"\n  Iter    dq (internals)          dx (cartesians)");
00106       fprintf(outfile,  "\n  ---- ----------------------  ");
00107       fprintf(outfile,"----------------------");  
00108   }
00109   else fprintf(outfile,"\n");
00110 
00111   if(print_lvl > RIDICULOUS_PRINT) 
00112       print_carts(1.0);
00113 
00114   while((!converged) && (loop<BT_LOOP) ) {
00115   
00116       /*compute A*/
00117       compute_G(); 
00118       compute_A(); 
00119       
00120       for(i=0;i<fnum_coords;++i) {
00121           dq[i] = c_new[i] - c_old[i];
00122       }
00123       if(print_lvl >= RIDICULOUS_PRINT) {
00124           for(i=0;i<fnum_coords;++i)
00125               fprintf(outfile,"\n dq[%d]=%.20lf",i+1,dq[i]);
00126       }
00127 
00128 
00129       /*compute dx = A dq */
00130       for(i=0;i<3*num_entries;++i) 
00131           dx[i]=0;
00132       for(i=0;i<3*num_entries;++i) {
00133           for(j=0;j<fnum_coords;++j) {
00134               dx[i] += A[i][j] * dq[j];
00135           }
00136       }
00137 
00138       pos=0;
00139       dx_sum = 0.0;
00140       int hack=0;
00141       for(i=0;i<3*num_entries;++i) {
00142           /* hack to keep proper orientation */
00143           if(fabs(carts[i])>ANOTHER_ZERO) 
00144               carts[i] += dx[i];
00145           else if(fabs(dx[i])>ANOTHER_ZERO)
00146               hack = 1;
00147           dx_sum += sqrt(dx[i]*dx[i]);
00148       }
00149       dx_sum /= (3*num_entries);
00150       if(hack)
00151           if(print_lvl >= RIDICULOUS_PRINT)
00152               fprintf(outfile,
00153                       "\n  WARNING: Using hack to keep proper orientation"); 
00154 
00155       if(print_lvl >= RIDICULOUS_PRINT) {
00156           for(i=0;i<3*num_entries;++i) 
00157               fprintf(outfile,"\n dx[%d]=%.20lf",i+1,dx[i]);
00158       }
00159       
00160       cart_to_internal(&c_old);
00161 
00162       compute_B();
00163     
00164       conv=0;
00165       for(i=0;i<fnum_coords;++i) {
00166           conv += sqrt((c_new[i] - c_old[i])*
00167                        (c_new[i] - c_old[i]));
00168       }
00169       conv /= fnum_coords;
00170       if(print_lvl>NORMAL_PRINT) {
00171           fprintf(outfile,"\n  %3d  %.20lf  %.20lf",loop+1,conv,dx_sum);
00172       }   
00173       if( (conv<BT_CONV) && (dx_sum<BT_CONV) )
00174           converged = 1;
00175       ++loop;
00176   }
00177   
00178   if(!converged) { 
00179       fprintf(outfile,
00180               "\n  Check for angles near 180.0 degrees, they're bad");
00181       fprintf(outfile,
00182               "\n  You may need to use z-matrix coordinates avoiding 180.0");
00183       fprintf(outfile," degree angles\n");
00184       punt("Back transformation to cartesians has failed");
00185   }
00186   else
00187       if(print_lvl > NORMAL_PRINT)
00188           fprintf(outfile,"\n  Back transformation to cartesians completed\n");
00189 
00190   
00191 
00192   return;
00193 }

double * internals::B_row_bond ( double *  carr,
int  atom1,
int  atom2 
) [protected]

Computes a row of B corrresponding to a simple bonding coordinate.

Parameters:
*carr full cartesian coordinate array
atom1 reference atom 1
atom2 atom bonded to 1

Definition at line 24 of file internals_B.cc.

References init_array(), and psi::extrema::coord_base_carts::num_entries.

Referenced by psi::extrema::zmat::compute_B().

00024                                                                     {
00025 
00026   int i;
00027   double *row, *unit12;
00028 
00029   row = init_array(num_entries*3);
00030 
00031   unit12 = unit_vec(c_arr, atom1, atom2);
00032 
00033   for(i=0;i<3;++i) {
00034       row[3*atom1+i] = -unit12[i];
00035       row[3*atom2+i] =  unit12[i];
00036     }
00037   return row;
00038 }

double * internals::B_row_angle ( double *  c_arr,
int  atom1,
int  atom3,
int  atom2 
) [protected]

Computes a row of B corresponding to a simple bending coordinate.

Parameters:
*c_arr full cartesian coordinate array
atom1 reference atom 1
atom3 atom bonded to 1
atom2 atom defining angle 1-3-2

Definition at line 51 of file internals_B.cc.

References init_array(), and psi::extrema::coord_base_carts::num_entries.

Referenced by psi::extrema::zmat::compute_B().

00052                                              {
00053 
00054   int i;
00055   double *row, *unit31, *unit32, r31, r32, cosine, sine;
00056  
00057   row = init_array(num_entries*3);
00058 
00059   unit31 = unit_vec(c_arr,atom3, atom1);
00060   unit32 = unit_vec(c_arr,atom3, atom2);
00061   cosine = dot_pdt(unit31,unit32);
00062   sine = sin(acos(cosine));
00063   r31 = vec_norm(c_arr,atom3,atom1);
00064   r32 = vec_norm(c_arr,atom3,atom2);
00065 
00066   for(i=0;i<3;++i) {
00067       row[3*atom1+i] = ( cosine * unit31[i] - unit32[i] ) / ( r31 * sine );
00068       row[3*atom2+i] = ( cosine * unit32[i] - unit31[i] ) / ( r32 * sine );
00069       row[3*atom3+i] = ((r31 - r32 * cosine) * unit31[i] + (r32 - r31 * cosine)
00070                         * unit32[i]) / (r31 * r32 * sine); 
00071     }
00072 
00073   return row;
00074 } 

double * internals::B_row_tors ( double *  c_arr,
int  atom1,
int  atom2,
int  atom3,
int  atom4 
) [protected]

Computes row of B matrix corresponding to a simple torsion coordinate.

Parameters:
*c_arr full cartesian array
atom1 reference atom 1
atom2 atom bonded to 1
atom3 atom defining angle 1-2-3
atom4 atom defining torsion 1-2-3-4

Definition at line 89 of file internals_B.cc.

References init_array(), and psi::extrema::coord_base_carts::num_entries.

Referenced by psi::extrema::zmat::compute_B().

00090                                                        {
00091 
00092   int i;
00093   double *row, *unit12, *unit23, *unit43, *unit32, *cross_12_23, *cross_43_32,
00094       cosine_an2, sine_an2, sine2_an2, 
00095       cosine_an3, sine_an3, sine2_an3, r12, r23;
00096 
00097   row = init_array(num_entries*3);
00098 
00099   unit12 = unit_vec(c_arr,atom1, atom2);
00100   unit23 = unit_vec(c_arr,atom2, atom3);
00101   unit43 = unit_vec(c_arr,atom4, atom3);
00102   unit32 = unit_vec(c_arr,atom3, atom2);
00103   cross_12_23 = cross_pdt(unit12, unit23);
00104   cross_43_32 = cross_pdt(unit43, unit32);
00105   cosine_an2 = dot_pdt( unit_vec(c_arr, atom2,atom1), 
00106                         unit_vec(c_arr, atom2,atom3));
00107   sine_an2 = sin(acos( cosine_an2 ));
00108   sine2_an2 = sine_an2 * sine_an2;
00109   cosine_an3 = dot_pdt( unit_vec(c_arr, atom3,atom2), 
00110                         unit_vec(c_arr, atom3,atom4));
00111   sine_an3 = sin(acos( cosine_an3 ));
00112   sine2_an3 = sine_an3 * sine_an3;
00113   r12 = vec_norm(c_arr,atom1,atom2);
00114   r23 = vec_norm(c_arr,atom2,atom3);
00115 
00116   for(i=0;i<3;++i) {
00117       row[3*atom1+i] = - cross_12_23[i] / ( r12 * sine2_an2 );
00118       row[3*atom2+i] = (r23 - r12 * cosine_an2) * cross_12_23[i] / 
00119           ( r23 * r12 * sine2_an2 ) +
00120           cosine_an3 * cross_43_32[i] / ( r23 * sine2_an3 );
00121     }
00122 
00123   
00124   /* entries for atom3 same as those for atom2 with permutation of 
00125      1 with 4 and 2 with 3, entries for atom4 same as those for 
00126      atom1 with same permutations*/
00127 
00128   /* I permute everything but the angles here */
00129   unit12 = unit_vec(c_arr,atom4, atom3);
00130   unit23 = unit_vec(c_arr,atom3, atom2);
00131   unit43 = unit_vec(c_arr,atom1, atom2);
00132   unit32 = unit_vec(c_arr,atom2, atom3);
00133   cross_12_23 = cross_pdt(unit12, unit23);
00134   cross_43_32 = cross_pdt(unit43, unit32);
00135   r12 = vec_norm(c_arr,atom4,atom3);
00136   r23 = vec_norm(c_arr,atom3,atom2);
00137 
00138   /* I permute the angles here */
00139   for(i=0;i<3;++i) {
00140       row[3*atom4+i] = - cross_12_23[i] / ( r12 * sine2_an3 );
00141       row[3*atom3+i] = (r23 - r12 * cosine_an3) * cross_12_23[i] / 
00142           ( r23 * r12 * sine2_an3 ) +
00143           cosine_an2 * cross_43_32[i] / ( r23 * sine2_an2 );
00144   }
00145   
00146   return row;
00147 } 


Field Documentation

int psi::extrema::internals::fnum_coords [protected]

full number of coordinates

Definition at line 32 of file extrema/Internals/internals.h.

Referenced by psi::extrema::zmat::back_transform(), back_transform(), compute_A(), psi::extrema::zmat::compute_B(), compute_G(), psi::extrema::deloc::deloc(), psi::extrema::zmat::grad_trans(), grad_trans(), psi::extrema::zmat::initial_Hi(), mem_alloc(), psi::extrema::zmat::newton_step(), print_B(), psi::extrema::zmat::zmat(), and ~internals().

int psi::extrema::internals::B_dim [protected]

fnum_coords for zmat, num_coords for deloc

Definition at line 32 of file extrema/Internals/internals.h.

double** psi::extrema::internals::full_geom [protected]

full geometry (including dummy atoms)

Definition at line 34 of file extrema/Internals/internals.h.

Referenced by mem_alloc(), and ~internals().

double * psi::extrema::internals::fcoords [protected]

full set of internal coordinate values

Definition at line 34 of file extrema/Internals/internals.h.

Referenced by mem_alloc(), and psi::extrema::zmat::zmat().

double * psi::extrema::internals::fcoords_old [protected]

full set of previous internal coordinate values

Definition at line 34 of file extrema/Internals/internals.h.

Referenced by mem_alloc(), and psi::extrema::zmat::zmat().

double * psi::extrema::internals::fgrads [protected]

full set of internal coordinate gradients

Definition at line 34 of file extrema/Internals/internals.h.

Referenced by psi::extrema::zmat::grad_trans(), grad_trans(), mem_alloc(), and psi::extrema::deloc::optimize().

double ** psi::extrema::internals::B [protected]

the B matrix

Definition at line 34 of file extrema/Internals/internals.h.

Referenced by compute_A(), psi::extrema::zmat::compute_B(), compute_G(), psi::extrema::deloc::deloc(), grad_trans(), mem_alloc(), print_B(), and ~internals().

double ** psi::extrema::internals::B_red [protected]

reduced dimension B matrix (no equivalent coords)

Definition at line 34 of file extrema/Internals/internals.h.

Referenced by psi::extrema::zmat::compute_B(), mem_alloc(), and ~internals().

double ** psi::extrema::internals::G [protected]

the G=B.B^t matrix

Definition at line 34 of file extrema/Internals/internals.h.

Referenced by compute_A(), compute_G(), psi::extrema::deloc::deloc(), mem_alloc(), and ~internals().

double ** psi::extrema::internals::A [protected]

the A=u.B^t.G matrix

Definition at line 34 of file extrema/Internals/internals.h.

Referenced by back_transform(), compute_A(), mem_alloc(), and ~internals().

double ** psi::extrema::internals::u [protected]

the u matrix (triplets of inverse atomic masses)

Definition at line 34 of file extrema/Internals/internals.h.

Referenced by compute_A(), psi::extrema::zmat::compute_B(), compute_G(), psi::extrema::deloc::deloc(), mem_alloc(), and ~internals().


The documentation for this class was generated from the following files:
Generated on Wed Feb 13 16:36:15 2008 for PSI by  doxygen 1.5.4