      subroutine efidp(m,n,icontr,inddyy,esp,efix,efiy,efiz,ccx,ccy,ccz,
     &                 g,density,zeta,axbx,ayby,azbz,px,py,pz,pxax,pyay,
     &                 pzaz,base,ss0,ss1,ss2,ss3,ps0,ps1,ps2,dscond,ds0,
     &                 ds1,fscond,fs0,dpcond,c,c1,pxcx,pycy,pzcz,u,
     &                 twozeta,ss4,cx,cy,cz,cx1,cy1,cz1,saxs0,says0,
     &                 sazs0,saxs1,says1,sazs1,saxs2,says2,sazs2,saxs3,
     &                 says3,sazs3,paxs0,pays0,pazs0,paxs1,pays1,pazs1,
     &                 paxs2,pays2,pazs2,daxs0,days0,dazs0,daxs1,days1,
     &                 dazs1,faxs0,fays0,fazs0,daxscond,dayscond,
     &                 dazscond,daxpcond,daypcond,dazpcond,faxscond,
     &                 fayscond,fazscond)

c  this subroutine written by roger gallant and alain st-amant of the
c  department of chemistry, university of ottawa, ottawa, canada.
c  all rights reserved.  this is part of the DeFT project.

      implicit real*8(a-h,o-z)

      dimension icontr(*),inddyy(m,*)

      dimension g(*),density(*),zeta(*),axbx(*),ayby(*),azbz(*),px(*),
     &          py(*),pz(*),pxax(*),pyay(*),pzaz(*),base(*),ss0(*),
     &          ss1(*),ss2(*),ss3(*),ss4(*),ps0(n,*),ps1(n,*),ps2(n,*),
     &          dscond(m,*),ds0(n,*),ds1(n,*),fscond(m,*),fs0(n,*),
     &          dpcond(m,*),c(*),c1(n,*),pxcx(*),pycy(*),pzcz(*),u(*),
     &          twozeta(*),cx(*),cy(*),cz(*),cx1(n,*),cy1(n,*),cz1(n,*),
     &          saxs0(*),says0(*),sazs0(*),saxs1(*),says1(*),sazs1(*),
     &          saxs2(*),says2(*),sazs2(*),saxs3(*),says3(*),sazs3(*),
     &          paxs0(n,*),pays0(n,*),pazs0(n,*),paxs1(n,*),pays1(n,*),
     &          pazs1(n,*),paxs2(n,*),pays2(n,*),pazs2(n,*),daxs0(n,*),
     &          days0(n,*),dazs0(n,*),daxs1(n,*),days1(n,*),dazs1(n,*),
     &          faxs0(n,*),fays0(n,*),fazs0(n,*),daxscond(m,*),
     &          dayscond(m,*),dazscond(m,*),faxscond(m,*),fayscond(m,*),
     &          fazscond(m,*),daxpcond(m,*),daypcond(m,*),dazpcond(m,*)

      data zero,two,three,pi/0.0d0,2.0d0,3.0d0,3.1415926535898d0/

      do 1001 i=1,n
 1001 twozeta(i)=two*zeta(i)

      twopi=two*pi

      do 1002 i=1,n
      pxcx(i)=px(i)-ccx
      pycy(i)=py(i)-ccy
      pzcz(i)=pz(i)-ccz
 1002 u(i)=zeta(i)*(pxcx(i)*pxcx(i)+pycy(i)*pycy(i)+pzcz(i)*pzcz(i))

      do 1003 i=1,n
      call augg(5,u(i),g)
      ss0(i)=base(i)*g(1)*twopi/zeta(i)
      ss1(i)=base(i)*g(2)*twopi/zeta(i)
      ss2(i)=base(i)*g(3)*twopi/zeta(i)
      ss3(i)=base(i)*g(4)*twopi/zeta(i)
 1003 ss4(i)=base(i)*g(5)*twopi/zeta(i)

      do 1004 i=1,n

      ps0(i,1)=pxax(i)*ss0(i)-pxcx(i)*ss1(i)
      ps0(i,2)=pyay(i)*ss0(i)-pycy(i)*ss1(i)
      ps0(i,3)=pzaz(i)*ss0(i)-pzcz(i)*ss1(i)

      ps1(i,1)=pxax(i)*ss1(i)-pxcx(i)*ss2(i)
      ps1(i,2)=pyay(i)*ss1(i)-pycy(i)*ss2(i)
      ps1(i,3)=pzaz(i)*ss1(i)-pzcz(i)*ss2(i)

      ps2(i,1)=pxax(i)*ss2(i)-pxcx(i)*ss3(i)
      ps2(i,2)=pyay(i)*ss2(i)-pycy(i)*ss3(i)
      ps2(i,3)=pzaz(i)*ss2(i)-pzcz(i)*ss3(i)

 1004 continue

      do 1005 i=1,n
 1005 c(i)=(ss0(i)-ss1(i))/twozeta(i)

      do 1006 i=1,n
      ds0(i,1)=pxax(i)*ps0(i,1)-pxcx(i)*ps1(i,1)+c(i)
      ds0(i,2)=pyay(i)*ps0(i,1)-pycy(i)*ps1(i,1)
      ds0(i,3)=pzaz(i)*ps0(i,1)-pzcz(i)*ps1(i,1)
      ds0(i,4)=pyay(i)*ps0(i,2)-pycy(i)*ps1(i,2)+c(i)
      ds0(i,5)=pzaz(i)*ps0(i,2)-pzcz(i)*ps1(i,2)
 1006 ds0(i,6)=pzaz(i)*ps0(i,3)-pzcz(i)*ps1(i,3)+c(i)

      do 1007 i=1,n
 1007 c(i)=(ss1(i)-ss2(i))/twozeta(i)

      do 1008 i=1,n
      ds1(i,1)=pxax(i)*ps1(i,1)-pxcx(i)*ps2(i,1)+c(i)
      ds1(i,2)=pyay(i)*ps1(i,1)-pycy(i)*ps2(i,1)
      ds1(i,3)=pzaz(i)*ps1(i,1)-pzcz(i)*ps2(i,1)
      ds1(i,4)=pyay(i)*ps1(i,2)-pycy(i)*ps2(i,2)+c(i)
      ds1(i,5)=pzaz(i)*ps1(i,2)-pzcz(i)*ps2(i,2)
 1008 ds1(i,6)=pzaz(i)*ps1(i,3)-pzcz(i)*ps2(i,3)+c(i)

      do 1009 j=1,3
      do 1009 i=1,n
 1009 c1(i,j)=(ps0(i,j)-ps1(i,j))/twozeta(i)

      do 1010 i=1,n
      fs0(i,01)=pxax(i)*ds0(i,1)-pxcx(i)*ds1(i,1)+c1(i,1)*two
      fs0(i,02)=pyay(i)*ds0(i,1)-pycy(i)*ds1(i,1)
      fs0(i,03)=pzaz(i)*ds0(i,1)-pzcz(i)*ds1(i,1)
      fs0(i,04)=pyay(i)*ds0(i,2)-pycy(i)*ds1(i,2)+c1(i,1)
      fs0(i,05)=pzaz(i)*ds0(i,2)-pzcz(i)*ds1(i,2)
      fs0(i,06)=pzaz(i)*ds0(i,3)-pzcz(i)*ds1(i,3)+c1(i,1)
      fs0(i,07)=pyay(i)*ds0(i,4)-pycy(i)*ds1(i,4)+c1(i,2)*two
      fs0(i,08)=pzaz(i)*ds0(i,4)-pzcz(i)*ds1(i,4)
      fs0(i,09)=pzaz(i)*ds0(i,5)-pzcz(i)*ds1(i,5)+c1(i,2)
 1010 fs0(i,10)=pzaz(i)*ds0(i,6)-pzcz(i)*ds1(i,6)+c1(i,3)*two

      do 1011 i=1,n

      saxs0(i)=twozeta(i)*pxcx(i)*ss1(i)
      says0(i)=twozeta(i)*pycy(i)*ss1(i)
      sazs0(i)=twozeta(i)*pzcz(i)*ss1(i)

      saxs1(i)=twozeta(i)*pxcx(i)*ss2(i)
      says1(i)=twozeta(i)*pycy(i)*ss2(i)
      sazs1(i)=twozeta(i)*pzcz(i)*ss2(i)

      saxs2(i)=twozeta(i)*pxcx(i)*ss3(i)
      says2(i)=twozeta(i)*pycy(i)*ss3(i)
      sazs2(i)=twozeta(i)*pzcz(i)*ss3(i)

      saxs3(i)=twozeta(i)*pxcx(i)*ss4(i)
      says3(i)=twozeta(i)*pycy(i)*ss4(i)
      sazs3(i)=twozeta(i)*pzcz(i)*ss4(i)

 1011 continue

      do 1012 i=1,n

      paxs0(i,1)=pxax(i)*saxs0(i)-pxcx(i)*saxs1(i)+ss1(i)
      paxs0(i,2)=pyay(i)*saxs0(i)-pycy(i)*saxs1(i)
      paxs0(i,3)=pzaz(i)*saxs0(i)-pzcz(i)*saxs1(i)
      pays0(i,1)=pxax(i)*says0(i)-pxcx(i)*says1(i)
      pays0(i,2)=pyay(i)*says0(i)-pycy(i)*says1(i)+ss1(i)
      pays0(i,3)=pzaz(i)*says0(i)-pzcz(i)*says1(i)
      pazs0(i,1)=pxax(i)*sazs0(i)-pxcx(i)*sazs1(i)
      pazs0(i,2)=pyay(i)*sazs0(i)-pycy(i)*sazs1(i)
      pazs0(i,3)=pzaz(i)*sazs0(i)-pzcz(i)*sazs1(i)+ss1(i)

      paxs1(i,1)=pxax(i)*saxs1(i)-pxcx(i)*saxs2(i)+ss2(i)
      paxs1(i,2)=pyay(i)*saxs1(i)-pycy(i)*saxs2(i)
      paxs1(i,3)=pzaz(i)*saxs1(i)-pzcz(i)*saxs2(i)
      pays1(i,1)=pxax(i)*says1(i)-pxcx(i)*says2(i)
      pays1(i,2)=pyay(i)*says1(i)-pycy(i)*says2(i)+ss2(i)
      pays1(i,3)=pzaz(i)*says1(i)-pzcz(i)*says2(i)
      pazs1(i,1)=pxax(i)*sazs1(i)-pxcx(i)*sazs2(i)
      pazs1(i,2)=pyay(i)*sazs1(i)-pycy(i)*sazs2(i)
      pazs1(i,3)=pzaz(i)*sazs1(i)-pzcz(i)*sazs2(i)+ss2(i)

      paxs2(i,1)=pxax(i)*saxs2(i)-pxcx(i)*saxs3(i)+ss3(i)
      paxs2(i,2)=pyay(i)*saxs2(i)-pycy(i)*saxs3(i)
      paxs2(i,3)=pzaz(i)*saxs2(i)-pzcz(i)*saxs3(i)
      pays2(i,1)=pxax(i)*says2(i)-pxcx(i)*says3(i)
      pays2(i,2)=pyay(i)*says2(i)-pycy(i)*says3(i)+ss3(i)
      pays2(i,3)=pzaz(i)*says2(i)-pzcz(i)*says3(i)
      pazs2(i,1)=pxax(i)*sazs2(i)-pxcx(i)*sazs3(i)
      pazs2(i,2)=pyay(i)*sazs2(i)-pycy(i)*sazs3(i)
      pazs2(i,3)=pzaz(i)*sazs2(i)-pzcz(i)*sazs3(i)+ss3(i)

 1012 continue

      do 1013 i=1,n
      cx(i)=(saxs0(i)-saxs1(i))/twozeta(i)
      cy(i)=(says0(i)-says1(i))/twozeta(i)
 1013 cz(i)=(sazs0(i)-sazs1(i))/twozeta(i)

      do 1014 i=1,n
      daxs0(i,1)=pxax(i)*paxs0(i,1)-pxcx(i)*paxs1(i,1)+cx(i)+ps1(i,1)
      daxs0(i,2)=pyay(i)*paxs0(i,1)-pycy(i)*paxs1(i,1)
      daxs0(i,3)=pzaz(i)*paxs0(i,1)-pzcz(i)*paxs1(i,1)
      daxs0(i,4)=pyay(i)*paxs0(i,2)-pycy(i)*paxs1(i,2)+cx(i)
      daxs0(i,5)=pzaz(i)*paxs0(i,2)-pzcz(i)*paxs1(i,2)
      daxs0(i,6)=pzaz(i)*paxs0(i,3)-pzcz(i)*paxs1(i,3)+cx(i)
      days0(i,1)=pxax(i)*pays0(i,1)-pxcx(i)*pays1(i,1)+cy(i)
      days0(i,2)=pyay(i)*pays0(i,1)-pycy(i)*pays1(i,1)      +ps1(i,1)
      days0(i,3)=pzaz(i)*pays0(i,1)-pzcz(i)*pays1(i,1)
      days0(i,4)=pyay(i)*pays0(i,2)-pycy(i)*pays1(i,2)+cy(i)+ps1(i,2)
      days0(i,5)=pzaz(i)*pays0(i,2)-pzcz(i)*pays1(i,2)
      days0(i,6)=pzaz(i)*pays0(i,3)-pzcz(i)*pays1(i,3)+cy(i)
      dazs0(i,1)=pxax(i)*pazs0(i,1)-pxcx(i)*pazs1(i,1)+cz(i)
      dazs0(i,2)=pyay(i)*pazs0(i,1)-pycy(i)*pazs1(i,1)
      dazs0(i,3)=pzaz(i)*pazs0(i,1)-pzcz(i)*pazs1(i,1)      +ps1(i,1)
      dazs0(i,4)=pyay(i)*pazs0(i,2)-pycy(i)*pazs1(i,2)+cz(i)
      dazs0(i,5)=pzaz(i)*pazs0(i,2)-pzcz(i)*pazs1(i,2)      +ps1(i,2)
 1014 dazs0(i,6)=pzaz(i)*pazs0(i,3)-pzcz(i)*pazs1(i,3)+cz(i)+ps1(i,3)

      do 1015 i=1,n
      cx(i)=(saxs1(i)-saxs2(i))/twozeta(i)
      cy(i)=(says1(i)-says2(i))/twozeta(i)
 1015 cz(i)=(sazs1(i)-sazs2(i))/twozeta(i)

      do 1016 i=1,n
      daxs1(i,1)=pxax(i)*paxs1(i,1)-pxcx(i)*paxs2(i,1)+cx(i)+ps2(i,1)
      daxs1(i,2)=pyay(i)*paxs1(i,1)-pycy(i)*paxs2(i,1)
      daxs1(i,3)=pzaz(i)*paxs1(i,1)-pzcz(i)*paxs2(i,1)
      daxs1(i,4)=pyay(i)*paxs1(i,2)-pycy(i)*paxs2(i,2)+cx(i)
      daxs1(i,5)=pzaz(i)*paxs1(i,2)-pzcz(i)*paxs2(i,2)
      daxs1(i,6)=pzaz(i)*paxs1(i,3)-pzcz(i)*paxs2(i,3)+cx(i)
      days1(i,1)=pxax(i)*pays1(i,1)-pxcx(i)*pays2(i,1)+cy(i)
      days1(i,2)=pyay(i)*pays1(i,1)-pycy(i)*pays2(i,1)      +ps2(i,1)
      days1(i,3)=pzaz(i)*pays1(i,1)-pzcz(i)*pays2(i,1)
      days1(i,4)=pyay(i)*pays1(i,2)-pycy(i)*pays2(i,2)+cy(i)+ps2(i,2)
      days1(i,5)=pzaz(i)*pays1(i,2)-pzcz(i)*pays2(i,2)
      days1(i,6)=pzaz(i)*pays1(i,3)-pzcz(i)*pays2(i,3)+cy(i)
      dazs1(i,1)=pxax(i)*pazs1(i,1)-pxcx(i)*pazs2(i,1)+cz(i)
      dazs1(i,2)=pyay(i)*pazs1(i,1)-pycy(i)*pazs2(i,1)
      dazs1(i,3)=pzaz(i)*pazs1(i,1)-pzcz(i)*pazs2(i,1)      +ps2(i,1)
      dazs1(i,4)=pyay(i)*pazs1(i,2)-pycy(i)*pazs2(i,2)+cz(i)
      dazs1(i,5)=pzaz(i)*pazs1(i,2)-pzcz(i)*pazs2(i,2)      +ps2(i,2)
 1016 dazs1(i,6)=pzaz(i)*pazs1(i,3)-pzcz(i)*pazs2(i,3)+cz(i)+ps2(i,3)
   
      do 1017 j=1,3
      do 1017 i=1,n
      cx1(i,j)=(paxs0(i,j)-paxs1(i,j))/twozeta(i)
      cy1(i,j)=(pays0(i,j)-pays1(i,j))/twozeta(i)
 1017 cz1(i,j)=(pazs0(i,j)-pazs1(i,j))/twozeta(i)

      do 1018 i=1,n
      faxs0(i,01)=pxax(i)*daxs0(i,1)-pxcx(i)*daxs1(i,1)+cx1(i,1)*two
     &                                                 +ds1(i,1)
      faxs0(i,02)=pyay(i)*daxs0(i,1)-pycy(i)*daxs1(i,1)
      faxs0(i,03)=pzaz(i)*daxs0(i,1)-pzcz(i)*daxs1(i,1)
      faxs0(i,04)=pyay(i)*daxs0(i,2)-pycy(i)*daxs1(i,2)+cx1(i,1)
      faxs0(i,05)=pzaz(i)*daxs0(i,2)-pzcz(i)*daxs1(i,2)
      faxs0(i,06)=pzaz(i)*daxs0(i,3)-pzcz(i)*daxs1(i,3)+cx1(i,1)
      faxs0(i,07)=pyay(i)*daxs0(i,4)-pycy(i)*daxs1(i,4)+cx1(i,2)*two
      faxs0(i,08)=pzaz(i)*daxs0(i,4)-pzcz(i)*daxs1(i,4)
      faxs0(i,09)=pzaz(i)*daxs0(i,5)-pzcz(i)*daxs1(i,5)+cx1(i,2)
      faxs0(i,10)=pzaz(i)*daxs0(i,6)-pzcz(i)*daxs1(i,6)+cx1(i,3)*two
      fays0(i,01)=pxax(i)*days0(i,1)-pxcx(i)*days1(i,1)+cy1(i,1)*two
      fays0(i,02)=pyay(i)*days0(i,1)-pycy(i)*days1(i,1)
     &                                                 +ds1(i,1)
      fays0(i,03)=pzaz(i)*days0(i,1)-pzcz(i)*days1(i,1)
      fays0(i,04)=pyay(i)*days0(i,2)-pycy(i)*days1(i,2)+cy1(i,1)
     &                                                 +ds1(i,2)
      fays0(i,05)=pzaz(i)*days0(i,2)-pzcz(i)*days1(i,2)
      fays0(i,06)=pzaz(i)*days0(i,3)-pzcz(i)*days1(i,3)+cy1(i,1)
      fays0(i,07)=pyay(i)*days0(i,4)-pycy(i)*days1(i,4)+cy1(i,2)*two
     &                                                 +ds1(i,4)
      fays0(i,08)=pzaz(i)*days0(i,4)-pzcz(i)*days1(i,4)
      fays0(i,09)=pzaz(i)*days0(i,5)-pzcz(i)*days1(i,5)+cy1(i,2)
      fays0(i,10)=pzaz(i)*days0(i,6)-pzcz(i)*days1(i,6)+cy1(i,3)*two
      fazs0(i,01)=pxax(i)*dazs0(i,1)-pxcx(i)*dazs1(i,1)+cz1(i,1)*two
      fazs0(i,02)=pyay(i)*dazs0(i,1)-pycy(i)*dazs1(i,1)
      fazs0(i,03)=pzaz(i)*dazs0(i,1)-pzcz(i)*dazs1(i,1)
     &                                                 +ds1(i,1)
      fazs0(i,04)=pyay(i)*dazs0(i,2)-pycy(i)*dazs1(i,2)+cz1(i,1)
      fazs0(i,05)=pzaz(i)*dazs0(i,2)-pzcz(i)*dazs1(i,2)
     &                                                 +ds1(i,2)
      fazs0(i,06)=pzaz(i)*dazs0(i,3)-pzcz(i)*dazs1(i,3)+cz1(i,1)    
     &                                                 +ds1(i,3)
      fazs0(i,07)=pyay(i)*dazs0(i,4)-pycy(i)*dazs1(i,4)+cz1(i,2)*two
      fazs0(i,08)=pzaz(i)*dazs0(i,4)-pzcz(i)*dazs1(i,4)
     &                                                 +ds1(i,4)
      fazs0(i,09)=pzaz(i)*dazs0(i,5)-pzcz(i)*dazs1(i,5)+cz1(i,2)
     &                                                 +ds1(i,5)
 1018 fazs0(i,10)=pzaz(i)*dazs0(i,6)-pzcz(i)*dazs1(i,6)+cz1(i,3)*two
     &                                                 +ds1(i,6)

      do 1019 j=1,6
      do 1019 i=1,m
 1019 dscond(i,j)=zero

      do 1020 j=1,10
      do 1020 i=1,m
 1020 fscond(i,j)=zero

      do 1021 j=1,6
      do 1021 i=1,n
 1021 dscond(icontr(i),j)=dscond(icontr(i),j)+ds0(i,j)

      do 1022 j=1,10
      do 1022 i=1,n
 1022 fscond(icontr(i),j)=fscond(icontr(i),j)+fs0(i,j)

      do 1023 i=1,m
      dpcond(i,01)=axbx(i)*dscond(i,1)+fscond(i,01)
      dpcond(i,02)=ayby(i)*dscond(i,1)+fscond(i,02)
      dpcond(i,03)=azbz(i)*dscond(i,1)+fscond(i,03)
      dpcond(i,04)=axbx(i)*dscond(i,2)+fscond(i,02)
      dpcond(i,05)=ayby(i)*dscond(i,2)+fscond(i,04)
      dpcond(i,06)=azbz(i)*dscond(i,2)+fscond(i,05)
      dpcond(i,07)=axbx(i)*dscond(i,3)+fscond(i,03)
      dpcond(i,08)=ayby(i)*dscond(i,3)+fscond(i,05)
      dpcond(i,09)=azbz(i)*dscond(i,3)+fscond(i,06)
      dpcond(i,10)=axbx(i)*dscond(i,4)+fscond(i,04)
      dpcond(i,11)=ayby(i)*dscond(i,4)+fscond(i,07)
      dpcond(i,12)=azbz(i)*dscond(i,4)+fscond(i,08)
      dpcond(i,13)=axbx(i)*dscond(i,5)+fscond(i,05)
      dpcond(i,14)=ayby(i)*dscond(i,5)+fscond(i,08)
      dpcond(i,15)=azbz(i)*dscond(i,5)+fscond(i,09)
      dpcond(i,16)=axbx(i)*dscond(i,6)+fscond(i,06)
      dpcond(i,17)=ayby(i)*dscond(i,6)+fscond(i,09)
 1023 dpcond(i,18)=azbz(i)*dscond(i,6)+fscond(i,10)

      do 1024 j=1,6
      do 1024 i=1,m
      daxscond(i,j)=zero
      dayscond(i,j)=zero
 1024 dazscond(i,j)=zero

      do 1025 j=1,10
      do 1025 i=1,m
      faxscond(i,j)=zero
      fayscond(i,j)=zero
 1025 fazscond(i,j)=zero

      do 1026 j=1,6
      do 1026 i=1,n
      daxscond(icontr(i),j)=daxscond(icontr(i),j)+daxs0(i,j)
      dayscond(icontr(i),j)=dayscond(icontr(i),j)+days0(i,j)
 1026 dazscond(icontr(i),j)=dazscond(icontr(i),j)+dazs0(i,j)

      do 1027 j=1,10
      do 1027 i=1,n
      faxscond(icontr(i),j)=faxscond(icontr(i),j)+faxs0(i,j)
      fayscond(icontr(i),j)=fayscond(icontr(i),j)+fays0(i,j)
 1027 fazscond(icontr(i),j)=fazscond(icontr(i),j)+fazs0(i,j)

      do 1028 i=1,m
      daxpcond(i,01)=axbx(i)*daxscond(i,1)+faxscond(i,01)
      daxpcond(i,02)=ayby(i)*daxscond(i,1)+faxscond(i,02)
      daxpcond(i,03)=azbz(i)*daxscond(i,1)+faxscond(i,03)
      daxpcond(i,04)=axbx(i)*daxscond(i,2)+faxscond(i,02)
      daxpcond(i,05)=ayby(i)*daxscond(i,2)+faxscond(i,04)
      daxpcond(i,06)=azbz(i)*daxscond(i,2)+faxscond(i,05)
      daxpcond(i,07)=axbx(i)*daxscond(i,3)+faxscond(i,03)
      daxpcond(i,08)=ayby(i)*daxscond(i,3)+faxscond(i,05)
      daxpcond(i,09)=azbz(i)*daxscond(i,3)+faxscond(i,06)
      daxpcond(i,10)=axbx(i)*daxscond(i,4)+faxscond(i,04)
      daxpcond(i,11)=ayby(i)*daxscond(i,4)+faxscond(i,07)
      daxpcond(i,12)=azbz(i)*daxscond(i,4)+faxscond(i,08)
      daxpcond(i,13)=axbx(i)*daxscond(i,5)+faxscond(i,05)
      daxpcond(i,14)=ayby(i)*daxscond(i,5)+faxscond(i,08)
      daxpcond(i,15)=azbz(i)*daxscond(i,5)+faxscond(i,09)
      daxpcond(i,16)=axbx(i)*daxscond(i,6)+faxscond(i,06)
      daxpcond(i,17)=ayby(i)*daxscond(i,6)+faxscond(i,09)
      daxpcond(i,18)=azbz(i)*daxscond(i,6)+faxscond(i,10)
      daypcond(i,01)=axbx(i)*dayscond(i,1)+fayscond(i,01)
      daypcond(i,02)=ayby(i)*dayscond(i,1)+fayscond(i,02)
      daypcond(i,03)=azbz(i)*dayscond(i,1)+fayscond(i,03)
      daypcond(i,04)=axbx(i)*dayscond(i,2)+fayscond(i,02)
      daypcond(i,05)=ayby(i)*dayscond(i,2)+fayscond(i,04)
      daypcond(i,06)=azbz(i)*dayscond(i,2)+fayscond(i,05)
      daypcond(i,07)=axbx(i)*dayscond(i,3)+fayscond(i,03)
      daypcond(i,08)=ayby(i)*dayscond(i,3)+fayscond(i,05)
      daypcond(i,09)=azbz(i)*dayscond(i,3)+fayscond(i,06)
      daypcond(i,10)=axbx(i)*dayscond(i,4)+fayscond(i,04)
      daypcond(i,11)=ayby(i)*dayscond(i,4)+fayscond(i,07)
      daypcond(i,12)=azbz(i)*dayscond(i,4)+fayscond(i,08)
      daypcond(i,13)=axbx(i)*dayscond(i,5)+fayscond(i,05)
      daypcond(i,14)=ayby(i)*dayscond(i,5)+fayscond(i,08)
      daypcond(i,15)=azbz(i)*dayscond(i,5)+fayscond(i,09)
      daypcond(i,16)=axbx(i)*dayscond(i,6)+fayscond(i,06)
      daypcond(i,17)=ayby(i)*dayscond(i,6)+fayscond(i,09)
      daypcond(i,18)=azbz(i)*dayscond(i,6)+fayscond(i,10)
      dazpcond(i,01)=axbx(i)*dazscond(i,1)+fazscond(i,01)
      dazpcond(i,02)=ayby(i)*dazscond(i,1)+fazscond(i,02)
      dazpcond(i,03)=azbz(i)*dazscond(i,1)+fazscond(i,03)
      dazpcond(i,04)=axbx(i)*dazscond(i,2)+fazscond(i,02)
      dazpcond(i,05)=ayby(i)*dazscond(i,2)+fazscond(i,04)
      dazpcond(i,06)=azbz(i)*dazscond(i,2)+fazscond(i,05)
      dazpcond(i,07)=axbx(i)*dazscond(i,3)+fazscond(i,03)
      dazpcond(i,08)=ayby(i)*dazscond(i,3)+fazscond(i,05)
      dazpcond(i,09)=azbz(i)*dazscond(i,3)+fazscond(i,06)
      dazpcond(i,10)=axbx(i)*dazscond(i,4)+fazscond(i,04)
      dazpcond(i,11)=ayby(i)*dazscond(i,4)+fazscond(i,07)
      dazpcond(i,12)=azbz(i)*dazscond(i,4)+fazscond(i,08)
      dazpcond(i,13)=axbx(i)*dazscond(i,5)+fazscond(i,05)
      dazpcond(i,14)=ayby(i)*dazscond(i,5)+fazscond(i,08)
      dazpcond(i,15)=azbz(i)*dazscond(i,5)+fazscond(i,09)
      dazpcond(i,16)=axbx(i)*dazscond(i,6)+fazscond(i,06)
      dazpcond(i,17)=ayby(i)*dazscond(i,6)+fazscond(i,09)
 1028 dazpcond(i,18)=azbz(i)*dazscond(i,6)+fazscond(i,10)

      do 1029 j=1,18
      do 1029 i=1,m
 1029 esp=esp-density(inddyy(i,j))*dpcond(i,j)

      do 1030 j=1,18
      do 1030 i=1,m
      efix=efix-density(inddyy(i,j))*daxpcond(i,j)
      efiy=efiy-density(inddyy(i,j))*daypcond(i,j)
 1030 efiz=efiz-density(inddyy(i,j))*dazpcond(i,j)

      return 
      end 
