#include "../include/dprec.fh" module pupildata use memory_module, only: ixStack=>ix, realStack=>x, ihStack=>ih logical pupactive integer puperror,iPup,jPup,bs1,bs2,iresPup,l_puptmp integer pupLevelData ! Data to pass throught PUPIL Interface integer pupStep integer pupQZchange ! is != 0 if quantum zone has changed integer pupqatoms ! total number of PUPIL quantum atoms integer pupnumnb14 ! total number of initial nonboded 14 pairs integer pupnbonh ! total number of initial boded H pairs integer pupnbona ! total number of initial boded pairs integer pupntheth ! total number of initial angled H tern. integer pupntheta ! total number of initial angled tern. integer pupnphih ! total number of initial dihed. H quat. integer pupnphia ! total number of initial dihed. quat. integer, allocatable :: pupnb14 (:) ! initial nonbonded 14 pair list 3,x since it contains ii,jj,ic0 integer, allocatable :: pupbonh (:) ! initial boded H pairs integer, allocatable :: pupbona (:) ! initial boded pairs integer, allocatable :: puptheth (:) ! initial angled H tern. integer, allocatable :: puptheta (:) ! initial angled tern. integer, allocatable :: pupphih (:) ! initial dihed. H quat. integer, allocatable :: pupphia (:) ! initial dihed. quat. integer, allocatable :: pupmask (:) ! mask: 0 -> classic 1-> quantum integer, allocatable :: pupqlist (:) ! atom number for quantum list integer, allocatable :: pupatm (:) integer, allocatable :: pupres (:) ! residue pointers !integer, allocatable :: ixStack (:) ! Captioning the ix pointer _REAL_, allocatable :: pupchg (:) ! Initial charges over each atoms in the system _REAL_, allocatable :: qcdata (:) _REAL_, allocatable :: qfpup (:) _REAL_, allocatable :: qcell (:) !_REAL_, allocatable :: realStack(:) ! Captioning the x pointer !character(len=4), dimension(:), allocatable :: ihStack ! Captioning the ih pointer character(len=10), dimension(:),allocatable :: keyres ! to keep the MM key RESIDUE_NAME character(len=20),dimension(:), allocatable :: keyMM ! to keep the MM key particles RESIDUE_NAME.ATOM_NAME character(len=20) strAux _REAL_ qmEnergy end module pupildata ! *********************************************************************** ! *********************************************************************** subroutine get_atomic_number_pupil(name,rmass,atomic_number) ! ! Joan Torras (February 2014) ! This subroutine assigns atomic number based upon the first and second ! letter of the atom symbol and checking with the atomic mass ! ! name : characters containing the atomic name ! rmass : atomic mass which has to fit with the given atom name ! atomic_number: atomic number which has to fit with the given atom name ! implicit none character(len=4), intent(in) :: name double precision, intent(in) :: rmass integer, intent(out) :: atomic_number logical :: check_wmass integer :: j, k, coincidence character(len=2), dimension(109):: elemnt=(/ & 'H ','HE', & 'LI','BE','B ','C ','N ','O ','F ','NE', & 'NA','MG','AL','SI','P ','S ','CL','AR', & 'K ','CA', & 'SC','TI','V ','CR','MN','FE','CO','NI','CU','ZN', & 'GA','GE','AS','SE','BR','KR', & 'RB','SR', & 'Y ','ZR','NB','MO','TC','RU','RH','PD','AG','CD', & 'IN','SN','SB','TE','I ','XE', & 'CS','BA', & 'LA','CE','PR','ND','PM','SM','EU', & 'GD','TB','DY','HO','ER','TM','YB', & 'LU','HF','TA','W ','RE','OS','IR','PT','AU','HG', & 'TL','PB','BI','PO','AT','RN', & 'FR','RA', & 'AC','TH','PA','U ','NP','PU','AM', & 'CM','BK','CF','ES','FM','MD','NO', & 'LR','RF','DB','SG','BH','HS','MT'/) ! ! compare values in name to those in elemnt and assign atomic ! numbers accordingly ! j = 1 coincidence = 0 do while ( (j .le. 109) .and. (coincidence .eq. 0) ) if(name(1:1) .eq. elemnt(j)(1:1)) then if(elemnt(j)(2:2) .eq. " ") then if(name(2:2) .eq. " ") then atomic_number = j coincidence = coincidence + 1 else do k=j+1,109 if((name(1:1) .eq. elemnt(k)(1:1)) .and. & (name(2:2) .eq. elemnt(k)(2:2)))then atomic_number = k coincidence = coincidence + 1 endif enddo if (coincidence .eq. 0) then atomic_number = j coincidence = coincidence + 1 endif endif else if(name(2:2) .eq. elemnt(j)(2:2)) then atomic_number = j coincidence = coincidence + 1 endif endif end if j = j + 1 enddo ! Checking coincidence between atomic_number and atomic mass if (.not. check_wmass(atomic_number,rmass)) then coincidence = 0 endif ! If not found so far lets shearch by atomic mass if(coincidence .eq. 0) then j = 1 do while ( (j .le. 109) .and. (coincidence .eq. 0) ) if (check_wmass(j,rmass)) then atomic_number = j coincidence = 1 endif j = j + 1 enddo endif if(coincidence .eq. 0) then write(6,*) 'PUPIL: Unable to correctly identify element ', name call mexit(6,1) else !write(6,*) 'FOUND: ',elemnt(atomic_number),' related to ',name endif return end subroutine get_atomic_number_pupil ! *********************************************************************** ! *********************************************************************** logical function check_wmass (atomic_number, rmass) ! ! Joan Torras (February 2014) ! ! This subroutine is checking if the pair made of atomic_number and rmass have ! coincidence with the values from periodic table ! ! return a logical value ==> Are the pair of given values: mass, and atomic ! number, in coincidence with the periodic table values with a maximum difference ! of 0.5 a.u. of mass ? implicit none integer, intent(in) :: atomic_number double precision, intent(in) :: rmass integer :: j, k, coincidence double precision, dimension(109):: wmass=(/ & 1.0079, 4.0026,6.941 ,9.0122 ,10.811 , 12.0107,14.0067 , 15.9994, 18.9984, 20.1797, & 22.9897, 24.305 ,26.9815 ,28.0855 ,30.9738 , 32.065 ,35.453 , 39.948 , 39.0983, 40.078 , & 44.9559, 47.867 ,50.9415 ,51.9961 ,54.938 , 55.845 ,58.9332 , 58.6934, 63.546 , 65.39 , & 69.723 , 72.64 ,74.9216 ,78.96 ,79.904 , 83.8 ,85.4678 , 87.62 , 88.9059, 91.224 , & 92.9064, 95.94 ,98.0 ,101.07 ,102.9055,106.42 ,107.8682,112.411 ,114.818 ,118.71 , & 121.76 ,127.6 ,126.9045,131.293 ,132.9055,137.327 ,138.9055,140.116 ,140.9077,144.24 , & 145.0 ,150.36 ,151.964 ,157.25 ,158.9253,162.5 ,164.9303,167.259 ,168.9342,173.04 , & 174.967 ,178.49 ,180.9479,183.84 ,186.207 ,190.23 ,192.217 ,195.078 ,196.9665,200.59 , & 204.3833,207.2 ,208.9804,209.0 ,210.0 ,222.0 ,223.0 ,226.0 ,227.0 ,232.0381, & 231.0359,238.0289,237.0 ,244.0 ,243.0 ,247.0 ,247.0 ,251.0 ,252.0 ,257.0 , & 258.0 ,259.0 ,262.0 ,261.0 ,262.0 ,266.0 ,264.0 ,277.0 ,268.0/) check_wmass = abs(rmass-wmass(atomic_number)) .le. 0.5 !write (6,*) check_wmass, atomic_number, rmass, wmass(atomic_number) end function check_wmass ! *********************************************************************** ! *********************************************************************** subroutine deleting_qm_atoms() use pupildata, x=>realStack , ix=>ixStack, ih=>ihStack use parms, only: req implicit none #include "../include/memory.h" #include "extra_pts.h" #include "../include/md.h" ! integer i ! Initializing the list of structures for a new QM zone call init_extra_pts( & ix(iibh),ix(ijbh),ix(iicbh), & ix(iiba),ix(ijba),ix(iicba), & ix(i24),ix(i26),ix(i28),ix(i30), & ix(i32),ix(i34),ix(i36),ix(i38), & ix(i40),ix(i42),ix(i44),ix(i46),ix(i48), & ix(i50),ix(i52),ix(i54),ix(i56),ix(i58), & ih(m06),ix,x,ix(i08),ix(i10), & nspm,ix(i70),x(l75),tmass,tmassinv,x(lmass),x(lwinv),req) if(nbonh.gt.0) call setbon(nbonh,ix(iibh),ix(ijbh),ix(iicbh), & ix(ibellygp)) ! remove bonds between QM atoms from list if(nbona.gt.0) call setbon(nbona,ix(iiba),ix(ijba),ix(iicba), & ix(ibellygp)) ! remove bonds between QM atoms from list if(ntheth.gt.0) call setang(ntheth,ix(i24),ix(i26),ix(i28),ix(i30), & ix(ibellygp)) ! remove angles between QM atoms from list if(ntheta.gt.0) call setang(ntheta,ix(i32),ix(i34),ix(i36),ix(i38),& ix(ibellygp)) ! remove angles between QM atoms from list if(nphih.gt.0) call setdih(nphih,ix(i40),ix(i42),ix(i44),ix(i46), & ix(i48), ix(ibellygp)) ! remove dihedrals between QM atoms from list if(nphia.gt.0) call setdih(nphia,ix(i50),ix(i52),ix(i54),ix(i56), & ix(i58), ix(ibellygp)) ! remove dihedrals between QM atoms from list !if(pupQZChange .ne. 0) then ! write(6,*) 'NUMBER OF H BONDS = ',nbonh ! write(6,*) 'NUMBER OF BONDS = ',nbona ! write(6,*) 'NUMBER OF H ANGLES = ',ntheth ! write(6,*) 'NUMBER OF ANGLES = ',ntheta ! write(6,*) 'NUMBER OF H DIHEDRALS = ',nphih ! write(6,*) 'NUMBER OF DIHEDRALS = ',nphia ! write(6,*) 'NUMBER OF NONBONDED 14 PAIRS = ',numnb14 ! do i=1,numnb14 ! write(6,*) ix(inb_14+(i-1)*3),ix(inb_14+(i-1)*3+1) ! enddo !endif ! write(6,*) 'nbonh',nbonh,'iibh',iibh,'nbona',nbona,'ntheth',ntheth ! write(6,*) 'ntheta',ntheta,'nphih',nphih,'nphia',nphia ! write(6,*) 'i26',i26,'i56',i56 ! write(6,*) 'IGB',igb return end subroutine deleting_qm_atoms ! *********************************************************************** ! *********************************************************************** subroutine write_energies(str,e) implicit none _REAL_ e(25) character(len=4) str integer i,k do i=1,5 write(6,"(a4,2x,'ENERG',5(2x,i2,2x,d15.8))") & str,( (i-1)*5+k,e((i-1)*5+k), k=1,5 ) enddo write(6,*) " ----------------------------------" end subroutine write_energies !***********************************************************************