1 : C ITGL2L SOURCE CHAT 05/01/13 00:44:39 5004 2 : CCC 3 : C ********************************************************************** 4 : CCC 5 : SUBROUTINE INTEGRA_LS_2LEVELS(XTRI,XINI,NDIMX,LAM,NMODEL,TOLREL, 6 : . jmax,nescri,ues,nnumer,deltax,kerre,jtot) 7 : 8 : IMPLICIT INTEGER(I-N) 9 : integer nmodel,nitmax,nescri,nnumer,ues,j,kerre, 10 : . ndimv,jmax,i,ii,k,n,ndimx,uesb,jl,jtot 11 : real*8 xtri(ndimx),xini(ndimx),lam,tolrel,deltax, 12 : . x(5),dx(5),r(5),err,normainf,norm_l2_2,xnew(5),rnew(5), 13 : . cotaere,lam0,lam1,lam2,f0,fp0,f1,f2,f3,s1,s2,aa,bb,cc, 14 : . sig(4),vecm(4),err0,dla,dlb,ddll,beta 15 : logical escr,last 16 : integer level,kmax 17 : common /linesearch/ level,kmax 18 : integer augla 19 : real*8 c 20 : common /auglagrang1/ augla 21 : common /auglagrang2/ c 22 : 23 : c muestra tensiones 24 : c if (nescri.eq.1) then 24 : c call der_enerelas_dpral(xini,x,nmodel) 25 : c write(ues,*)' Tensiones iniciales' 26 : c write(ues,*)(x(i),i=1,3) 27 : c if (ndimx.eq.4) write(ues,*)' Variable interna inicial',x(4) 28 : c endif 29 : 30 : c muestra iteraciones de line-search 31 : uesb=299 32 : escr=.false. 33 : if (nescri.eq.1) escr=.true. 34 : 35 : c parametros generales 36 : kerre=0 37 : n=ndimx+1 38 : ndimv=1 38 : call equv(x,xini,ndimx) 39 : x(n)=lam 40 : beta=1.D-4 41 : j=0 42 : 43 : jtot=-1 44 : last=.false. 45 : jl=0 46 : 47 : C CALCULOS PARA LAMBDA; SUPONEMOS SIG_INI=SIG_TRIAL 48 : 2 continue 49 : if ((augla.eq.1).and.(jl.eq.0)) then 50 : jl=1 51 : goto 10 52 : endif 53 : jl=jl+1 54 : C A; B 54 : call actualiza_lambda(x,ndimx,nmodel,lam,dla,dlb) 55 : ddll=-dla/dlb 56 : lam=lam+ddll 57 : if (lam.lt.0) then 58 : write(ues,*)' GAMMA negative !',lam,' => ',0.D0 59 : write(ues,*)' dll:',ddll 60 : lam=0 61 : endif 62 : C CHECK DL => GOTO END 63 : jtot=jtot+1 64 : c if (nescri.eq.1) write(ues,999)jtot,jl,lam,ddll,dla 65 : if (last) goto 100 66 : if (abs(ddll).lt.tolrel) last=.true. 67 : if (jl.gt.50+1) then 68 : write(ues,*)' ERROR: loop gamma no ha convergido en', 50 69 : kerre=1 70 : goto 100 71 : endif 72 : 73 : C COMPUTE SIGMA_J FROM LAMBDA_J WITH LINE-SEARCH 74 : 10 continue 75 : j=j+1 75 : call residuo_2levels(x,r,ndimx,ndimx,nmodel,xtri,sig,vecm,lam) 76 : f0=norm_l2_2(r,ndimx) 77 : fp0=-2.D0*f0 77 : call direction_2levels 78 : . (x,r,dx,ndimx,ndimx,nmodel,sig,vecm,nnumer,deltax,lam) 79 : err0=normainf(dx,ndimx) 80 : lam0=1.D0 80 : call addvec(xnew,x,lam0,dx,ndimx) 80 : call residuo_2levels 81 : . (xnew,rnew,ndimx,ndimx,nmodel,xtri,sig,vecm,lam) 82 : f1=norm_l2_2(rnew,ndimx) 83 : cotaere=f0+beta*lam0*fp0 84 : lam1=lam0 85 : if (escr) write(uesb,*)' more ',0,lam0 86 : if (f1.lt.cotaere) goto 80 87 : if ((f0.lt.tolrel).and.(f1.lt.tolrel)) goto 80 88 : if (kmax.le.0) goto 70 89 : 90 : k=1 91 : lam1=(-fp0)/2.D0/(f1-f0-fp0) 92 : if (lam1.lt.0.1D0) lam1=0.1D0 92 : call addvec(xnew,x,lam1,dx,ndimx) 92 : call residuo_2levels(xnew,rnew,n,ndimx,nmodel,xtri,sig,vecm,lam) 93 : f2=norm_l2_2(rnew,ndimx) 94 : cotaere=f0+beta*lam1*fp0 95 : if (escr) write(uesb,*)' more ',k,lam1 96 : if (f2.lt.cotaere) goto 80 97 : if (kmax.le.1) goto 70 98 : 99 : 20 continue 100 : k=k+1 101 : lam2=(-fp0*lam1**2)/2.D0/(f2-f0-fp0*lam1) 102 : if (lam2.lt.(0.1D0*lam1)) lam2=0.1D0*lam1 103 : if (lam2.gt.(0.5D0*lam1)) lam2=0.5D0*lam1 103 : call addvec(xnew,x,lam2,dx,ndimx) 103 : call residuo_2levels(xnew,rnew,n,ndimx,nmodel,xtri,sig,vecm,lam) 104 : f3=norm_l2_2(rnew,ndimx) 105 : cotaere=f0+beta*lam2*fp0 106 : f1=f2 107 : f2=f3 108 : lam0=lam1 109 : lam1=lam2 110 : if (escr) write(uesb,*)' more ',k,lam2 111 : if (f3.lt.cotaere) goto 80 112 : if (k.lt.kmax) goto 20 113 : 70 continue 114 : if (escr) write(uesb,*)' more line search',kmax 115 : 80 continue 116 : err=err0*lam1/normainf(xnew,ndimx) 116 : call equv(x,xnew,ndimx) 117 : c MUESTRA ITERACIONES EN SIGMA_J 118 : c if (nescri.eq.1) then 118 : c call Invari_t(x,3,aa) 118 : c call Invari_q(x,3,bb) 119 : c write(ues,999)jtot,j,(x(i),i=1,3),lam,err 120 : c write(ues,999)j,aa,bb,(x(i),i=4,n),err 121 : c 999 format(2I4,' ',6(E10.4,1x)) 121 : c call der_enerelas_dpral(x,xnew,nmodel) 121 : c call Invari_t(xnew,3,cc) 121 : c call Invari_q(xnew,3,bb) 121 : c call Invari_p(xnew,3,aa) 122 : c write(ues,1999)j,aa,bb,cc*45.D0/atan(1.D0),err,rnew(n) 123 : c1999 format(I3,' ',5(E14.8,1x)) 124 : c endif 125 : if (err.lt.tolrel) goto 99 126 : if (j.lt.jmax) goto 10 127 : kerre=1 128 : write(ues,*)' ERROR: Loop interno sigmas no ha convergido' 129 : 99 continue 130 : c FIN ITERACIONES EN SIGMA_J 131 : jtot=jtot+j 132 : j=0 133 : goto 2 134 : 135 : c RESULTADOS FINALES 136 : 100 continue 136 : call equv(xini,x,ndimx) 137 : 138 : c MUESTRA RESULTADOS 139 : c if (nescri.eq.1) then 140 : c write(ues,*)' Deform elas finales' 141 : c write(ues,*)(x(i),i=1,3) 142 : c write(ues,*)' Deform plas finales' 143 : c write(ues,*)(xtri(i)-x(i),i=1,3) 144 : c if (ndimx.eq.4) write(ues,*)' Variable interna final',x(4) 144 : c call der_enerelas_dpral(x,xini,nmodel) 145 : c write(ues,*)' Tensiones finales' 146 : c write(ues,*)(xini(i),i=1,3) 147 : c if (ndimx.eq.4) write(ues,*)' Variable interna inicial',x(4) 147 : c call equv(xini,x,ndimx) 148 : c endif 149 : end 150 : 151 : 152 : 153 :
© Cast3M 2003 - All rights reserved.
Disclaimer