1C$Id$ 2************************************************************************ 3* * 4 subroutine ecp_cstrans (l,ncart,ntr,l_hi,l_lo,nsph,cart,ldc, 5 & sph,lds,csco,lcsco,lstart,indicator,position,mode) 6* * 7* Transform from cartesian monomials x^i y^j z^k to real solid * 8* spherical tensor functions r^n G_lm(theta,phi). These are * 9* sqrt(2pi(1+delta_m,0)/2l+1) times the real spherical harmonics. * 10* * 11* l (inp) - angular momentum of shell, l = i+j+k * 12* ncart (inp) - number of cartesian components = (l+1)*(l+2)/2 * 13* ntr (inp) - number of quantities to transform * 14* l_hi - maximum angular momentum of spherical quantities * 15* l_lo - minimum angular momentum of spherical quantities * 16* nsph (out) - number of spherical components * 17* cart (i/o) - array of cartesian quantities. * 18* ldc (inp) - leading dimension of cart * 19* sph (i/o) - array of spherical quantities * 20* lds (inp) - leading dimension of sph * 21* csco (inp) - transformation arrays * 22* lcsco (inp) - length of transformation arrays * 23* lstart - pointers to beginning of transformation for given l * 24* indicator (inp) - +1 = cartesian to spherical, * 25* -1 = spherical to cartesian * 26* position (inp) - +1 = post-multiply, array is ntr*(ncart or nsph) * 27* -1 = premultiply, array is (ncart or nsph)*ntr * 28* mode (inp) - +1 = transform with regular matrix, * 29* -1 = transform with inverse matrix * 30* * 31* Written by K. G. Dyall * 32* * 33************************************************************************ 34 implicit none 35#include "errquit.fh" 36 integer l,ncart,ntr,l_hi,l_lo,nsph,ldc,lds,lcsco,lstart(0:1,0:l), 37 & indicator,position,mode 38 integer i,j,offset_start,offset_end 39 double precision cart(ldc,ncart),sph(lds,ncart) 40 double precision csco(lcsco),zero,one 41 parameter (zero = 0.0d00, one = 1.0d00) 42* 43* Test validity of input parameters 44* 45 if (abs(indicator) .ne. 1) call errquit ( 46 & 'Illegal value of indicator in ecp_cstrans',99, BASIS_ERR) 47 if (abs(position) .ne. 1) call errquit ( 48 & 'Illegal value of position in ecp_cstrans',99, BASIS_ERR) 49 if (abs(mode) .ne. 1) call errquit ( 50 & 'Illegal value of mode in ecp_cstrans',99, BASIS_ERR) 51* 52C write (6,*) 'l_hi,l_lo',l_hi,l_lo 53 offset_start = (l_hi+1)*(l_hi+2)/2 54 offset_end = (l_lo-1)*(l_lo)/2 55C write (6,*) 'offset_start,offset_end',offset_start,offset_end 56 nsph = offset_start-offset_end 57* 58 if (position .gt. 0) then 59 if (lds .lt. ntr) call errquit ( 60 & ' lds < ntr in ecp_cstrans',99, BASIS_ERR) 61 if (ldc .lt. ntr) call errquit ( 62 & ' ldc < ntr in ecp_cstrans',99, BASIS_ERR) 63 else 64 if (lds .lt. nsph) call errquit ( 65 & ' lds < nsph in ecp_cstrans',99, BASIS_ERR) 66 if (ldc .lt. ncart) call errquit ( 67 & ' ldc < ncart in ecp_cstrans',99, BASIS_ERR) 68 end if 69* 70 i = lstart(0,l)+ncart*(ncart-offset_start) 71 j = lstart(1,l)+ncart-offset_start 72C write (6,*) lstart(0,l),lstart(1,l) 73C write (6,*) 'l,i,j',l,i,j 74C call ecp_matpr(csco(i),1,ncart,1,ncart,1,ncart,1,nsph, 75C & 'Cart to sph','F',81,5) 76C call ecp_matpr(csco(j),1,ncart,1,ncart,1,nsph,1,ncart, 77C & 'Sph to cart','F',81,5) 78* 79* Cartesian to spherical 80* 81 if (indicator .gt. 0) then 82 if (position .gt. 0) then 83 if (mode .gt. 0) then 84 call dgemm ('N','N',ntr,nsph,ncart,one,cart,ldc, 85 & csco(i),ncart,zero,sph,lds) 86 else 87 call dgemm ('N','T',ntr,nsph,ncart,one,cart,ntr, 88 & csco(j),ncart,zero,sph,lds) 89 end if 90 else 91 if (mode .gt. 0) then 92 call dgemm ('T','N',nsph,ntr,ncart,one,csco(i), 93 & ncart,cart,ldc,zero,sph,lds) 94 else 95 call dgemm ('N','N',nsph,ntr,ncart,one,csco(j), 96 & ncart,cart,ldc,zero,sph,lds) 97 end if 98 end if 99* 100* Spherical to cartesian 101* 102 else 103 if (position .gt. 0) then 104 if (mode .gt. 0) then 105 call dgemm ('N','N',ntr,ncart,nsph,one,sph,lds, 106 & csco(j),ncart,zero,cart,ldc) 107 else 108 call dgemm ('N','T',ntr,ncart,nsph,one,sph,lds, 109 & csco(i),ncart,zero,cart,ldc) 110 end if 111 else 112 if (mode .gt. 0) then 113 call dgemm ('T','N',ncart,ntr,nsph,one,csco(j), 114 & ncart,sph,lds,zero,cart,ldc) 115 else 116 call dgemm ('N','N',ncart,ntr,nsph,one,csco(i), 117 & ncart,sph,lds,zero,cart,ldc) 118 end if 119 end if 120 end if 121* 122 return 123 end 124