Esope source files

   1 : C ITGA1J    SOURCE    CHAT      05/01/13    00:44:30     5004           
   2 : CCC
   3 : C **********************************************************************
   4 : CCC
   5 :       subroutine Integra1J2 (xtri,x,ndimx,ndimv,lam,ddefpl,ndims,
   6 :      .                       tolrel,nitmax,nescri,ues,kerre,iiter)
   7 :       IMPLICIT INTEGER(I-N)
   8 :       IMPLICIT REAL*8 (A-H,O-Z)
   9 :       integer ndims,nitmax,nescri,ndimx,ues
  10 :       real*8  xtri(*),x(*),ddefpl(ndims),lam,tolrel
  11 :       integer i,j,iiter,kerre,ndimv
  12 :       real*8  dl,lres,siginv(3),auxr1,vtLUiw
  13 :       real*8  xres(7),dx(7),Amat(49),Gmat(49),vecn(7),vecm(7)
  14 :       kerre=0
  15 :         do i=1,7
  16 :          xres(i)=0.D0
  17 :          dx(i)=0.D0
  18 :          vecn(i)=0.D0
  19 :          vecm(i)=0.D0
  20 :         enddo
  21 :         do i=1,49
  22 :          amat(i)=0.D0
  23 :          gmat(i)=0.D0
  24 :         enddo
  25 :       iiter=-1
  25 :       call MatGenHookinv(Gmat,ndimx,ndims)
  26 :         dl=1.D0
  27 : 10    continue
  28 :       iiter=iiter+1
  28 :       call yielddJ2(x,ndims,x(ndimx),ndimv,lres)
  28 :       call vflsigJ2(x,ndims,x(ndimx),ndimv,vecm)
  28 :       call vflvarJ2(x,ndims,x(ndimx),ndimv,vecm(ndimx))
  29 :       do i=1,ndimx
  30 :          xres(i)=lam*vecm(i)
  31 :          do j=1,ndimx
  32 :             xres(i)=xres(i)+Gmat((i-1)*ndimx+j)*(x(j)-xtri(j))
  33 :          enddo
  34 :       enddo
  35 :       auxmax1=abs(dx(1))
  36 :       auxmax2=abs(x(1))
  37 :       do i=2,ndimx
  38 :        if (abs(dx(i)).gt.auxmax1) auxmax1=abs(dx(i))
  39 :        if (abs(x(i)).gt.auxmax2) auxmax2=abs(x(i))
  40 :       enddo
  41 :       auxr1=max(auxmax1,abs(dl))/max(auxmax2,1.D0)
  42 :       if (nescri.eq.1) then
  43 :          write(ues,'(I3,X,8(1pE8.2,1x))')
  44 :      .        iiter,(abs(xres(i)),i=1,ndimx),abs(lres),auxr1    
  45 :       endif
  46 :       if (auxr1.lt.tolrel) then
  47 : c CONVERGIDO
  48 :          do i=1,ndims
  49 :             ddefpl(i)=0.D0
  50 :             do j=1,ndims
  51 :                ddefpl(i)=ddefpl(i)+Gmat((i-1)*ndimx+j)*(xtri(j)-x(j))
  52 :             enddo
  53 :          enddo
  54 :          return
  55 :       endif
  56 : c NO CONVERGIDO
  57 :       if (iiter.eq.nitmax) then
  58 :          kerre=1
  59 :          return
  60 :       endif
  60 :       call vyisigJ2(x,ndims,x(ndimx),ndimv,vecn)
  60 :       call vyivarJ2(x,ndims,x(ndimx),ndimv,vecn(ndimx))
  60 :       call HessFlsigJ2(x,ndims,x(ndimx),ndimv,Amat,ndimx)
  61 :       do i=1,ndimx*ndimx                                                   
  62 :          Amat(i)=Gmat(i)+lam*Amat(i)
  63 :       enddo
  63 :       call DescLU(Amat,ndimx)
  64 :       dl=(lres-vtLUiw(vecn,Amat,xres,ndimx))
  65 :      .         /vtLUiw(vecn,Amat,vecm,ndimx)
  66 :       do i=1,ndimx
  67 :          xres(i)=-xres(i)-dl*vecm(i)
  68 :       enddo
  68 :       call LUiw(Amat,xres,dx,ndimx)
  69 :       lam=lam+dl
  70 :       do i=1,ndimx
  71 :          x(i)=x(i)+dx(i)
  72 :       enddo                              
  73 :       if (x(ndimx).lt.0.D0) x(ndimx)=0.D0
  74 :       go to 10
  75 :       end
  76 :  
  77 :  
  78 :