X-Git-Url: http://mmka.chem.univ.gda.pl/gitweb/?a=blobdiff_plain;f=source%2Fwham%2Fsrc-NEWSC-NEWCORR%2Fodlodc.f;fp=source%2Fwham%2Fsrc-NEWSC-NEWCORR%2Fodlodc.f;h=c18ac72ace8dd62830468aef58f6ea7cc9dcbe8a;hb=d101c97dea752458d76055fdbae49c26fff03c1f;hp=0000000000000000000000000000000000000000;hpb=325eda160c9ad2982501e091ca40606a29043712;p=unres.git

diff --git a/source/wham/src-NEWSC-NEWCORR/odlodc.f b/source/wham/src-NEWSC-NEWCORR/odlodc.f
new file mode 100644
index 0000000..c18ac72
--- /dev/null
+++ b/source/wham/src-NEWSC-NEWCORR/odlodc.f
@@ -0,0 +1,55 @@
+      subroutine odlodc(r1,r2,a,b,uu,vv,aa,bb,dd)
+      implicit real*8 (a-h,o-z)
+      dimension r1(3),r2(3),a(3),b(3),x(3),y(3)
+      odl(u,v) = (r1(1)-r2(1))**2+(r1(2)-r2(2))**2+(r1(3)-r2(3))**2
+     & + 2*ar*u - 2*br*v - 2*ab*u*v + aa*u**2 + bb*v**2
+c      print *,"r1",(r1(i),i=1,3)
+c      print *,"r2",(r2(i),i=1,3)
+c      print *,"a",(a(i),i=1,3)
+c      print *,"b",(b(i),i=1,3)
+      aa = a(1)**2+a(2)**2+a(3)**2
+      bb = b(1)**2+b(2)**2+b(3)**2
+      ab = a(1)*b(1)+a(2)*b(2)+a(3)*b(3) 
+      ar = a(1)*(r1(1)-r2(1))+a(2)*(r1(2)-r2(2))+a(3)*(r1(3)-r2(3))
+      br = b(1)*(r1(1)-r2(1))+b(2)*(r1(2)-r2(2))+b(3)*(r1(3)-r2(3))
+      det = aa*bb-ab**2
+c      print *,'aa',aa,' bb',bb,' ab',ab,' ar',ar,' br',br,' det',det
+      uu = (-ar*bb+br*ab)/det
+      vv = (br*aa-ar*ab)/det
+c      print *,u,v
+      uu=dmin1(uu,1.0d0)
+      uu=dmax1(uu,0.0d0)
+      vv=dmin1(vv,1.0d0)
+      vv=dmax1(vv,0.0d0)
+      dd1 = odl(uu,vv)
+      dd2 = odl(0.0d0,0.0d0)
+      dd3 = odl(0.0d0,1.0d0)
+      dd4 = odl(1.0d0,0.0d0)
+      dd5 = odl(1.0d0,1.0d0)
+      dd = dsqrt(dmin1(dd1,dd2,dd3,dd4,dd5))
+      if (dd.eq.dd2) then
+        uu=0.0d0
+        vv=0.0d0
+      else if (dd.eq.dd3) then
+        uu=0.0d0
+        vv=1.0d0
+      else if (dd.eq.dd4) then
+        uu=1.0d0
+        vv=0.0d0
+      else if (dd.eq.dd5) then
+        uu=1.0d0
+        vv=1.0d0
+      endif 
+c Control check
+c      do i=1,3
+c        x(i)=r1(i)+u*a(i)
+c        y(i)=r2(i)+v*b(i)
+c      enddo
+c      dd1 = (x(1)-y(1))**2+(x(2)-y(2))**2+(x(3)-y(3))**2
+c      dd1 = dsqrt(dd1)
+      aa = dsqrt(aa)
+      bb = dsqrt(bb)
+c      write (8,*) uu,vv,dd,dd1
+c      print *,dd,dd1
+      return
+      end