hdiff output

r11830/sff2.c 2017-01-21 10:40:24.828348695 +0000 r11829/sff2.c 2017-01-21 10:40:25.096362741 +0000
1278:      *            e = pk(ic) * (1.0+phase*cos(pn(ic)*phi)1278:      *            e = pk(ic) * (1.0+phase*cos(pn(ic)*phi)
1279:      *            where phase = 1.0 or -1.0, and pn = 1,2,3,4, or 61279:      *            where phase = 1.0 or -1.0, and pn = 1,2,3,4, or 6
1280: 1280: 
1281:      *     ----- energy terms for dihedrals are expressed in terms of1281:      *     ----- energy terms for dihedrals are expressed in terms of
1282:      *           cosphi in order to eliminate problems for planar angles ----1282:      *           cosphi in order to eliminate problems for planar angles ----
1283:      */1283:      */
1284: 1284: 
1285:     /*  sf344> try to modify the potential to have a different form, 1285:     /*  sf344> try to modify the potential to have a different form, 
1286:                consistent with Amber:1286:                consistent with Amber:
1287:                   e = pk(ic) * (1.0 + cos(phase)*cos(pn(ic)*phi) + sin(phase)*sin(pn(ic)*phi)1287:                   e = pk(ic) * (1.0 + cos(phase)*cos(pn(ic)*phi) + sin(phase)*sin(pn(ic)*phi)
1288:       sf344, 06-10-2008:  
1289:         the (only) problem with this form is that when abs(cos(phi))=1, we are dividing by zero in the  
1290:         derivative and second derivative expression for the term with sin(phi). Small angle  
1291:         approximation solves the problem (see below). It remains to be seen what kind of  
1292:         discontinuities will this cause for parameters that have sin(phase)>0 (most parametrisations 
1293:         have it explicitly set to zero, I guess because of this problem). 
1294:  
1295: 1288: 
1296:      */1289:      */
1297: 1290: 
1298: 1291: 
1299:     ct2 = ct*ct;1292:     ct2 = ct*ct;
1300:     st2 = 1.0-ct2;1293:     st2 = 1.0-ct2;
1301:     st = sqrt(st2)*s/fabs(s) ;1294:     st = sqrt(st2)*s/fabs(s) ;
1302:  
1303:     if (st >= 0 && st < 1.e-3) 
1304:       st = 1.e-3; 
1305:     else if (st < 0 && st > -1.e-3) 
1306:       st = -1.e-3; 
1307:     if (st2 >= 0 && st2 < 1.e-3) 
1308:       st2 = 1.e-3; 
1309:     else if (st2 < 0 && st2 > -1.e-3) 
1310:       st2 = -1.e-3; 
1311:  
1312:   multi_term:1295:   multi_term:
1313:     iper = (int) fabs(Pn[atyp]) + 0.0001;1296:     iper = (int) fabs(Pn[atyp]) + 0.0001;
1314:     assert(iper != (fabs(Pn[atyp]) - 0.0001));1297:     assert(iper != (fabs(Pn[atyp]) - 0.0001));
1315:     assert(iper >= 1 && iper <= 6);1298:     assert(iper >= 1 && iper <= 6);
1316:     switch (iper) {1299:     switch (iper) {
1317: 1300: 
1318:     case 1:1301:     case 1:
1319:       e = ct;1302:       e = ct;
1320:       df = 1.0;1303:       df = 1.0;
1321:       ddf = 0.0;1304:       ddf = 0.0;


legend
Lines Added 
Lines changed
 Lines Removed

hdiff - version: 2.1.0