Sources Esope | Cast3M

jac0j2
C JAC0J2    SOURCE    PV        22/04/25    21:15:08     11344          
CCC
C **********************************************************************
CCC
      SUBROUTINE JAC0J2 (KMAT,NDIMK,X,NDIMS,NDIMV,LAM,
     .                   nescri,ues,kerre)
      IMPLICIT INTEGER(I-N)
      integer ndims,nescri,ues,ndimv,i,j,k,kerre,ndimk
      real*8  x(ndims),lam,void(1),Kmat(ndimk,ndimk),g
      real*8  Amat(36),Gmat(36),vecn(6),vecm(6),Avecm(6),vecntA(6),
     .        vecEm(6),Hmat(36)
      kerre=0
      void(1) =0.D0
        do i=1,6
         Avecm(i)=0.D0
         vecntA(i)=0.D0
         vecn(i)=0.D0
         vecm(i)=0.D0
         vecEm(i)=0.D0
        enddo
        do i=1,36
         amat(i)=0.D0
         gmat(i)=0.D0
         hmat(i)=0.D0
        enddo
c m
      call vflsigJ2(x,ndims,void,ndimv,vecm)
c n
      call vyisigJ2(x,ndims,void,ndimv,vecn)
c H
      call HessFlsigJ2(x,ndims,void,ndimv,Hmat,ndims)
c G=E
      call MatGenHook(Gmat,ndims,ndims)
c Em = E*m
      do i=1,ndims
       do j=1,ndims
         vecEm(i)=vecEm(i)+Gmat(ndims*(j-1)+i)*vecm(j)
       enddo
      enddo
c A = I+l*E*H
      do i=1,ndims
       do j=1,ndims
        do k=1,ndims
          Amat(ndims*(j-1)+i)=Amat(ndims*(j-1)+i)+
     .                        Gmat(ndims*(k-1)+i)*Hmat(ndims*(j-1)+k)
        enddo
        Amat(ndims*(j-1)+i)=lam*Amat(ndims*(j-1)+i)
       enddo
       Amat(ndims*(i-1)+i)=1.D0+Amat(ndims*(i-1)+i)
      enddo
c G = A-1
      call DescLU(Amat,ndims)
      call LUinv(Amat,Gmat,ndims)
c g = nt*A-1*Em
      g=0.D0
      do i=1,ndims
       do j=1,ndims
          g=g+vecn(i)*Gmat(i+ndims*(j-1))*vecEm(j)
       enddo
      enddo
c G*Em; nt*G
      do i=1,ndims
       do j=1,ndims
         Avecm(i)=Avecm(i)+Gmat(i+ndims*(j-1))*vecEm(j)
         vecntA(i)=vecntA(i)+vecn(j)*Gmat(j+ndims*(i-1))
       enddo
      enddo
c K = G-(G*m * nt*G)/g
      do i=1,ndims
       do j=1,ndims
          Kmat(i,j)=Gmat((j-1)*ndims+i)-Avecm(i)*vecntA(j)/g
       enddo
      enddo
      return
      end