1 subroutine intcartderiv
2 implicit real*8 (a-h,o-z)
11 include 'COMMON.INTERACT'
12 include 'COMMON.DERIV'
13 include 'COMMON.IOUNITS'
14 include 'COMMON.LOCAL'
15 double precision dcostheta(3,2,maxres),
16 & dcosphi(3,3,maxres),dsinphi(3,3,maxres),
17 & dcosalpha(3,3,maxres),dcosomega(3,3,maxres),
18 & dsinomega(3,3,maxres),vo1(3),vo2(3),vo3(3),
19 & dummy(3),vp1(3),vp2(3),vp3(3),vpp1(3),n(3)
21 #if defined(MPI) && defined(PARINTDER)
22 if (nfgtasks.gt.1 .and. me.eq.king)
23 & call MPI_Bcast(8,1,MPI_INTEGER,king,FG_COMM,IERROR)
28 c write (iout,*) "iphi1_start",iphi1_start," iphi1_end",iphi1_end
29 c Derivatives of theta's
30 #if defined(MPI) && defined(PARINTDER)
31 c We need dtheta(:,:,i-1) to compute dphi(:,:,i)
32 do i=max0(ithet_start-1,3),ithet_end
37 sint=sqrt(1-cost*cost)
39 dcostheta(j,1,i)=-(dc_norm(j,i-1)+cost*dc_norm(j,i-2))/
41 dtheta(j,1,i)=-1/sint*dcostheta(j,1,i)
42 dcostheta(j,2,i)=-(dc_norm(j,i-2)+cost*dc_norm(j,i-1))/
44 dtheta(j,2,i)=-1/sint*dcostheta(j,2,i)
48 #if defined(MPI) && defined(PARINTDER)
49 c We need dtheta(:,:,i-1) to compute dphi(:,:,i)
50 do i=max0(ithet_start-1,3),ithet_end
54 if (itype(i-1).ne.10) then
55 cost1=dcos(omicron(1,i))
56 sint1=sqrt(1-cost1*cost1)
57 cost2=dcos(omicron(2,i))
58 sint2=sqrt(1-cost2*cost2)
60 CC Calculate derivative over first omicron (Cai-2,Cai-1,SCi-1)
61 dcosomicron(j,1,1,i)=-(-dc_norm(j,i-1+nres)+
62 & cost1*dc_norm(j,i-2))/
64 domicron(j,1,1,i)=-1/sint1*dcosomicron(j,1,1,i)
65 dcosomicron(j,1,2,i)=-(dc_norm(j,i-2)
66 & +cost1*(-dc_norm(j,i-1+nres)))/
68 domicron(j,1,2,i)=-1/sint1*dcosomicron(j,1,2,i)
69 CC Calculate derivative over second omicron Sci-1,Cai-1 Cai
70 CC Looks messy but better than if in loop
71 dcosomicron(j,2,1,i)=-(-dc_norm(j,i-1+nres)
72 & +cost2*dc_norm(j,i-1))/
74 domicron(j,2,1,i)=-1/sint2*dcosomicron(j,2,1,i)
75 dcosomicron(j,2,2,i)=-(dc_norm(j,i-1)
76 & +cost1*(-dc_norm(j,i-1+nres)))/
78 domicron(j,2,2,i)=-1/sint2*dcosomicron(j,2,2,i)
86 c If phi is 0 or 180 degrees, then the formulas
87 c have to be derived by power series expansion of the
88 c conventional formulas around 0 and 180.
90 do i=iphi1_start,iphi1_end
94 c the conventional case
96 sint1=dsin(theta(i-1))
99 cost1=dcos(theta(i-1))
101 scalp=scalar(dc_norm(1,i-3),dc_norm(1,i-1))
102 fac0=1.0d0/(sint1*sint)
105 fac3=cosg*cost1/(sint1*sint1)
106 fac4=cosg*cost/(sint*sint)
107 c Obtaining the gamma derivatives from sine derivative
108 if (phi(i).gt.-pi4.and.phi(i).le.pi4.or.
109 & phi(i).gt.pi34.and.phi(i).le.pi.or.
110 & phi(i).gt.-pi.and.phi(i).le.-pi34) then
111 call vecpr(dc_norm(1,i-1),dc_norm(1,i-2),vp1)
112 call vecpr(dc_norm(1,i-3),dc_norm(1,i-1),vp2)
113 call vecpr(dc_norm(1,i-3),dc_norm(1,i-2),vp3)
118 dsinphi(j,1,i)=-sing*ctgt1*dtheta(j,1,i-1)
119 & -(fac0*vp1(j)+sing*dc_norm(j,i-3))*vbld_inv(i-2)
120 dphi(j,1,i)=cosg_inv*dsinphi(j,1,i)
122 & -sing*(ctgt1*dtheta(j,2,i-1)+ctgt*dtheta(j,1,i))
123 & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
124 dphi(j,2,i)=cosg_inv*dsinphi(j,2,i)
125 c Bug fixed 3/24/05 (AL)
126 dsinphi(j,3,i)=-sing*ctgt*dtheta(j,2,i)
127 & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i)
128 c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
129 dphi(j,3,i)=cosg_inv*dsinphi(j,3,i)
131 c Obtaining the gamma derivatives from cosine derivative
134 dcosphi(j,1,i)=fac1*dcostheta(j,1,i-1)+fac3*
135 & dcostheta(j,1,i-1)-fac0*(dc_norm(j,i-1)-scalp*
136 & dc_norm(j,i-3))/vbld(i-2)
137 dphi(j,1,i)=-1/sing*dcosphi(j,1,i)
138 dcosphi(j,2,i)=fac1*dcostheta(j,2,i-1)+fac2*
139 & dcostheta(j,1,i)+fac3*dcostheta(j,2,i-1)+fac4*
141 dphi(j,2,i)=-1/sing*dcosphi(j,2,i)
142 dcosphi(j,3,i)=fac2*dcostheta(j,2,i)+fac4*
143 & dcostheta(j,2,i)-fac0*(dc_norm(j,i-3)-scalp*
144 & dc_norm(j,i-1))/vbld(i)
145 dphi(j,3,i)=-1/sing*dcosphi(j,3,i)
150 Calculate derivative of Tauangle
152 do i=iphi1_start-1,iphi1_end
156 cc dtauangle(j,intertyp,dervityp,residue number)
157 cc INTERTYP=1 SC...Ca...Ca..Ca
158 c the conventional case
160 sint1=dsin(omicron(2,i-1))
161 sing=dsin(tauangle(1,i))
163 cost1=dcos(omicron(2,i-1))
164 cosg=dcos(tauangle(1,i))
166 dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres)
168 scalp=scalar(dc_norm2(1,i-2+nres),dc_norm(1,i-1))
169 fac0=1.0d0/(sint1*sint)
172 fac3=cosg*cost1/(sint1*sint1)
173 fac4=cosg*cost/(sint*sint)
174 c Obtaining the gamma derivatives from sine derivative
175 if (tauangle(1,i).gt.-pi4.and.tauangle(1,i).le.pi4.or.
176 & tauangle(1,i).gt.pi34.and.tauangle(1,i).le.pi.or.
177 & tauangle(1,i).gt.-pi.and.tauangle(1,i).le.-pi34) then
179 dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres)
181 call vecpr(dc_norm(1,i-1),dc_norm(1,i-2),vp1)
182 call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-1),vp2)
183 call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-2),vp3)
188 dsintau(j,1,1,i)=-sing*ctgt1*domicron(j,2,1,i-1)
189 &-(fac0*vp1(j)+sing*(-dc_norm(j,i-2+nres)))
190 & *vbld_inv(i-2+nres)
191 dtauangle(j,1,1,i)=cosg_inv*dsintau(j,1,1,i)
193 & -sing*(ctgt1*domicron(j,2,2,i-1)+ctgt*dtheta(j,1,i))
194 & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
195 dtauangle(j,1,2,i)=cosg_inv*dsintau(j,1,2,i)
196 c Bug fixed 3/24/05 (AL)
197 dsintau(j,1,3,i)=-sing*ctgt*dtheta(j,2,i)
198 & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i)
199 c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
200 dtauangle(j,1,3,i)=cosg_inv*dsinphi(j,1,3,i)
202 c Obtaining the gamma derivatives from cosine derivative
205 dcostau(j,1,1,i)=fac1*dcosomicron(j,2,1,i-1)+fac3*
206 & dcosomicron(j,2,1,i-1)-fac0*(dc_norm(j,i-1)-scalp*
207 & (-dc_norm(j,i-2+nres)))/vbld(i-2+nres)
208 dtauangle(j,1,1,i)=-1/sing*dcostau(j,1,1,i)
209 dcostau(j,1,2,i)=fac1*dcosomicron(j,2,2,i-1)+fac2*
210 & dcostheta(j,1,i)+fac3*dcosomicron(j,2,2,i-1)+fac4*
212 dtauangle(j,1,2,i)=-1/sing*dcostau(j,1,2,i)
213 dcostau(j,1,3,i)=fac2*dcostheta(j,2,i)+fac4*
214 & dcostheta(j,2,i)-fac0*(-dc_norm(j,i-2+nres)-scalp*
215 & dc_norm(j,i-1))/vbld(i)
216 dtauangle(j,1,3,i)=-1/sing*dcostau(j,1,3,i)
220 CC Second case Ca...Ca...Ca...SC
222 do i=iphi1_start,iphi1_end
226 c the conventional case
227 sint=dsin(omicron(1,i))
228 sint1=dsin(theta(i-1))
229 sing=dsin(tauangle(2,i)
230 cost=dcos(omicron(1,i)
231 cost1=dcos(theta(i-1))
232 cosg=dcos(tauangle(2,i))
234 dc_norm2(j,i-1+nres)=-dc_norm(j,i-1+nres)
236 scalp=scalar(dc_norm(1,i-3),dc_norm2(1,i-1+nres))
237 fac0=1.0d0/(sint1*sint)
240 fac3=cosg*cost1/(sint1*sint1)
241 fac4=cosg*cost/(sint*sint)
242 c Obtaining the gamma derivatives from sine derivative
243 if (tauangle(2,i).gt.-pi4.and.tauangle(2,i).le.pi4.or.
244 & tauangle(2,i).gt.pi34.and.tauangle(2,i).le.pi.or.
245 & tauangle(2,i).gt.-pi.and.tauangle(2,i).le.-pi34) then
246 call vecpr(dc_norm2(1,i-1+nres),dc_norm(1,i-2),vp1)
247 call vecpr(dc_norm(1,i-3),dc_norm2(1,i-1+nres),vp2)
248 call vecpr(dc_norm(1,i-3),dc_norm(1,i-2),vp3)
253 dsintau(j,2,1,i)=-sing*ctgt1*dtheta(j,1,i-1)
254 & -(fac0*vp1(j)+sing*dc_norm(j,i-3))*vbld_inv(i-2)
255 dtauangle(j,2,1,i)=cosg_inv*dsintau(j,2,1,i)
257 & -sing*(ctgt1*dtheta(j,2,i-1)+ctgt*domicron(j,1,1,i))
258 & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
259 dphi(j,2,i)=cosg_inv*dsintau(j,2,2,i)
260 c Bug fixed 3/24/05 (AL)
261 dsintau(j,2,3,i)=-sing*ctgt*domicron(j,1,2,i)
262 & +(fac0*vp3(j)-sing*dc_norm2(j,i-1+nres))*vbld_inv(i)
263 c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
264 dtauangle(j,2,3,i)=cosg_inv*dsinphi(j,2,3,i)
266 c Obtaining the gamma derivatives from cosine derivative
269 dcostau(j,2,1,i)=fac1*dcostheta(j,1,i-1)+fac3*
270 & dcostheta(j,1,i-1)-fac0*(dc_norm2(j,i-1+nres)-scalp*
271 & dc_norm(j,i-3))/vbld(i-2)
272 dtauangle(j,2,1,i)=-1/sing*dcostau(j,2,1,i)
273 dcostau(j,2,2,i)=fac1*dcostheta(j,2,i-1)+fac2*
274 & dcosomicron(j,1,1,i)+fac3*dcostheta(j,2,i-1)+fac4*
275 & dcosomicron(j,1,1,i)
276 dtauanlge(j,2,2,i)=-1/sing*dcostau(j,2,2,i)
277 dcostau(j,2,3,i)=fac2*dcosomicron(j,1,2,i)+fac4*
278 & dcostheta(j,1,2,i)-fac0*(dc_norm(j,i-3)-scalp*
279 & dc_norm2(j,i-1+nres))/vbld(i-1+nres)
280 dtauanlge(j,2,3,i)=-1/sing*dcosphi(j,3,i)
286 CCC third case SC...Ca...Ca...SC
289 do i=iphi1_start-1,iphi1_end
293 c the conventional case
294 sint=dsin(omicron(1,i))
295 sint1=dsin(omicron(2,i-1))
296 sing=dsin(tauangle(3,i))
297 cost=dcos(omicron(1,i))
298 cost1=dcos(omicron(2,i-1))
299 cosg=dcos(tauangle(3,i))
301 dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres)
302 dc_norm2(j,i-1+nres)=-dc_norm(j,i-2+nres)
304 scalp=scalar(dc_norm2(1,i-2+nres),dc_norm2(1,i-1+nres))
305 fac0=1.0d0/(sint1*sint)
308 fac3=cosg*cost1/(sint1*sint1)
309 fac4=cosg*cost/(sint*sint)
310 c Obtaining the gamma derivatives from sine derivative
311 if (tauangle(3,i).gt.-pi4.and.tauangle(3,i).le.pi4.or.
312 & tauangle(3,i).gt.pi34.and.tauangle(3,i).le.pi.or.
313 & tauangle(3,i).gt.-pi.and.tauangle(3,i).le.-pi34) then
314 call vecpr(dc_norm2(1,i-1+nres),dc_norm(1,i-2),vp1)
315 call vecpr(dc_norm2(1,i-2+nres),dc_norm2(1,i-1+nres),vp2)
316 call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-2),vp3)
321 dsintau(j,3,1,i)=-sing*ctgt1*domicron(j,2,1,i-1)
322 & -(fac0*vp1(j)+sing*dc_norm2(j,i-2+nres))
323 & *vbld_inv(i-2+nres)
324 dtauangle(j,3,1,i)=cosg_inv*dsintau(j,3,1,i)
326 & -sing*(ctgt1*domicron(j,2,2,i-1)+ctgt*domicron(j,1,1,i))
327 & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1+nres)
328 dtauangle(j,3,2,i)=cosg_inv*dsintau(j,3,2,i)
329 c Bug fixed 3/24/05 (AL)
330 dsintau(j,3,3,i)=-sing*ctgt*domicron(j,1,2,i)
331 & +(fac0*vp3(j)-sing*dc_norm2(j,i-1+nres))
332 & *vbld_inv(i-1+nres)
333 c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
334 dphi(j,3,3,i)=cosg_inv*dsintau(j,3,3,i)
336 c Obtaining the gamma derivatives from cosine derivative
339 dcostau(j,3,1,i)=fac1*dcosomicron(j,2,1,i-1)+fac3*
340 & dcostheta(j,1,i-1)-fac0*(dc_norm2(j,i-1+nres)-scalp*
341 & dc_norm2(j,i-2+nres))/vbld(i-2+nres)
342 dtauangle(j,3,1,i)=-1/sing*dcostau(j,3,1,i)
343 dcostau(j,3,2,i)=fac1*dcosomicron(j,2,2,i-1)+fac2*
344 & dcosomicron(j,1,1,i)+fac3*dcosomicron(j,2,2,i-1)+fac4*
345 & dcosomicron(j,1,1,i)
346 dtauangle(j,3,2,i)=-1/sing*dcostau(j,3,2,i)
347 dcostau(j,3,3,i)=fac2*dcosomicron(j,1,2,i)+fac4*
348 & dcostau(j,3,2,i)-fac0*(dc_norm2(j,i-2+nres)-scalp*
349 & dc_norm2(j,i-1+nres))/vbld(i-1+nres)
350 dtauangle(j,3,3,i)=-1/sing*dcostau(j,3,3,i)
355 c Derivatives of side-chain angles alpha and omega
356 #if defined(MPI) && defined(PARINTDER)
357 do i=ibond_start,ibond_end
361 if(itype(i).ne.10) then
362 fac5=1.0d0/dsqrt(2*(1+dcos(theta(i+1))))
366 fac9=fac5/vbld(i+nres)
367 scala1=scalar(dc_norm(1,i-1),dc_norm(1,i+nres))
368 scala2=scalar(dc_norm(1,i),dc_norm(1,i+nres))
369 cosa=dsqrt(0.5d0/(1.0d0+dcos(theta(i+1))))*(
370 & scalar(dC_norm(1,i),dC_norm(1,i+nres))
371 & -scalar(dC_norm(1,i-1),dC_norm(1,i+nres)))
372 sina=sqrt(1-cosa*cosa)
375 dcosalpha(j,1,i)=fac6*(scala1*dc_norm(j,i-1)-
376 & dc_norm(j,i+nres))-cosa*fac7*dcostheta(j,1,i+1)
377 dalpha(j,1,i)=-1/sina*dcosalpha(j,1,i)
378 dcosalpha(j,2,i)=fac8*(dc_norm(j,i+nres)-
379 & scala2*dc_norm(j,i))-cosa*fac7*dcostheta(j,2,i+1)
380 dalpha(j,2,i)=-1/sina*dcosalpha(j,2,i)
381 dcosalpha(j,3,i)=(fac9*(dc_norm(j,i)-
382 & dc_norm(j,i-1))-(cosa*dc_norm(j,i+nres))/
384 dalpha(j,3,i)=-1/sina*dcosalpha(j,3,i)
386 c obtaining the derivatives of omega from sines
387 if(omeg(i).gt.-pi4.and.omeg(i).le.pi4.or.
388 & omeg(i).gt.pi34.and.omeg(i).le.pi.or.
389 & omeg(i).gt.-pi.and.omeg(i).le.-pi34) then
390 fac15=dcos(theta(i+1))/(dsin(theta(i+1))*
392 fac16=dcos(alph(i))/(dsin(alph(i))*dsin(alph(i)))
393 fac17=1.0d0/(dsin(theta(i+1))*dsin(alph(i)))
394 call vecpr(dc_norm(1,i+nres),dc_norm(1,i),vo1)
395 call vecpr(dc_norm(1,i+nres),dc_norm(1,i-1),vo2)
396 call vecpr(dc_norm(1,i),dc_norm(1,i-1),vo3)
397 coso_inv=1.0d0/dcos(omeg(i))
399 dsinomega(j,1,i)=sino*(fac15*dcostheta(j,1,i+1)
400 & +fac16*dcosalpha(j,1,i))-fac17/vbld(i)*vo1(j)-(
401 & sino*dc_norm(j,i-1))/vbld(i)
402 domega(j,1,i)=coso_inv*dsinomega(j,1,i)
403 dsinomega(j,2,i)=sino*(fac15*dcostheta(j,2,i+1)
404 & +fac16*dcosalpha(j,2,i))+fac17/vbld(i+1)*vo2(j)
405 & -sino*dc_norm(j,i)/vbld(i+1)
406 domega(j,2,i)=coso_inv*dsinomega(j,2,i)
407 dsinomega(j,3,i)=sino*fac16*dcosalpha(j,3,i)-
408 & fac17/vbld(i+nres)*vo3(j)-sino*dc_norm(j,i+nres)/
410 domega(j,3,i)=coso_inv*dsinomega(j,3,i)
413 c obtaining the derivatives of omega from cosines
414 fac10=sqrt(0.5d0*(1-dcos(theta(i+1))))
415 fac11=sqrt(0.5d0*(1+dcos(theta(i+1))))
420 dcosomega(j,1,i)=(-(0.25d0*cosa/fac11*
421 & dcostheta(j,1,i+1)+fac11*dcosalpha(j,1,i))*fac12+
422 & (0.25d0/fac10*sina*dcostheta(j,1,i+1)+cosa/sina*
423 & fac10*dcosalpha(j,1,i))*(scala2-fac11*cosa))/fac13
424 domega(j,1,i)=-1/sino*dcosomega(j,1,i)
425 dcosomega(j,2,i)=(((dc_norm(j,i+nres)-scala2*
426 & dc_norm(j,i))/vbld(i+1)-0.25d0*cosa/fac11*
427 & dcostheta(j,2,i+1)-fac11*dcosalpha(j,2,i))*fac12+
428 & (scala2-fac11*cosa)*(0.25d0*sina/fac10*
429 & dcostheta(j,2,i+1)+fac10*cosa/sina*dcosalpha(j,2,i)
431 domega(j,2,i)=-1/sino*dcosomega(j,2,i)
432 dcosomega(j,3,i)=1/fac10*((1/vbld(i+nres)*(dc_norm(j,i)-
433 & scala2*dc_norm(j,i+nres))-fac11*dcosalpha(j,3,i))*sina+
434 & (scala2-fac11*cosa)*(cosa/sina*dcosalpha(j,3,i)))/fac14
435 domega(j,3,i)=-1/sino*dcosomega(j,3,i)
441 #if defined(MPI) && defined(PARINTDER)
442 if (nfgtasks.gt.1) then
444 cd write (iout,*) "Gather dtheta"
446 write (iout,*) "dtheta before gather"
448 write (iout,'(i3,3(3f8.5,3x))') i,((dtheta(j,k,i),k=1,3),j=1,2)
451 call MPI_Gatherv(dtheta(1,1,ithet_start),ithet_count(fg_rank),
452 & MPI_THET,dtheta(1,1,1),ithet_count(0),ithet_displ(0),MPI_THET,
453 & king,FG_COMM,IERROR)
455 cd write (iout,*) "Gather dphi"
457 write (iout,*) "dphi before gather"
459 write (iout,'(i3,3(3f8.5,3x))') i,((dphi(j,k,i),k=1,3),j=1,3)
462 call MPI_Gatherv(dphi(1,1,iphi1_start),iphi1_count(fg_rank),
463 & MPI_GAM,dphi(1,1,1),iphi1_count(0),iphi1_displ(0),MPI_GAM,
464 & king,FG_COMM,IERROR)
465 cd write (iout,*) "Gather dalpha"
468 call MPI_Gatherv(dalpha(1,1,ibond_start),ibond_count(fg_rank),
469 & MPI_GAM,dalpha(1,1,1),ibond_count(0),ibond_displ(0),MPI_GAM,
470 & king,FG_COMM,IERROR)
471 cd write (iout,*) "Gather domega"
473 call MPI_Gatherv(domega(1,1,ibond_start),ibond_count(fg_rank),
474 & MPI_GAM,domega(1,1,1),ibond_count(0),ibond_displ(0),MPI_GAM,
475 & king,FG_COMM,IERROR)
480 write (iout,*) "dtheta after gather"
482 write (iout,'(i3,3(3f8.5,3x))') i,((dtheta(j,k,i),j=1,3),j=1,2)
484 write (iout,*) "dphi after gather"
486 write (iout,'(i3,3(3f8.5,3x))') i,((dphi(j,k,i),j=1,3),k=1,3)
492 subroutine checkintcartgrad
493 implicit real*8 (a-h,o-z)
498 include 'COMMON.CHAIN'
501 include 'COMMON.INTERACT'
502 include 'COMMON.DERIV'
503 include 'COMMON.IOUNITS'
504 include 'COMMON.SETUP'
505 double precision dthetanum(3,2,maxres),dphinum(3,3,maxres)
506 & ,dalphanum(3,3,maxres), domeganum(3,3,maxres)
507 double precision theta_s(maxres),phi_s(maxres),alph_s(maxres),
508 & omeg_s(maxres),dc_norm_s(3)
509 double precision aincr /1.0d-5/
517 c Check theta gradient
519 & "Analytical (upper) and numerical (lower) gradient of theta"
526 call int_from_cart1(.false.)
527 dthetanum(j,1,i)=(theta(i)-theta_s(i))/aincr
530 dc(j,i-1)=dc(j,i-1)+aincr
532 dthetanum(j,2,i)=(theta(i)-theta_s(i))/aincr
535 write (iout,'(i5,3f10.5,5x,3f10.5)') i,(dtheta(j,1,i),j=1,3),
536 & (dtheta(j,2,i),j=1,3)
537 write (iout,'(5x,3f10.5,5x,3f10.5)') (dthetanum(j,1,i),j=1,3),
538 & (dthetanum(j,2,i),j=1,3)
539 write (iout,'(5x,3f10.5,5x,3f10.5)')
540 & (dthetanum(j,1,i)/dtheta(j,1,i),j=1,3),
541 & (dthetanum(j,2,i)/dtheta(j,2,i),j=1,3)
544 c Check gamma gradient
546 & "Analytical (upper) and numerical (lower) gradient of gamma"
552 dphinum(j,1,i)=(phi(i)-phi_s(i))/aincr
557 dphinum(j,2,i)=(phi(i)-phi_s(i))/aincr
560 dc(j,i-1)=dc(j,i-1)+aincr
562 dphinum(j,3,i)=(phi(i)-phi_s(i))/aincr
565 write (iout,'(i5,3(3f10.5,5x))') i,(dphi(j,1,i),j=1,3),
566 & (dphi(j,2,i),j=1,3),(dphi(j,3,i),j=1,3)
567 write (iout,'(5x,3(3f10.5,5x))') (dphinum(j,1,i),j=1,3),
568 & (dphinum(j,2,i),j=1,3),(dphinum(j,3,i),j=1,3)
569 write (iout,'(5x,3(3f10.5,5x))')
570 & (dphinum(j,1,i)/dphi(j,1,i),j=1,3),
571 & (dphinum(j,2,i)/dphi(j,2,i),j=1,3),
572 & (dphinum(j,3,i)/dphi(j,3,i),j=1,3)
575 c Check alpha gradient
577 & "Analytical (upper) and numerical (lower) gradient of alpha"
579 if(itype(i).ne.10) then
584 dalphanum(j,1,i)=(alph(i)-alph_s(i))
590 dalphanum(j,2,i)=(alph(i)-alph_s(i))
594 dc(j,i+nres)=dc(j,i+nres)+aincr
596 dalphanum(j,3,i)=(alph(i)-alph_s(i))
601 write (iout,'(i5,3(3f10.5,5x))') i,(dalpha(j,1,i),j=1,3),
602 & (dalpha(j,2,i),j=1,3),(dalpha(j,3,i),j=1,3)
603 write (iout,'(5x,3(3f10.5,5x))') (dalphanum(j,1,i),j=1,3),
604 & (dalphanum(j,2,i),j=1,3),(dalphanum(j,3,i),j=1,3)
605 write (iout,'(5x,3(3f10.5,5x))')
606 & (dalphanum(j,1,i)/dalpha(j,1,i),j=1,3),
607 & (dalphanum(j,2,i)/dalpha(j,2,i),j=1,3),
608 & (dalphanum(j,3,i)/dalpha(j,3,i),j=1,3)
611 c Check omega gradient
613 & "Analytical (upper) and numerical (lower) gradient of omega"
615 if(itype(i).ne.10) then
620 domeganum(j,1,i)=(omeg(i)-omeg_s(i))
626 domeganum(j,2,i)=(omeg(i)-omeg_s(i))
630 dc(j,i+nres)=dc(j,i+nres)+aincr
632 domeganum(j,3,i)=(omeg(i)-omeg_s(i))
637 write (iout,'(i5,3(3f10.5,5x))') i,(domega(j,1,i),j=1,3),
638 & (domega(j,2,i),j=1,3),(domega(j,3,i),j=1,3)
639 write (iout,'(5x,3(3f10.5,5x))') (domeganum(j,1,i),j=1,3),
640 & (domeganum(j,2,i),j=1,3),(domeganum(j,3,i),j=1,3)
641 write (iout,'(5x,3(3f10.5,5x))')
642 & (domeganum(j,1,i)/domega(j,1,i),j=1,3),
643 & (domeganum(j,2,i)/domega(j,2,i),j=1,3),
644 & (domeganum(j,3,i)/domega(j,3,i),j=1,3)
650 subroutine chainbuild_cart
651 implicit real*8 (a-h,o-z)
656 include 'COMMON.SETUP'
657 include 'COMMON.CHAIN'
658 include 'COMMON.LOCAL'
659 include 'COMMON.TIME1'
660 include 'COMMON.IOUNITS'
663 if (nfgtasks.gt.1) then
664 c write (iout,*) "BCAST in chainbuild_cart"
666 c Broadcast the order to build the chain and compute internal coordinates
667 c to the slaves. The slaves receive the order in ERGASTULUM.
669 c write (iout,*) "CHAINBUILD_CART: DC before BCAST"
671 c write (iout,'(i3,3f10.5,5x,3f10.5)') i,(dc(j,i),j=1,3),
672 c & (dc(j,i+nres),j=1,3)
675 & call MPI_Bcast(7,1,MPI_INTEGER,king,FG_COMM,IERROR)
676 time_bcast7=time_bcast7+MPI_Wtime()-time00
678 call MPI_Bcast(dc(1,0),6*(nres+1),MPI_DOUBLE_PRECISION,
680 c write (iout,*) "CHAINBUILD_CART: DC after BCAST"
682 c write (iout,'(i3,3f10.5,5x,3f10.5)') i,(dc(j,i),j=1,3),
683 c & (dc(j,i+nres),j=1,3)
685 c write (iout,*) "End BCAST in chainbuild_cart"
687 time_bcast=time_bcast+MPI_Wtime()-time00
688 time_bcastc=time_bcastc+MPI_Wtime()-time01
696 c(j,i)=c(j,i-1)+dc(j,i-1)
701 c(j,i+nres)=c(j,i)+dc(j,i+nres)
704 c write (iout,*) "CHAINBUILD_CART"
706 call int_from_cart1(.false.)