Difference for src/GridForceGrid.C from version 1.2 to 1.3

version 1.2version 1.3
Line 13
Line 13
 #include "common.h" #include "common.h"
  
 #define MIN_DEBUG_LEVEL 3 #define MIN_DEBUG_LEVEL 3
 //#define DEBUGM #define DEBUGM
 #include "Debug.h" #include "Debug.h"
  
  
Line 196
Line 196
      for (int i1 = 0; i1 < 3; i1++) {      for (int i1 = 0; i1 < 3; i1++) {
  if (cross(Avec[i0], Kvec[i1]) == 0) {  if (cross(Avec[i0], Kvec[i1]) == 0) {
      cont[i1] = TRUE;      cont[i1] = TRUE;
      offset[i1] = simParams->gridforceVOffset[i0];      offset[i1] = simParams->gridforceVOffset[i0] * factor;
      gap[i1] = (inv * (Avec[i0] - Kvec[i1])).length(); // want in grid-point units (normal = 1)      gap[i1] = (inv * (Avec[i0] - Kvec[i1])).length(); // want in grid-point units (normal = 1)
            
      if (gap[i1] < 0) {      if (gap[i1] < 0) {
Line 277
Line 277
   
  if (cont[i0] && fabs(offset[i0] - (p_avg[i0]-n_avg[i0])) > modThresh) {  if (cont[i0] && fabs(offset[i0] - (p_avg[i0]-n_avg[i0])) > modThresh) {
      iout << iWARN << "GRID FORCE POTENTIAL DIFFERENCE IN K" << i0      iout << iWARN << "GRID FORCE POTENTIAL DIFFERENCE IN K" << i0
   << " DIRECTION IS GREATER THAN " << modThresh << " KCAL/MOL*E\n" << endi;   << " DIRECTION IS " << offset[i0] - (p_avg[i0]-n_avg[i0]) << " KCAL/MOL*E\n" << endi;
  }  }
  DebugM(3, "n_avg[" << i0 << "] = " << n_avg[i0] << "\n" << endi);  DebugM(3, "n_avg[" << i0 << "] = " << n_avg[i0] << "\n" << endi);
  DebugM(3, "p_avg[" << i0 << "] = " << p_avg[i0] << "\n" << endi);  DebugM(3, "p_avg[" << i0 << "] = " << p_avg[i0] << "\n" << endi);
Line 378
Line 378
  int inds2[3];  int inds2[3];
  int ind2 = 0;  int ind2 = 0;
   
   float voff = 0.0;
  int bit = 1; // bit = 2^i1 in the below loop  int bit = 1; // bit = 2^i1 in the below loop
  for (int i1 = 0; i1 < 3; i1++) {  for (int i1 = 0; i1 < 3; i1++) {
      inds2[i1] = (inds[i1] + ((i0 & bit) ? 1 : 0)) % k[i1];      inds2[i1] = (inds[i1] + ((i0 & bit) ? 1 : 0)) % k[i1];
      ind2 += inds2[i1] * dk[i1];      ind2 += inds2[i1] * dk[i1];
       
       // Deal with voltage offsets
       if (cont[i1] && inds[i1] == (k[i1]-1) && inds2[i1] == 0) {
   voff += offset[i1];
   DebugM(3, "offset[" << i1 << "] = " << offset[i1] << "\n" << endi);
       }
       
      bit <<= 1; // i.e. multiply by 2      bit <<= 1; // i.e. multiply by 2
  }  }
   
Line 402
Line 410
  Bool zero_derivs = FALSE;  Bool zero_derivs = FALSE;
  int dk_hi[3];  int dk_hi[3];
  int dk_lo[3];  int dk_lo[3];
   float voffs[3];
  for (int i1 = 0; i1 < 3; i1++) {  for (int i1 = 0; i1 < 3; i1++) {
      if (inds2[i1] == 0) {      if (inds2[i1] == 0) {
  if (cont[i1]) {  if (cont[i1]) {
      dk_hi[i1] = dk[i1];      dk_hi[i1] = dk[i1];
      dk_lo[i1] = -(k[i1]-1) * dk[i1];      dk_lo[i1] = -(k[i1]-1) * dk[i1];
       voffs[i1] = offset[i1];
  }  }
  else zero_derivs = TRUE;  else zero_derivs = TRUE;
      }      }
Line 414
Line 424
  if (cont[i1]) {  if (cont[i1]) {
      dk_hi[i1] = -(k[i1]-1) * dk[i1];      dk_hi[i1] = -(k[i1]-1) * dk[i1];
      dk_lo[i1] = dk[i1];      dk_lo[i1] = dk[i1];
       voffs[i1] = offset[i1];
  }  }
  else zero_derivs = TRUE;  else zero_derivs = TRUE;
      }      }
      else {      else {
  dk_hi[i1] = dk[i1];  dk_hi[i1] = dk[i1];
  dk_lo[i1] = dk[i1];  dk_lo[i1] = dk[i1];
   voffs[i1] = 0.0;
      }      }
  }  }
   
Line 433
Line 445
  box->loc = dg;  box->loc = dg;
   
  // V  // V
  box->b[i0] = grid[ind2];  box->b[i0] = grid[ind2] + voff;
   
  // First derivatives  if (zero_derivs) {
  //  dV/dx      box->b[8+i0] = 0.0;
  box->b[8+i0] = (zero_derivs) ? 0.0 :      box->b[16+i0] = 0.0;
      0.5 * (grid[ind2 + dk_hi[0]] - grid[ind2 - dk_lo[0]]);      box->b[24+i0] = 0.0;
  //  dV/dy      box->b[32+i0] = 0.0;
  box->b[16+i0] = (zero_derivs) ? 0.0 :      box->b[40+i0] = 0.0;
      0.5 * (grid[ind2 + dk_hi[1]] - grid[ind2 - dk_lo[1]]);      box->b[48+i0] = 0.0;
  //  dV/dz      box->b[56+i0] = 0.0;
  box->b[24+i0] = (zero_derivs) ? 0.0 :  } else {
      0.5 * (grid[ind2 + dk_hi[2]] - grid[ind2 - dk_lo[2]]);      box->b[8+i0] = 0.5 * (grid[ind2 + dk_hi[0]] - grid[ind2 - dk_lo[0]] + voffs[0]); //  dV/dx
        box->b[16+i0] = 0.5 * (grid[ind2 + dk_hi[1]] - grid[ind2 - dk_lo[1]] + voffs[1]); //  dV/dy
  // Mixed second derivatives      box->b[24+i0] = 0.5 * (grid[ind2 + dk_hi[2]] - grid[ind2 - dk_lo[2]] + voffs[2]); //  dV/dz
  //  d2V/dxdy      box->b[32+i0] = 0.25 * (grid[ind2 + dk_hi[0] + dk_hi[1]] - grid[ind2 - dk_lo[0] + dk_hi[1]]
  box->b[32+i0] = (zero_derivs) ? 0.0 :      - grid[ind2 + dk_hi[0] - dk_lo[1]] + grid[ind2 - dk_lo[0] - dk_lo[1]]); //  d2V/dxdy
      0.25 * (grid[ind2 + dk_hi[0] + dk_hi[1]] - grid[ind2 - dk_lo[0] + dk_hi[1]]      box->b[40+i0] = 0.25 * (grid[ind2 + dk_hi[0] + dk_hi[2]] - grid[ind2 - dk_lo[0] + dk_hi[2]]
      - grid[ind2 + dk_hi[0] - dk_lo[1]] + grid[ind2 - dk_lo[0] - dk_lo[1]]);      - grid[ind2 + dk_hi[0] - dk_lo[2]] + grid[ind2 - dk_lo[0] - dk_lo[2]]); //  d2V/dxdz
  //  d2V/dxdz      box->b[48+i0] = 0.25 * (grid[ind2 + dk_hi[1] + dk_hi[2]] - grid[ind2 - dk_lo[1] + dk_hi[2]]
  box->b[40+i0] = (zero_derivs) ? 0.0 :      - grid[ind2 + dk_hi[1] - dk_lo[2]] + grid[ind2 - dk_lo[1] - dk_lo[2]]); //  d2V/dydz
      0.25 * (grid[ind2 + dk_hi[0] + dk_hi[2]] - grid[ind2 - dk_lo[0] + dk_hi[2]] 
      - grid[ind2 + dk_hi[0] - dk_lo[2]] + grid[ind2 - dk_lo[0] - dk_lo[2]]); 
  //  d2V/dydz 
  box->b[48+i0] = (zero_derivs) ? 0.0 : 
      0.25 * (grid[ind2 + dk_hi[1] + dk_hi[2]] - grid[ind2 - dk_lo[1] + dk_hi[2]] 
      - grid[ind2 + dk_hi[1] - dk_lo[2]] + grid[ind2 - dk_lo[1] - dk_lo[2]]); 
   
  // d3V/dxdydz      box->b[56+i0] = // d3V/dxdydz
  box->b[56+i0] = (zero_derivs) ? 0.0 : 
      0.125 * (grid[ind2 + dk_hi[0] + dk_hi[1] + dk_hi[2]] - grid[ind2 + dk_hi[0] + dk_hi[1] - dk_lo[2]]      0.125 * (grid[ind2 + dk_hi[0] + dk_hi[1] + dk_hi[2]] - grid[ind2 + dk_hi[0] + dk_hi[1] - dk_lo[2]]
       - grid[ind2 + dk_hi[0] - dk_lo[1] + dk_hi[2]] - grid[ind2 - dk_lo[0] + dk_hi[1] + dk_hi[2]]       - grid[ind2 + dk_hi[0] - dk_lo[1] + dk_hi[2]] - grid[ind2 - dk_lo[0] + dk_hi[1] + dk_hi[2]]
       + grid[ind2 + dk_hi[0] - dk_lo[1] - dk_lo[2]] + grid[ind2 - dk_lo[0] + dk_hi[1] - dk_lo[2]]       + grid[ind2 + dk_hi[0] - dk_lo[1] - dk_lo[2]] + grid[ind2 - dk_lo[0] + dk_hi[1] - dk_lo[2]]
       + grid[ind2 - dk_lo[0] - dk_lo[1] + dk_hi[2]] - grid[ind2 - dk_lo[0] - dk_lo[1] - dk_lo[2]]);       + grid[ind2 - dk_lo[0] - dk_lo[1] + dk_hi[2]] - grid[ind2 - dk_lo[0] - dk_lo[1] - dk_lo[2]]);
   }
   
  DebugM(2, "V = " << box->b[i0] << "\n");  DebugM(2, "V = " << box->b[i0] << "\n");
   

Legend:
Removed in v.1.2 
changed lines
 Added in v.1.3


Made by using version 1.53 of cvs2html