********************************************************************************
*                                                                              *
*  2DHF version 1-2003                                                         *
*  Copyright (C) 1996  Jacek Kobus, Leif Laaksonen, Dage Sundholm              *
*                                                                              *
*  This software may be used and distributed according to the terms            *
*  of the GNU General Public License, see README and COPYING.                  *
*                                                                              *
********************************************************************************
c ### difni ###
c
c    This routine calculates 
c
c     (\frac{\partial^2}{\partial\ni^2} + 
c           d(\ni,\mu) \frac{\partial}{\partial \ni}) f(\ni,\mu)
c
c     The function is imersed in the array f(nni+8,nmut+8)
c     To allow for the differentiation over ni variable in the form
c     of matrix times vector the array f has to be transposed first.
c     Results of differentiation are stored directly column-wise in 
c     another matrix which means doing back tranposition before 
c     returning to a calling routine.
c
c     This routine differentiates a function over ni variable.
c     the function is imersed in the array f(nni+6,nmut+6)
c
      subroutine difni (n,f,fd,wtran1,wtran2)
      implicit integer*4 (i-n)
      implicit real*8 (a-h,o-z)
      include 'commons8.inc'
      dimension f(nni+8,n+8),fd(nni+8,n+8),
     &          wtran1(n+8,nni+8),wtran2(n+8,nni+8)

c     To allow for the differentiation over ni variable in the form
c     of matrix times vector the array f has to be transposed first.
c     Results of differentiation will be stored directly
c     column-wise in another matrix which means doing back tranposition
c     before returning to calling routine.  

      call gmtran (f,wtran1,nni+8,n+8)

      do j=1,nni
         call mxv (wtran1(1,j),n+8,dni(1,j),9,wtran2(1,j+4))
      enddo

      call gmtran(wtran2,fd,n+8,nni+8)

      return
      end