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 include 'COMMON.SCCOR'
16 double precision dcostheta(3,2,maxres),
17 & dcosphi(3,3,maxres),dsinphi(3,3,maxres),
18 & dcosalpha(3,3,maxres),dcosomega(3,3,maxres),
19 & dsinomega(3,3,maxres),vo1(3),vo2(3),vo3(3),
20 & dummy(3),vp1(3),vp2(3),vp3(3),vpp1(3),n(3)
22 #if defined(MPI) && defined(PARINTDER)
23 if (nfgtasks.gt.1 .and. me.eq.king)
24 & call MPI_Bcast(8,1,MPI_INTEGER,king,FG_COMM,IERROR)
29 c write (iout,*) "iphi1_start",iphi1_start," iphi1_end",iphi1_end
39 c Derivatives of thetas
40 #if defined(MPI) && defined(PARINTDER)
41 c We need dtheta(:,:,i-1) to compute dphi(:,:,i)
42 do i=max0(ithet_start-1,3),ithet_end
47 sint=sqrt(1-cost*cost)
49 dcostheta(j,1,i)=-(dc_norm(j,i-1)+cost*dc_norm(j,i-2))/
51 c if (itype(i-1).ne.ntyp1)
52 dtheta(j,1,i)=-dcostheta(j,1,i)/sint
53 dcostheta(j,2,i)=-(dc_norm(j,i-2)+cost*dc_norm(j,i-1))/
55 c if (itype(i-1).ne.ntyp1)
56 dtheta(j,2,i)=-dcostheta(j,2,i)/sint
59 #if defined(MPI) && defined(PARINTDER)
60 c We need dtheta(:,:,i-1) to compute dphi(:,:,i)
61 do i=max0(ithet_start-1,3),ithet_end
65 if ((itype(i-1).ne.10).and.(itype(i-1).ne.ntyp1)) then
66 cost1=dcos(omicron(1,i))
67 sint1=sqrt(1-cost1*cost1)
68 cost2=dcos(omicron(2,i))
69 sint2=sqrt(1-cost2*cost2)
71 CC Calculate derivative over first omicron (Cai-2,Cai-1,SCi-1)
72 dcosomicron(j,1,1,i)=-(dc_norm(j,i-1+nres)+
73 & cost1*dc_norm(j,i-2))/
75 domicron(j,1,1,i)=-1/sint1*dcosomicron(j,1,1,i)
76 dcosomicron(j,1,2,i)=-(dc_norm(j,i-2)
77 & +cost1*(dc_norm(j,i-1+nres)))/
79 domicron(j,1,2,i)=-1/sint1*dcosomicron(j,1,2,i)
80 CC Calculate derivative over second omicron Sci-1,Cai-1 Cai
81 CC Looks messy but better than if in loop
82 dcosomicron(j,2,1,i)=-(-dc_norm(j,i-1+nres)
83 & +cost2*dc_norm(j,i-1))/
85 domicron(j,2,1,i)=-1/sint2*dcosomicron(j,2,1,i)
86 dcosomicron(j,2,2,i)=-(dc_norm(j,i-1)
87 & +cost2*(-dc_norm(j,i-1+nres)))/
89 c write(iout,*) "vbld", i,itype(i),vbld(i-1+nres)
90 domicron(j,2,2,i)=-1/sint2*dcosomicron(j,2,2,i)
96 c If phi is 0 or 180 degrees, then the formulas
97 c have to be derived by power series expansion of the
98 c conventional formulas around 0 and 180.
100 do i=iphi1_start,iphi1_end
104 c if (itype(i-2).eq.ntyp1.or. itype(i-1).eq.ntyp1
105 c & .or. itype(i).eq.ntyp1 .or. itype(i-3).eq.ntyp1) cycle
106 c the conventional case
108 sint1=dsin(theta(i-1))
111 cost1=dcos(theta(i-1))
113 scalp=scalar(dc_norm(1,i-3),dc_norm(1,i-1))
114 fac0=1.0d0/(sint1*sint)
117 fac3=cosg*cost1/(sint1*sint1)
118 fac4=cosg*cost/(sint*sint)
119 c Obtaining the gamma derivatives from sine derivative
120 if (phi(i).gt.-pi4.and.phi(i).le.pi4.or.
121 & phi(i).gt.pi34.and.phi(i).le.pi.or.
122 & phi(i).ge.-pi.and.phi(i).le.-pi34) then
123 call vecpr(dc_norm(1,i-1),dc_norm(1,i-2),vp1)
124 call vecpr(dc_norm(1,i-3),dc_norm(1,i-1),vp2)
125 call vecpr(dc_norm(1,i-3),dc_norm(1,i-2),vp3)
130 c if (itype(i-1).ne.ntyp1 .and. itype(i-2).ne.ntyp1) then
131 dsinphi(j,1,i)=-sing*ctgt1*dtheta(j,1,i-1)
132 & -(fac0*vp1(j)+sing*dc_norm(j,i-3))*vbld_inv(i-2)
133 dphi(j,1,i)=cosg_inv*dsinphi(j,1,i)
135 & -sing*(ctgt1*dtheta(j,2,i-1)+ctgt*dtheta(j,1,i))
136 & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
137 dphi(j,2,i)=cosg_inv*dsinphi(j,2,i)
138 dsinphi(j,3,i)=-sing*ctgt*dtheta(j,2,i)
139 & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i)
140 c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
141 dphi(j,3,i)=cosg_inv*dsinphi(j,3,i)
143 c Bug fixed 3/24/05 (AL)
145 c Obtaining the gamma derivatives from cosine derivative
148 c if (itype(i-1).ne.ntyp1 .and. itype(i-2).ne.ntyp1) then
149 dcosphi(j,1,i)=fac1*dcostheta(j,1,i-1)+fac3*
150 & dcostheta(j,1,i-1)-fac0*(dc_norm(j,i-1)-scalp*
151 & dc_norm(j,i-3))/vbld(i-2)
152 dphi(j,1,i)=-1/sing*dcosphi(j,1,i)
153 dcosphi(j,2,i)=fac1*dcostheta(j,2,i-1)+fac2*
154 & dcostheta(j,1,i)+fac3*dcostheta(j,2,i-1)+fac4*
156 dphi(j,2,i)=-1/sing*dcosphi(j,2,i)
157 dcosphi(j,3,i)=fac2*dcostheta(j,2,i)+fac4*
158 & dcostheta(j,2,i)-fac0*(dc_norm(j,i-3)-scalp*
159 & dc_norm(j,i-1))/vbld(i)
160 dphi(j,3,i)=-1/sing*dcosphi(j,3,i)
165 Calculate derivative of Tauangle
168 dc_norm2(j,i+nres)=-dc_norm(j,i+nres)
172 do i=itau_start,itau_end
176 if ((itype(i-2).eq.ntyp1).or.(itype(i-2).eq.10)) cycle
177 c if ((itype(i-2).eq.ntyp1).or.(itype(i-2).eq.10).or.
178 c & (itype(i-1).eq.ntyp1).or.(itype(i).eq.ntyp1)) cycle
179 cc dtauangle(j,intertyp,dervityp,residue number)
180 cc INTERTYP=1 SC...Ca...Ca..Ca
181 c the conventional case
183 sint1=dsin(omicron(2,i-1))
184 sing=dsin(tauangle(1,i))
186 cost1=dcos(omicron(2,i-1))
187 cosg=dcos(tauangle(1,i))
189 dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres)
190 cc write(iout,*) dc_norm2(j,i-2+nres),"dcnorm"
192 scalp=scalar(dc_norm2(1,i-2+nres),dc_norm(1,i-1))
193 fac0=1.0d0/(sint1*sint)
196 fac3=cosg*cost1/(sint1*sint1)
197 fac4=cosg*cost/(sint*sint)
198 cc write(iout,*) "faki",fac0,fac1,fac2,fac3,fac4
199 c Obtaining the gamma derivatives from sine derivative
200 if (tauangle(1,i).gt.-pi4.and.tauangle(1,i).le.pi4.or.
201 & tauangle(1,i).gt.pi34.and.tauangle(1,i).le.pi.or.
202 & tauangle(1,i).ge.-pi.and.tauangle(1,i).le.-pi34) then
203 call vecpr(dc_norm(1,i-1),dc_norm(1,i-2),vp1)
204 call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-1),vp2)
205 call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-2),vp3)
210 dsintau(j,1,1,i)=-sing*ctgt1*domicron(j,2,2,i-1)
211 &-(fac0*vp1(j)+sing*(dc_norm2(j,i-2+nres)))
212 & *vbld_inv(i-2+nres)
213 dtauangle(j,1,1,i)=cosg_inv*dsintau(j,1,1,i)
215 & -sing*(ctgt1*domicron(j,2,1,i-1)+ctgt*dtheta(j,1,i))
216 & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
217 c write(iout,*) "dsintau", dsintau(j,1,2,i)
218 dtauangle(j,1,2,i)=cosg_inv*dsintau(j,1,2,i)
219 c Bug fixed 3/24/05 (AL)
220 dsintau(j,1,3,i)=-sing*ctgt*dtheta(j,2,i)
221 & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i)
222 c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
223 dtauangle(j,1,3,i)=cosg_inv*dsintau(j,1,3,i)
225 c Obtaining the gamma derivatives from cosine derivative
228 dcostau(j,1,1,i)=fac1*dcosomicron(j,2,2,i-1)+fac3*
229 & dcosomicron(j,2,2,i-1)-fac0*(dc_norm(j,i-1)-scalp*
230 & (dc_norm2(j,i-2+nres)))/vbld(i-2+nres)
231 dtauangle(j,1,1,i)=-1/sing*dcostau(j,1,1,i)
232 dcostau(j,1,2,i)=fac1*dcosomicron(j,2,1,i-1)+fac2*
233 & dcostheta(j,1,i)+fac3*dcosomicron(j,2,1,i-1)+fac4*
235 dtauangle(j,1,2,i)=-1/sing*dcostau(j,1,2,i)
236 dcostau(j,1,3,i)=fac2*dcostheta(j,2,i)+fac4*
237 & dcostheta(j,2,i)-fac0*(-dc_norm(j,i-2+nres)-scalp*
238 & dc_norm(j,i-1))/vbld(i)
239 dtauangle(j,1,3,i)=-1/sing*dcostau(j,1,3,i)
240 c write (iout,*) "else",i
244 c write(iout,*) "tu",i,k,(dtauangle(j,1,k,i),j=1,3)
247 CC Second case Ca...Ca...Ca...SC
249 do i=itau_start,itau_end
253 if ((itype(i-1).eq.ntyp1).or.(itype(i-1).eq.10).or.
254 & (itype(i-2).eq.ntyp1).or.(itype(i-3).eq.ntyp1)) cycle
255 c the conventional case
256 sint=dsin(omicron(1,i))
257 sint1=dsin(theta(i-1))
258 sing=dsin(tauangle(2,i))
259 cost=dcos(omicron(1,i))
260 cost1=dcos(theta(i-1))
261 cosg=dcos(tauangle(2,i))
263 c dc_norm2(j,i-1+nres)=-dc_norm(j,i-1+nres)
265 scalp=scalar(dc_norm(1,i-3),dc_norm(1,i-1+nres))
266 fac0=1.0d0/(sint1*sint)
269 fac3=cosg*cost1/(sint1*sint1)
270 fac4=cosg*cost/(sint*sint)
271 c Obtaining the gamma derivatives from sine derivative
272 if (tauangle(2,i).gt.-pi4.and.tauangle(2,i).le.pi4.or.
273 & tauangle(2,i).gt.pi34.and.tauangle(2,i).le.pi.or.
274 & tauangle(2,i).gt.-pi.and.tauangle(2,i).le.-pi34) then
275 call vecpr(dc_norm2(1,i-1+nres),dc_norm(1,i-2),vp1)
276 call vecpr(dc_norm(1,i-3),dc_norm(1,i-1+nres),vp2)
277 call vecpr(dc_norm(1,i-3),dc_norm(1,i-2),vp3)
282 dsintau(j,2,1,i)=-sing*ctgt1*dtheta(j,1,i-1)
283 & +(fac0*vp1(j)-sing*dc_norm(j,i-3))*vbld_inv(i-2)
284 c write(iout,*) i,j,dsintau(j,2,1,i),sing*ctgt1*dtheta(j,1,i-1),
285 c &fac0*vp1(j),sing*dc_norm(j,i-3),vbld_inv(i-2),"dsintau(2,1)"
286 dtauangle(j,2,1,i)=cosg_inv*dsintau(j,2,1,i)
288 & -sing*(ctgt1*dtheta(j,2,i-1)+ctgt*domicron(j,1,1,i))
289 & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
290 c write(iout,*) "sprawdzenie",i,j,sing*ctgt1*dtheta(j,2,i-1),
291 c & sing*ctgt*domicron(j,1,2,i),
292 c & (fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
293 dtauangle(j,2,2,i)=cosg_inv*dsintau(j,2,2,i)
294 c Bug fixed 3/24/05 (AL)
295 dsintau(j,2,3,i)=-sing*ctgt*domicron(j,1,2,i)
296 & +(fac0*vp3(j)-sing*dc_norm(j,i-1+nres))*vbld_inv(i-1+nres)
297 c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
298 dtauangle(j,2,3,i)=cosg_inv*dsintau(j,2,3,i)
300 c Obtaining the gamma derivatives from cosine derivative
303 dcostau(j,2,1,i)=fac1*dcostheta(j,1,i-1)+fac3*
304 & dcostheta(j,1,i-1)-fac0*(dc_norm(j,i-1+nres)-scalp*
305 & dc_norm(j,i-3))/vbld(i-2)
306 dtauangle(j,2,1,i)=-1/sing*dcostau(j,2,1,i)
307 dcostau(j,2,2,i)=fac1*dcostheta(j,2,i-1)+fac2*
308 & dcosomicron(j,1,1,i)+fac3*dcostheta(j,2,i-1)+fac4*
309 & dcosomicron(j,1,1,i)
310 dtauangle(j,2,2,i)=-1/sing*dcostau(j,2,2,i)
311 dcostau(j,2,3,i)=fac2*dcosomicron(j,1,2,i)+fac4*
312 & dcosomicron(j,1,2,i)-fac0*(dc_norm(j,i-3)-scalp*
313 & dc_norm(j,i-1+nres))/vbld(i-1+nres)
314 dtauangle(j,2,3,i)=-1/sing*dcostau(j,2,3,i)
315 c write(iout,*) i,j,"else", dtauangle(j,2,3,i)
320 CCC third case SC...Ca...Ca...SC
323 do i=itau_start,itau_end
327 c the conventional case
328 if ((itype(i-1).eq.ntyp1).or.(itype(i-1).eq.10).or.
329 &(itype(i-2).eq.ntyp1).or.(itype(i-2).eq.10)) cycle
330 sint=dsin(omicron(1,i))
331 sint1=dsin(omicron(2,i-1))
332 sing=dsin(tauangle(3,i))
333 cost=dcos(omicron(1,i))
334 cost1=dcos(omicron(2,i-1))
335 cosg=dcos(tauangle(3,i))
336 c write (iout,*) "i",i," omicron",omicron(1,i),omicron(2,i)
337 c write (iout,*) "i",i," tauangle",tauangle(1,i),tauangle(2,i),
340 dc_norm2(j,i-2+nres)=-dc_norm(j,i-2+nres)
341 c dc_norm2(j,i-1+nres)=-dc_norm(j,i-1+nres)
343 scalp=scalar(dc_norm2(1,i-2+nres),dc_norm(1,i-1+nres))
344 fac0=1.0d0/(sint1*sint)
347 fac3=cosg*cost1/(sint1*sint1)
348 fac4=cosg*cost/(sint*sint)
349 c Obtaining the gamma derivatives from sine derivative
350 if (tauangle(3,i).gt.-pi4.and.tauangle(3,i).le.pi4.or.
351 & tauangle(3,i).gt.pi34.and.tauangle(3,i).le.pi.or.
352 & tauangle(3,i).ge.-pi.and.tauangle(3,i).le.-pi34) then
353 call vecpr(dc_norm(1,i-1+nres),dc_norm(1,i-2),vp1)
354 call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-1+nres),vp2)
355 call vecpr(dc_norm2(1,i-2+nres),dc_norm(1,i-2),vp3)
360 dsintau(j,3,1,i)=-sing*ctgt1*domicron(j,2,2,i-1)
361 & -(fac0*vp1(j)-sing*dc_norm(j,i-2+nres))
362 & *vbld_inv(i-2+nres)
363 dtauangle(j,3,1,i)=cosg_inv*dsintau(j,3,1,i)
365 & -sing*(ctgt1*domicron(j,2,1,i-1)+ctgt*domicron(j,1,1,i))
366 & -(fac0*vp2(j)+sing*dc_norm(j,i-2))*vbld_inv(i-1)
367 dtauangle(j,3,2,i)=cosg_inv*dsintau(j,3,2,i)
368 c Bug fixed 3/24/05 (AL)
369 dsintau(j,3,3,i)=-sing*ctgt*domicron(j,1,2,i)
370 & +(fac0*vp3(j)-sing*dc_norm(j,i-1+nres))
371 & *vbld_inv(i-1+nres)
372 c & +(fac0*vp3(j)-sing*dc_norm(j,i-1))*vbld_inv(i-1)
373 dtauangle(j,3,3,i)=cosg_inv*dsintau(j,3,3,i)
375 c Obtaining the gamma derivatives from cosine derivative
378 dcostau(j,3,1,i)=fac1*dcosomicron(j,2,2,i-1)+fac3*
379 & dcosomicron(j,2,2,i-1)-fac0*(dc_norm(j,i-1+nres)-scalp*
380 & dc_norm2(j,i-2+nres))/vbld(i-2+nres)
381 dtauangle(j,3,1,i)=-1/sing*dcostau(j,3,1,i)
382 dcostau(j,3,2,i)=fac1*dcosomicron(j,2,1,i-1)+fac2*
383 & dcosomicron(j,1,1,i)+fac3*dcosomicron(j,2,1,i-1)+fac4*
384 & dcosomicron(j,1,1,i)
385 dtauangle(j,3,2,i)=-1/sing*dcostau(j,3,2,i)
386 dcostau(j,3,3,i)=fac2*dcosomicron(j,1,2,i)+fac4*
387 & dcosomicron(j,1,2,i)-fac0*(dc_norm2(j,i-2+nres)-scalp*
388 & dc_norm(j,i-1+nres))/vbld(i-1+nres)
389 dtauangle(j,3,3,i)=-1/sing*dcostau(j,3,3,i)
390 c write(iout,*) "else",i
396 c Derivatives of side-chain angles alpha and omega
397 #if defined(MPI) && defined(PARINTDER)
398 do i=ibond_start,ibond_end
402 if(itype(i).ne.10 .and. itype(i).ne.ntyp1) then
403 fac5=1.0d0/dsqrt(2*(1+dcos(theta(i+1))))
407 fac9=fac5/vbld(i+nres)
408 scala1=scalar(dc_norm(1,i-1),dc_norm(1,i+nres))
409 scala2=scalar(dc_norm(1,i),dc_norm(1,i+nres))
410 cosa=dsqrt(0.5d0/(1.0d0+dcos(theta(i+1))))*(
411 & scalar(dC_norm(1,i),dC_norm(1,i+nres))
412 & -scalar(dC_norm(1,i-1),dC_norm(1,i+nres)))
413 sina=sqrt(1-cosa*cosa)
415 c write (iout,*) "i",i," cosa",cosa," sina",sina," sino",sino
417 dcosalpha(j,1,i)=fac6*(scala1*dc_norm(j,i-1)-
418 & dc_norm(j,i+nres))-cosa*fac7*dcostheta(j,1,i+1)
419 dalpha(j,1,i)=-1/sina*dcosalpha(j,1,i)
420 dcosalpha(j,2,i)=fac8*(dc_norm(j,i+nres)-
421 & scala2*dc_norm(j,i))-cosa*fac7*dcostheta(j,2,i+1)
422 dalpha(j,2,i)=-1/sina*dcosalpha(j,2,i)
423 dcosalpha(j,3,i)=(fac9*(dc_norm(j,i)-
424 & dc_norm(j,i-1))-(cosa*dc_norm(j,i+nres))/
426 dalpha(j,3,i)=-1/sina*dcosalpha(j,3,i)
428 c obtaining the derivatives of omega from sines
429 if(omeg(i).gt.-pi4.and.omeg(i).le.pi4.or.
430 & omeg(i).gt.pi34.and.omeg(i).le.pi.or.
431 & omeg(i).ge.-pi.and.omeg(i).le.-pi34) then
432 fac15=dcos(theta(i+1))/(dsin(theta(i+1))*
434 fac16=dcos(alph(i))/(dsin(alph(i))*dsin(alph(i)))
435 fac17=1.0d0/(dsin(theta(i+1))*dsin(alph(i)))
436 call vecpr(dc_norm(1,i+nres),dc_norm(1,i),vo1)
437 call vecpr(dc_norm(1,i+nres),dc_norm(1,i-1),vo2)
438 call vecpr(dc_norm(1,i),dc_norm(1,i-1),vo3)
439 coso_inv=1.0d0/dcos(omeg(i))
441 dsinomega(j,1,i)=sino*(fac15*dcostheta(j,1,i+1)
442 & +fac16*dcosalpha(j,1,i))-fac17/vbld(i)*vo1(j)-(
443 & sino*dc_norm(j,i-1))/vbld(i)
444 domega(j,1,i)=coso_inv*dsinomega(j,1,i)
445 dsinomega(j,2,i)=sino*(fac15*dcostheta(j,2,i+1)
446 & +fac16*dcosalpha(j,2,i))+fac17/vbld(i+1)*vo2(j)
447 & -sino*dc_norm(j,i)/vbld(i+1)
448 domega(j,2,i)=coso_inv*dsinomega(j,2,i)
449 dsinomega(j,3,i)=sino*fac16*dcosalpha(j,3,i)-
450 & fac17/vbld(i+nres)*vo3(j)-sino*dc_norm(j,i+nres)/
452 domega(j,3,i)=coso_inv*dsinomega(j,3,i)
455 c obtaining the derivatives of omega from cosines
456 fac10=sqrt(0.5d0*(1-dcos(theta(i+1))))
457 fac11=sqrt(0.5d0*(1+dcos(theta(i+1))))
462 dcosomega(j,1,i)=(-(0.25d0*cosa/fac11*
463 & dcostheta(j,1,i+1)+fac11*dcosalpha(j,1,i))*fac12+
464 & (0.25d0/fac10*sina*dcostheta(j,1,i+1)+cosa/sina*
465 & fac10*dcosalpha(j,1,i))*(scala2-fac11*cosa))/fac13
466 domega(j,1,i)=-1/sino*dcosomega(j,1,i)
467 dcosomega(j,2,i)=(((dc_norm(j,i+nres)-scala2*
468 & dc_norm(j,i))/vbld(i+1)-0.25d0*cosa/fac11*
469 & dcostheta(j,2,i+1)-fac11*dcosalpha(j,2,i))*fac12+
470 & (scala2-fac11*cosa)*(0.25d0*sina/fac10*
471 & dcostheta(j,2,i+1)+fac10*cosa/sina*dcosalpha(j,2,i)
473 domega(j,2,i)=-1/sino*dcosomega(j,2,i)
474 dcosomega(j,3,i)=1/fac10*((1/vbld(i+nres)*(dc_norm(j,i)-
475 & scala2*dc_norm(j,i+nres))-fac11*dcosalpha(j,3,i))*sina+
476 & (scala2-fac11*cosa)*(cosa/sina*dcosalpha(j,3,i)))/fac14
477 domega(j,3,i)=-1/sino*dcosomega(j,3,i)
490 #if defined(MPI) && defined(PARINTDER)
491 if (nfgtasks.gt.1) then
493 cd write (iout,*) "Gather dtheta"
495 write (iout,*) "dtheta before gather"
497 write (iout,'(i3,3(3f8.5,3x))') i,((dtheta(j,k,i),k=1,3),j=1,2)
500 call MPI_Gatherv(dtheta(1,1,ithet_start),ithet_count(fg_rank),
501 & MPI_THET,dtheta(1,1,1),ithet_count(0),ithet_displ(0),MPI_THET,
502 & king,FG_COMM,IERROR)
504 cd write (iout,*) "Gather dphi"
506 write (iout,*) "dphi before gather"
508 write (iout,'(i3,3(3f8.5,3x))') i,((dphi(j,k,i),k=1,3),j=1,3)
511 call MPI_Gatherv(dphi(1,1,iphi1_start),iphi1_count(fg_rank),
512 & MPI_GAM,dphi(1,1,1),iphi1_count(0),iphi1_displ(0),MPI_GAM,
513 & king,FG_COMM,IERROR)
514 cd write (iout,*) "Gather dalpha"
517 call MPI_Gatherv(dalpha(1,1,ibond_start),ibond_count(fg_rank),
518 & MPI_GAM,dalpha(1,1,1),ibond_count(0),ibond_displ(0),MPI_GAM,
519 & king,FG_COMM,IERROR)
520 cd write (iout,*) "Gather domega"
522 call MPI_Gatherv(domega(1,1,ibond_start),ibond_count(fg_rank),
523 & MPI_GAM,domega(1,1,1),ibond_count(0),ibond_displ(0),MPI_GAM,
524 & king,FG_COMM,IERROR)
529 write (iout,*) "dtheta after gather"
531 write (iout,'(i3,3(3f8.5,3x))') i,((dtheta(j,k,i),j=1,3),k=1,2)
533 write (iout,*) "dphi after gather"
535 write (iout,'(i3,3(3f8.5,3x))') i,((dphi(j,k,i),j=1,3),k=1,3)
537 write (iout,*) "dalpha after gather"
539 write (iout,'(i3,3(3f8.5,3x))') i,((dalpha(j,k,i),j=1,3),k=1,3)
541 write (iout,*) "domega after gather"
543 write (iout,'(i3,3(3f8.5,3x))') i,((domega(j,k,i),j=1,3),k=1,3)
545 write (iout,*) "dtauangle after gather"
547 write (iout,'(i3,3(3f8.5,3x))') i,
548 & (((dtauangle(j,k,l,i),j=1,3),k=1,3),l=1,3)
554 subroutine checkintcartgrad
555 implicit real*8 (a-h,o-z)
560 include 'COMMON.CHAIN'
563 include 'COMMON.INTERACT'
564 include 'COMMON.DERIV'
565 include 'COMMON.IOUNITS'
566 include 'COMMON.SETUP'
567 double precision dthetanum(3,2,maxres),dphinum(3,3,maxres)
568 & ,dalphanum(3,3,maxres), domeganum(3,3,maxres)
569 double precision theta_s(maxres),phi_s(maxres),alph_s(maxres),
570 & omeg_s(maxres),dc_norm_s(3)
571 double precision aincr /1.0d-5/
579 c Check theta gradient
581 & "Analytical (upper) and numerical (lower) gradient of theta"
588 call int_from_cart1(.false.)
589 dthetanum(j,1,i)=(theta(i)-theta_s(i))/aincr
592 dc(j,i-1)=dc(j,i-1)+aincr
594 dthetanum(j,2,i)=(theta(i)-theta_s(i))/aincr
597 write (iout,'(i5,3f10.5,5x,3f10.5)') i,(dtheta(j,1,i),j=1,3),
598 & (dtheta(j,2,i),j=1,3)
599 write (iout,'(5x,3f10.5,5x,3f10.5)') (dthetanum(j,1,i),j=1,3),
600 & (dthetanum(j,2,i),j=1,3)
601 write (iout,'(5x,3f10.5,5x,3f10.5)')
602 & (dthetanum(j,1,i)/dtheta(j,1,i),j=1,3),
603 & (dthetanum(j,2,i)/dtheta(j,2,i),j=1,3)
606 c Check gamma gradient
608 & "Analytical (upper) and numerical (lower) gradient of gamma"
614 dphinum(j,1,i)=(phi(i)-phi_s(i))/aincr
619 dphinum(j,2,i)=(phi(i)-phi_s(i))/aincr
622 dc(j,i-1)=dc(j,i-1)+aincr
624 dphinum(j,3,i)=(phi(i)-phi_s(i))/aincr
627 write (iout,'(i5,3(3f10.5,5x))') i,(dphi(j,1,i),j=1,3),
628 & (dphi(j,2,i),j=1,3),(dphi(j,3,i),j=1,3)
629 write (iout,'(5x,3(3f10.5,5x))') (dphinum(j,1,i),j=1,3),
630 & (dphinum(j,2,i),j=1,3),(dphinum(j,3,i),j=1,3)
631 write (iout,'(5x,3(3f10.5,5x))')
632 & (dphinum(j,1,i)/dphi(j,1,i),j=1,3),
633 & (dphinum(j,2,i)/dphi(j,2,i),j=1,3),
634 & (dphinum(j,3,i)/dphi(j,3,i),j=1,3)
637 c Check alpha gradient
639 & "Analytical (upper) and numerical (lower) gradient of alpha"
641 if(itype(i).ne.10) then
646 dalphanum(j,1,i)=(alph(i)-alph_s(i))
652 dalphanum(j,2,i)=(alph(i)-alph_s(i))
656 dc(j,i+nres)=dc(j,i+nres)+aincr
658 dalphanum(j,3,i)=(alph(i)-alph_s(i))
663 write (iout,'(i5,3(3f10.5,5x))') i,(dalpha(j,1,i),j=1,3),
664 & (dalpha(j,2,i),j=1,3),(dalpha(j,3,i),j=1,3)
665 write (iout,'(5x,3(3f10.5,5x))') (dalphanum(j,1,i),j=1,3),
666 & (dalphanum(j,2,i),j=1,3),(dalphanum(j,3,i),j=1,3)
667 write (iout,'(5x,3(3f10.5,5x))')
668 & (dalphanum(j,1,i)/dalpha(j,1,i),j=1,3),
669 & (dalphanum(j,2,i)/dalpha(j,2,i),j=1,3),
670 & (dalphanum(j,3,i)/dalpha(j,3,i),j=1,3)
673 c Check omega gradient
675 & "Analytical (upper) and numerical (lower) gradient of omega"
677 if(itype(i).ne.10) then
682 domeganum(j,1,i)=(omeg(i)-omeg_s(i))
688 domeganum(j,2,i)=(omeg(i)-omeg_s(i))
692 dc(j,i+nres)=dc(j,i+nres)+aincr
694 domeganum(j,3,i)=(omeg(i)-omeg_s(i))
699 write (iout,'(i5,3(3f10.5,5x))') i,(domega(j,1,i),j=1,3),
700 & (domega(j,2,i),j=1,3),(domega(j,3,i),j=1,3)
701 write (iout,'(5x,3(3f10.5,5x))') (domeganum(j,1,i),j=1,3),
702 & (domeganum(j,2,i),j=1,3),(domeganum(j,3,i),j=1,3)
703 write (iout,'(5x,3(3f10.5,5x))')
704 & (domeganum(j,1,i)/domega(j,1,i),j=1,3),
705 & (domeganum(j,2,i)/domega(j,2,i),j=1,3),
706 & (domeganum(j,3,i)/domega(j,3,i),j=1,3)
711 c------------------------------------------------------------
712 subroutine chainbuild_cart
713 implicit real*8 (a-h,o-z)
718 include 'COMMON.SETUP'
719 include 'COMMON.CHAIN'
720 include 'COMMON.LOCAL'
721 include 'COMMON.TIME1'
722 include 'COMMON.IOUNITS'
725 if (nfgtasks.gt.1) then
726 c write (iout,*) "BCAST in chainbuild_cart"
728 c Broadcast the order to build the chain and compute internal coordinates
729 c to the slaves. The slaves receive the order in ERGASTULUM.
731 c write (iout,*) "CHAINBUILD_CART: DC before BCAST"
733 c write (iout,'(i3,3f10.5,5x,3f10.5)') i,(dc(j,i),j=1,3),
734 c & (dc(j,i+nres),j=1,3)
737 & call MPI_Bcast(7,1,MPI_INTEGER,king,FG_COMM,IERROR)
738 time_bcast7=time_bcast7+MPI_Wtime()-time00
740 call MPI_Bcast(dc(1,0),6*(nres+1),MPI_DOUBLE_PRECISION,
742 c write (iout,*) "CHAINBUILD_CART: DC after BCAST"
744 c write (iout,'(i3,3f10.5,5x,3f10.5)') i,(dc(j,i),j=1,3),
745 c & (dc(j,i+nres),j=1,3)
747 c write (iout,*) "End BCAST in chainbuild_cart"
749 time_bcast=time_bcast+MPI_Wtime()-time00
750 time_bcastc=time_bcastc+MPI_Wtime()-time01
759 c(j,i)=c(j,i-1)+dc(j,i-1)
764 c(j,i+nres)=c(j,i)+dc(j,i+nres)
768 c write (iout,*) "CHAINBUILD_CART"
770 call int_from_cart1(.false.)