| /*************************************************************************/ |
| /* */ |
| /* Copyright (c) 1994 Stanford University */ |
| /* */ |
| /* All rights reserved. */ |
| /* */ |
| /* Permission is given to use, copy, and modify this software for any */ |
| /* non-commercial purpose as long as this copyright notice is not */ |
| /* removed. All other uses, including redistribution in whole or in */ |
| /* part, are forbidden without prior written permission. */ |
| /* */ |
| /* This software is provided with absolutely no warranty and no */ |
| /* support. */ |
| /* */ |
| /*************************************************************************/ |
| |
| EXTERN_ENV |
| #include "mdvar.h" |
| #include "parameters.h" |
| #include "mddata.h" |
| #include "split.h" |
| #include "global.h" |
| |
| /* predicts new values for displacement and its five derivatives |
| * |
| * NOR1 : NOR1 = NORDER + 1 = 7 (for a sixth-order method) |
| */ |
| void PREDIC(double *C, long NOR1, long ProcID) |
| { |
| /* this routine calculates predicted F(X), F'(X), F''(X), ... */ |
| |
| long JIZ; |
| long JI; |
| long L; |
| double S; |
| long func, mol, dir, atom; |
| |
| JIZ=2; |
| /* .....loop over F(X), F'(X), F''(X), ..... */ |
| for (func = 0; func < NORDER; func++) { |
| for (mol = StartMol[ProcID]; mol < StartMol[ProcID+1]; mol++) |
| for ( dir = 0; dir < NDIR; dir++) |
| for ( atom = 0; atom < NATOM; atom++ ) { |
| JI = JIZ; |
| /* sum over Taylor Series */ |
| S = 0.0; |
| for ( L = func; L < NORDER; L++) { |
| S += C[JI] * VAR[mol].F[L+1][dir][atom]; |
| JI++; |
| } /* for */ |
| VAR[mol].F[func][dir][atom] += S; |
| } /* for atom */ |
| JIZ += NOR1; |
| } /* for func */ |
| } /* end of subroutine PREDIC */ |
| |
| /* corrects the predicted values, based on forces etc. computed in the interim |
| * |
| * PCC : the predictor-corrector constants |
| * NOR1: NORDER + 1 = 7 for a sixth-order method) |
| */ |
| |
| void CORREC(double *PCC, long NOR1, long ProcID) |
| { |
| /* |
| .....this routine calculates corrected F(X), F'(X), F"(X), .... |
| from corrected F(X) = predicted F(X) + PCC(1)*(FR-SD) |
| where SD is predicted accl. F"(X) and FR is computed |
| accl. (force/mass) at predicted position |
| */ |
| |
| double Y; |
| long mol, dir, atom, func; |
| |
| for (mol = StartMol[ProcID]; mol < StartMol[ProcID+1]; mol++) { |
| for (dir = 0; dir < NDIR; dir++) { |
| for (atom = 0; atom < NATOM; atom++) { |
| Y = VAR[mol].F[FORCES][dir][atom] - VAR[mol].F[ACC][dir][atom]; |
| for ( func = 0; func < NOR1; func++) |
| VAR[mol].F[func][dir][atom] += PCC[func] * Y; |
| } /* for atom */ |
| } /* for dir */ |
| } /* for mol */ |
| |
| } /* end of subroutine CORREC */ |