1c 2c $Id: hnd_dervxyz.F 19696 2010-10-29 16:53:42Z d3y133 $ 3c 4c Taken and modified from HONDO 5c 6 subroutine dim_dervxyz(mder) 7 implicit none 8c 9c Routine can do ddvxyz and dvxyz by setting mder to 2 or 1 10c 11c ----- gauss-hermite quadrature using minimum point formula ----- 12c 13#include "hnd_whermt.fh" 14c 15 double precision xint, yint, zint, t, x0, y0, z0, xi, yi, zi 16 double precision xj, yj, zj, cx, cy, cz, zero, ptx, pty, ptz 17 double precision dum, px, py, pz, ax, ay, az, bx, by, bz 18 integer ni, nj, npts, imin, imax, i, mder, ii, jj 19 common/hnd_xyzder/xint,yint,zint,t,x0,y0,z0,xi,yi,zi,xj,yj,zj 20 1 ,ni,nj,cx,cy,cz 21 data zero /0.0d+00/ 22 xint = zero 23 yint = zero 24 zint = zero 25 npts = (ni+nj+mder-2)/2+1 26 imin = hermin(npts) 27 imax = hermax(npts) 28 do 11 i = imin,imax 29 dum = h(i)*t 30 ptx = dum+x0 31 pty = dum+y0 32 ptz = dum+z0 33 px = ptx-cx 34 py = pty-cy 35 pz = ptz-cz 36 px = px**mder 37 py = py**mder 38 pz = pz**mder 39 ax = ptx-xi 40 ay = pty-yi 41 az = ptz-zi 42 bx = ptx-xj 43 by = pty-yj 44 bz = ptz-zj 45 do ii = 1, ni-1 46 px=px*ax 47 py=py*ay 48 pz=pz*az 49 enddo 50 do jj = 1, nj-1 51 px=px*bx 52 py=py*by 53 pz=pz*bz 54 enddo 55 dum = w(i) 56 xint = xint+dum*px 57 yint = yint+dum*py 58 zint = zint+dum*pz 59 11 continue 60 return 61 end 62