Esope source files

   1 : C ITGA21    SOURCE    CHAT      05/01/13    00:44:35     5004           
   2 : CCC
   3 : C **********************************************************************
   4 : CCC
   5 :       SUBROUTINE INTEGRA21 (XTRI,X,NDIMX,LAM,DDEFPL,NDIMS,
   6 :      .                    nmodel,tolrel,nitmax,nescri,ues,nnumer,
   7 :      .                    deltax,kerre,iiter)
   8 :       IMPLICIT INTEGER(I-N)
   9 :       
  10 :       integer ndims,nmodel,nitmax,nescri,nnumer,ndimx,ues,kdummy
  11 :       real*8    xtri(ndimx),x(ndimx),ddefpl(ndims),lam,tolrel
  12 :       integer i,j,iiter,kerre,npcon,npcap,kdu,ndimv
  13 :       real*8    dl,lres,siginv(3)
  14 :       real*8    auxr1,vtLUiw,auxmax1,auxmax2,deltax
  15 :       real*8    xres(8),dx(8),Amat(64),Gmat(64),vecn(8),vecm(8)
  16 : 
  17 :       kerre=0
  18 :       ndimv=2
  19 :       npcon=ndims+1
  20 :       npcap=ndims+2
  21 :       kdummy=nmodel
  22 :       if (nmodel.eq.22) kdummy=23
  23 :       do i=1,8
  24 :           xres(i)=0.D0
  25 :           dx(i)=0.D0
  26 :           vecn(i)=0.D0
  27 :           vecm(i)=0.D0
  28 :         enddo
  29 :       do i=1,64
  30 :           Amat(i)=0.D0
  31 :           Gmat(i)=0.D0
  32 :         enddo
  33 :       dl=1.D0
  34 :       iiter=-1
  34 :       call MatGenHookinv(Gmat,ndimx,ndims)
  35 : 10    continue
  36 :       iiter=iiter+1
  37 :       kdu=0
  37 :       call yielddMAC(x,ndims,x(npcon),ndimv,lres,kdummy)
  37 :       call VecFlMAC(x,ndimx,vecm,nmodel)
  38 :       do i=1,ndimx
  39 :          xres(i)=lam*vecm(i)
  40 :          do j=1,ndimx
  41 :             xres(i)=xres(i)+Gmat((i-1)*ndimx+j)*(x(j)-xtri(j))
  42 :          enddo
  43 :       enddo
  44 :       auxmax1=dabs(dx(1))
  45 :       do i=2,ndims+1
  46 :          if (dabs(dx(i)).gt.auxmax1) auxmax1=dabs(dx(i))
  47 :       enddo
  48 :       auxmax2=dabs(x(1))
  49 :       do i=2,ndimx
  50 :          if (dabs(x(i)).gt.auxmax2) auxmax2=dabs(x(i))
  51 :       enddo
  52 :         auxr1=max(auxmax1,dabs(dl))/max(auxmax2,dabs(lam))
  53 :       if (iiter.eq.0) auxr1=1.D0
  54 :       if (nescri.eq.1) then
  54 :          call InvariantesPQT(x,ndims,siginv)
  55 :          write(ues,998)iiter,(siginv(i),i=1,3),x(ndimx-1),
  56 :      .                 x(ndimx),auxr1
  57 : 998      format(I3,'Sig ',3(E10.4,1x),'V ',E10.4,1x,E10.4,' L ',E10.4)
  58 :       endif
  59 :       if ((auxr1.lt.tolrel).or.
  60 :      .   ((iiter.eq.nitmax).and.(auxr1.lt.sqrt(tolrel)))) then
  61 : c CONVERGIDO
  62 :          do i=1,ndims
  63 :             ddefpl(i)=0.D0
  64 :             do j=1,ndims
  65 :                ddefpl(i)=ddefpl(i)+Gmat((i-1)*ndimx+j)*(xtri(j)-x(j))
  66 :             enddo
  67 :          enddo
  68 :          return
  69 :       endif
  70 : c NO CONVERGIDO
  71 :       if (iiter.eq.nitmax) then
  72 :          kerre=1
  73 :          return
  74 :       endif
  74 :       call VecYiMAC(x,ndimx,vecn,nmodel)
  74 :       call HessMAC(x,ndimx,Amat,ndimx,nnumer,deltax,nmodel)
  75 :       do i=1,ndimx*ndimx                                                   
  76 :          Amat(i)=Gmat(i)+lam*Amat(i)
  77 :       enddo
  77 :       call DescLU(Amat,ndimx)
  78 :       dl=(lres-vtLUiw(vecn,Amat,xres,ndimx))
  79 :      .         /vtLUiw(vecn,Amat,vecm,ndimx)
  80 :       do i=1,ndimx
  81 :          xres(i)=-xres(i)-dl*vecm(i)
  82 :       enddo
  82 :       call LUiw(Amat,xres,dx,ndimx)
  83 :       lam=lam+dl
  84 :       do i=1,ndimx
  85 :          x(i)=x(i)+dx(i)
  86 :       enddo                              
  87 :       if (x(npcon).lt.0.D0) x(npcon)=0.D0
  88 :       if (x(npcap).lt.0.D0) x(npcap)=0.D0
  89 :       go to 10
  90 :       end
  91 :  
  92 :  
  93 :