[7137e5d] | 1 | real*8 function enyshe_simple()
|
---|
| 2 |
|
---|
| 3 | include "INCL.H"
|
---|
| 4 |
|
---|
| 5 | integer i,j,numberOfAtoms,indexOfMolecule, i14
|
---|
| 6 | integer indexOfFirstResidue, indexOfLastResidue
|
---|
| 7 | integer indexOfFirstAtom, indexOfLastAtom, atomTypeI, indexOfAtom
|
---|
| 8 | integer typeOfDihedral, classOfPotential
|
---|
| 9 | integer indexOfFirstDihedral, indexOfLastDihedral
|
---|
| 10 | real*8 totalEnergy, coulombEnergy, rij2, rij4, rij6, rij12, E0
|
---|
| 11 | real*8 vdWEnergy, vdWDiff, torsionEnergy,xatI, yatI, zatI, cgI, tmpE
|
---|
| 12 | real*8 hydrogenBondEnergy
|
---|
| 13 |
|
---|
| 14 | indexOfMolecule = 1
|
---|
| 15 | indexOfFirstResidue = irsml1(indexOfMolecule)
|
---|
| 16 | indexOfLastResidue = irsml2(indexOfMolecule)
|
---|
| 17 | indexOfFirstAtom = iatrs1(indexOfFirstResidue)
|
---|
| 18 | indexOfLastAtom = iatrs2(indexOfLastResidue)
|
---|
| 19 | numberOfAtoms = indexOfLastAtom - indexOfFirstAtom
|
---|
| 20 | !$acc region
|
---|
| 21 | do i = indexOfFirstAtom, indexOfLastAtom
|
---|
| 22 | xatI = xat(i)
|
---|
| 23 | yatI = yat(i)
|
---|
| 24 | zatI = zat(i)
|
---|
| 25 | cgI = cgat(i)
|
---|
| 26 | atomTypeI = ityat(i)
|
---|
| 27 |
|
---|
| 28 | do j = i + 1, indexOfLastAtom
|
---|
| 29 | rij2 = (xatI - xat(j)) * (xatI - xat(j)) + (yatI - yat(j)) * (yatI - yat(j)) &
|
---|
| 30 | + (zatI - zat(j)) * (zatI - zat(j))
|
---|
| 31 | rij4 = rij2 * rij2
|
---|
| 32 | rij6 = rij4 * rij2
|
---|
| 33 |
|
---|
| 34 | ! Calculate Coulomb energy
|
---|
| 35 | coulombEnergy = coulombEnergy + cgI * cgat(j) / sqrt(rij2)
|
---|
| 36 | ! Calculate van-der-Waals and hydrogen bond energy
|
---|
| 37 | vdWEnergy = vdWEnergy + aij(atomTypeI, ityat(j)) / (rij6 * rij6) - cij(atomTypeI, ityat(j)) / rij6
|
---|
| 38 | ! Calculate hydrogen bond energy
|
---|
| 39 | hydrogenBondEnergy = hydrogenBondEnergy + ahb(ityat(i),ityat(j))/(rij6*rij6) &
|
---|
| 40 | - chb(ityat(i),ityat(j))/(rij6*rij4)
|
---|
| 41 | end do
|
---|
| 42 |
|
---|
| 43 | do i14=i14at1(i),i14at2(i)
|
---|
| 44 | j = l14at(i14)
|
---|
| 45 | rij2 = (xat(i) - xat(j)) * (xat(i) - xat(j)) + (yat(i) - yat(j)) * (yat(i) - yat(j)) &
|
---|
| 46 | + (zat(i) - zat(j)) * (zat(i) - zat(j))
|
---|
| 47 | rij4 = rij2 * rij2
|
---|
| 48 | rij12 = rij4 * rij4 * rij4
|
---|
| 49 |
|
---|
| 50 | vdWDiff = vdWDiff - 0.5 * aij(ityat(i), ityat(j)) / (rij12)
|
---|
| 51 | end do
|
---|
| 52 |
|
---|
| 53 | end do
|
---|
| 54 |
|
---|
| 55 | indexOfFirstDihedral = ivrml1(indexOfMolecule)
|
---|
| 56 | indexOfLastDihedral = indexOfFirstDihedral + nvrml(indexOfMolecule)
|
---|
| 57 |
|
---|
| 58 | do i = indexOfFirstDihedral, indexOfLastDihedral
|
---|
| 59 | indexOfAtom = iatvr(i)
|
---|
| 60 | typeOfDihedral = ityvr(i)
|
---|
| 61 | classOfPotential = iclvr(i)
|
---|
| 62 | if (typeOfDihedral == 3) then
|
---|
| 63 | if (e0to(classOfPotential) == 0.0) then
|
---|
| 64 | tmpE = torsionEnergy + e0to(classOfPotential) * (1.0 + sgto(classOfPotential) &
|
---|
| 65 | * cos(toat(indexOfAtom) * rnto(classOfPotential)))
|
---|
| 66 | endif
|
---|
| 67 | endif
|
---|
| 68 |
|
---|
| 69 | torsionEnergy = tmpE
|
---|
| 70 | end do
|
---|
| 71 |
|
---|
| 72 | totalEnergy = coulombEnergy + vdWEnergy + vdWDiff + torsionEnergy
|
---|
| 73 |
|
---|
| 74 | enyshe_simple = totalEnergy
|
---|
| 75 | !$acc end region
|
---|
| 76 | return
|
---|
| 77 |
|
---|
| 78 | end function enyshe_simple
|
---|