source: init_molecule.f@ bd2278d

Last change on this file since bd2278d was bd2278d, checked in by baerbaer <baerbaer@…>, 16 years ago

Reformatting comments and continuation marks.

Fortran 90 and higher use ! to mark comments no matter where they are in the
code. The only valid continuation marker is &.
I also added the SMMP.kdevelop.filelist to the repository to make it easier
to use kdevelop.

git-svn-id: svn+ssh://svn.berlios.de/svnroot/repos/smmp/trunk@12 26dc1dd8-5c4e-0410-9ffe-d298b4865968

  • Property mode set to 100644
File size: 6.3 KB
Line 
1! **************************************************************
2!
3! This file contains the subroutines: init_molecule
4!
5! Copyright 2003-2005 Frank Eisenmenger, U.H.E. Hansmann,
6! Shura Hayryan, Chin-Ku
7! Copyright 2007 Frank Eisenmenger, U.H.E. Hansmann,
8! Jan H. Meinke, Sandipan Mohanty
9!
10! **************************************************************
11! FIXME: Data in varfile determines which molecule is changed.
12
13 subroutine init_molecule(iabin,grpn,grpc,seqfile,varfile)
14
15! ----------------------------------------------------------
16! PURPOSE: construct starting structure of molecule(s)
17!
18! iabin = 1 : ab Initio using sequence &
19! variables given in input files
20! iabin != 1 : sequence, variable information
21! from PDB-file
22!
23! grpn: N-terminal group
24! grpc: C-terminal group
25!
26! CALLS: addend,bldmol,c_alfa,getmol,iendst, mklist, nursvr,
27! pdbread,pdbvars,redseq,redvar,setmvs
28!
29! ----------------------------------------------------------
30
31 include 'INCL.H'
32 include 'INCP.H'
33
34!f2py character*80 optional, intent(in) :: seqfile = ' '
35!f2py character*80 optional, intent(in) :: varfile = ' '
36
37 character grpn*4,grpc*4
38 character navr*3, nars*4
39 character seqfile*80, varfile*80
40 integer ontlml
41 logical readFromStdin
42
43 ontlml = 1
44 readFromStdin = .false.
45
46 write (*,*) 'init_molecule: Solvent: ', itysol
47 if (iabin.eq.1) then
48
49! ----------------------------------------- get sequence for molecule(s)
50 lunseq=11
51 if (ntlml.gt.0) then
52 ontlml = ntlml + 1
53 endif
54 if (iendst(seqfile).le.1.or.seqfile.eq.' ') then
55 1 write (*,'(/,a,$)') ' file with SEQUENCE:'
56 seqfil=' '
57 read (*,'(a)',err=1) seqfil
58 readFromStdin = .true.
59 else
60 seqfil = seqfile
61 endif
62 call redseq
63
64 write (*,*) 'File with sequence is ', seqfil(1:iendst(seqfil))
65
66! --------------------------------- read & assemble data from libraries
67! initial coordinates, interaction lists
68
69 ntl = ntlml
70 do i=ontlml, ntl
71
72 call getmol(i) ! assemble data from libraries
73
74 do j=1,6 ! initialize global parameters
75 gbpr(j,i)=0.d0
76 enddo
77
78 call bldmol(i) ! co-ordinates
79
80 ntlml = i
81 call addend(i,grpn,grpc) ! modify ends
82 call setmvs(i) ! determine sets of moving atoms for given variables
83 call mklist(i) ! compile lists of interaction partners
84
85 enddo
86
87! --------------------------- Read the initial conformation if necessary
88 if(readFromStdin) then
89 write (*,'(a,$)') ' file with VARIABLES:'
90!
91 varfil=' '
92 read(*,'(a)',end=2,err=2) varfil
93 else
94 varfil = varfile
95 endif
96 l=iendst(varfil)
97 if (l.gt.0.and.varfil.ne.' ') then
98 write (*,'(1x,a,/)') varfil(1:l)
99 lunvar=13
100
101 call redvar ! get vars. and rebuild
102
103 endif
104
105 2 write(*,*) ' '
106
107! -------------------- get: nvr,idvr, vlvr, olvlvr
108 nvr = 0
109 do i=1,ivrml1(ntlml)+nvrml(ntlml)-1
110
111 if (.not.fxvr(i)) then
112 nvr=nvr+1
113 idvr(nvr)=i ! index of not fixed var.
114 endif
115
116 it=ityvr(i)
117
118 if (it.eq.3) then ! torsion
119 vlvr(i)=toat(iatvr(i))
120 elseif (it.eq.2) then ! b.angle
121 vlvr(i)=baat(iatvr(i))
122 elseif (it.eq.1) then ! b.length
123 vlvr(i)=blat(iatvr(i))
124 endif
125
126 olvlvr(i) = vlvr(i)
127 enddo
128
129 ireg = 0
130
131 else ! =========================== from PDB
132 if (iendst(seqfile).le.1) then
133 3 write (*,'(/,a,$)') ' PDB-file:'
134 seqfil=' '
135 read (*,'(a)',err=3) seqfil
136 else
137 seqfil = seqfile
138 endif
139 write (*,*) 'PDB structure ',seqfil(1:iendst(seqfil))
140 print *, 'calling readpdb with ',seqfile
141 call pdbread(seqfil,ier)
142
143 if (ier.ne.0) stop
144
145 call pdbvars()
146
147 ireg = 1
148
149 endif
150
151! -------------------------- set var. amplitudes for simulations
152
153 do i=1,ivrml1(ntlml)+nvrml(ntlml)-1
154
155 if (ityvr(i).eq.3.and..not.fxvr(i)) then ! torsion
156
157 navr = nmvr(i)
158
159 ir = nursvr(i)
160 nars = seq(ir)
161
162 if ( navr(1:2).eq.'om'
163
164 & .or.nars(1:3).eq.'arg'.and.(navr(1:2).eq.'x5'
165 & .or.navr(1:2).eq.'x6')
166
167 & .or.(nars(1:3).eq.'asn'.or.nars(1:3).eq.'asp')
168 & .and.navr(1:2).eq.'x3'
169
170 & .or.(nars(1:3).eq.'gln'.or.nars(1:3).eq.'glu')
171 & .and.navr(1:2).eq.'x4'
172
173 & ) then
174
175! axvr(i) = pi/9.d0 ! 20 deg.
176 axvr(i) = pi2 ! Trying out 360 deg. for these as well
177
178 else
179 axvr(i) = pi2 ! 360 deg.
180 endif
181
182 else
183 axvr(i) = 0.d0
184 endif
185
186 enddo ! vars.
187
188! --------------------- initialize solvation pars. if necessary
189
190 if (itysol.ne.0) then
191
192 i1=iatrs1(irsml1(1)) ! 1st atom of 1st molecule
193 i2=iatrs2(irsml2(ntlml)) ! last atom of last molecule
194
195 its = iabs(itysol)
196
197 do i=i1,i2 ! all atoms
198 it=ityat(i)
199 sigma(i)=coef_sl(its,it)
200 rvdw(i) =rad_vdw(its,it)
201
202 if (nmat(i)(1:1).ne.'h') rvdw(i)=rvdw(i)+rwater
203
204 enddo
205
206 endif
207! Initialize calpha array
208 do i=ontlml, ntlml
209 call c_alfa(i,1)
210 enddo
211
212! Initialize arrays used in the BGS update
213 call init_lund()
214 return
215 end
216
217
Note: See TracBrowser for help on using the repository browser.