C ITGA21 SOURCE CB215821 16/04/21 21:17:18 8920 CCC C ********************************************************************** CCC SUBROUTINE INTEGRA21 (XTRI,X,NDIMX,LAM,DDEFPL,NDIMS, . nmodel,tolrel,nitmax,nescri,ues,nnumer, . deltax,kerre,iiter) IMPLICIT INTEGER(I-N) integer ndims,nmodel,nitmax,nescri,nnumer,ndimx,ues,kdummy real*8 xtri(ndimx),x(ndimx),ddefpl(ndims),lam,tolrel integer i,j,iiter,kerre,npcon,npcap,kdu,ndimv real*8 dl,lres,siginv(3) real*8 auxr1,vtLUiw,auxmax1,auxmax2,deltax real*8 xres(8),dx(8),Amat(64),Gmat(64),vecn(8),vecm(8) kerre=0 ndimv=2 npcon=ndims+1 npcap=ndims+2 kdummy=nmodel if (nmodel.eq.22) kdummy=23 do i=1,8 xres(i)=0.D0 dx(i)=0.D0 vecn(i)=0.D0 vecm(i)=0.D0 enddo do i=1,64 Amat(i)=0.D0 Gmat(i)=0.D0 enddo dl=1.D0 iiter=-1 call MatGenHookinv(Gmat,ndimx,ndims) 10 continue iiter=iiter+1 kdu=0 call yielddMAC(x,ndims,x(npcon),ndimv,lres,kdummy) call VecFlMAC(x,ndimx,vecm,nmodel) do i=1,ndimx xres(i)=lam*vecm(i) do j=1,ndimx xres(i)=xres(i)+Gmat((i-1)*ndimx+j)*(x(j)-xtri(j)) enddo enddo auxmax1=ABS(dx(1)) do i=2,ndims+1 if (ABS(dx(i)).gt.auxmax1) auxmax1=ABS(dx(i)) enddo auxmax2=ABS(x(1)) do i=2,ndimx if (ABS(x(i)).gt.auxmax2) auxmax2=ABS(x(i)) enddo auxr1=max(auxmax1,ABS(dl))/max(auxmax2,ABS(lam)) if (iiter.eq.0) auxr1=1.D0 if (nescri.eq.1) then call InvariantesPQT(x,ndims,siginv) write(ues,998)iiter,(siginv(i),i=1,3),x(ndimx-1), . x(ndimx),auxr1 998 format(I3,'Sig ',3(E10.4,1x),'V ',E10.4,1x,E10.4,' L ',E10.4) endif if ((auxr1.lt.tolrel).or. . ((iiter.eq.nitmax).and.(auxr1.lt.sqrt(tolrel)))) then c CONVERGIDO do i=1,ndims ddefpl(i)=0.D0 do j=1,ndims ddefpl(i)=ddefpl(i)+Gmat((i-1)*ndimx+j)*(xtri(j)-x(j)) enddo enddo return endif c NO CONVERGIDO if (iiter.eq.nitmax) then kerre=1 return endif call VecYiMAC(x,ndimx,vecn,nmodel) call HessMAC(x,ndimx,Amat,ndimx,nnumer,deltax,nmodel) do i=1,ndimx*ndimx Amat(i)=Gmat(i)+lam*Amat(i) enddo call DescLU(Amat,ndimx) dl=(lres-vtLUiw(vecn,Amat,xres,ndimx)) . /vtLUiw(vecn,Amat,vecm,ndimx) do i=1,ndimx xres(i)=-xres(i)-dl*vecm(i) enddo call LUiw(Amat,xres,dx,ndimx) lam=lam+dl do i=1,ndimx x(i)=x(i)+dx(i) enddo if (x(npcon).lt.0.D0) x(npcon)=0.D0 if (x(npcap).lt.0.D0) x(npcap)=0.D0 go to 10 end