Télécharger chole1.eso

Retour à la liste

Numérotation des lignes :

chole1
  1. C CHOLE1 SOURCE PV090527 24/01/12 21:15:02 11821
  2. FUNCTION CHOLE1(ILIGF,iprell,imasql,VALF,DAAG,IPKNO,IPPVF,KHG,
  3. >IVPOF,KIDEP,KI1,KQ,idep,prec,nbo)
  4. IMPLICIT INTEGER(I-N)
  5. IMPLICIT REAL*8 (A-H,O-Z)
  6. -INC SMMATRI
  7. -INC SMRIGID
  8. -INC CCHOLE
  9. -INC CCREEL
  10. * DAAG contient l'inverse de la diagonale de facon a faire des multiplications et non des divisions
  11. DIMENSION ILIGF(*),VALF(*),DAAG(*),IPKNO(*),IPPVF(*),IVPOF(*)
  12. dimension imasql(*)
  13. * constante pour la somme compensee utilisee pour l'accumulation. CF ddot2
  14. parameter (c=2.D0**26+1.D0)
  15. xmatri=matric
  16. nbnnma=nbnnmc
  17. IPPKHG=IPPVF(KHG)
  18. KBAS=IPKNO(KIDEP)
  19. KHAU=IPKNO(KI1)
  20. KDIAG=KI1+1
  21. DNORM=ABS(VALF(KDIAG))*PREC
  22. KPREM=IVPOF(KHG)-IPPKHG
  23. ILIG=IPRELL+KHG-NBNNMA-1
  24. IECAR=KQ-IPRELL+1
  25. DO 30 NNJ=MAX(1,KIDEP+IECAR),KI1+IECAR
  26. KK=NNJ-IECAR
  27. ICOL=KQ+KK-NBNNMA
  28. NNJJ=IPPVF(NNJ+1)
  29. NJ=NNJJ-IPPVF(NNJ)
  30. LLOL=MIN(NJ,KK)-1
  31. LLON=MIN(LLOL-KK+KPREM+1,LLOL-NNJJ+IVPOF(NNJ)+1)
  32. LLON=MIN(LLON,NBNNMA-KQ-KK+LLOL+1)
  33. IF (LLON.GT.0.and.kk.ge.1) THEN
  34. IEC1=KK-LLOL-1
  35. IEC2=NNJJ-IPPKHG-KK
  36. ideq=1+iec1+idep-1
  37. if (llon.gt.masdim+1) then
  38. p=ddotpw(llon,VALF(1+iec1),VALF(1+iec1+iec2),
  39. > imasql(1),1+idep-1,nbo)
  40. else
  41. if (imasql(ideq/masdim+1).gt.0.or.
  42. > imasql((ideq+llon-1)/masdim+1).gt.0) then
  43. p=ddotpv(llon,VALF(1+iec1),VALF(1+iec1+iec2))
  44. nbo=nbo+llon
  45. else
  46. p=0.d0
  47. endif
  48. endif
  49. VALF(KK)=VALF(KK)-P
  50. IF (ABS(VALF(KK)).GT.DNORM) then
  51. KPREM=KK
  52. imasql((kk+idep-1)/masdim+1) =1
  53. * si on remonte, on tombe au terme diagonal ou apres, mais ce n'est qu'un seul terme
  54. imasql(kk/masdim+1) =1
  55. ELSE
  56. * annuler le terme car on l'ignorera par la suite
  57. valf(kk)=0.d0
  58. ENDIF
  59. ENDIF
  60. if (ilig.ge.1.and.icol.ge.1) then
  61. RE(ILIG,ICOL,1)=VALF(KK)
  62. RE(ICOL,ILIG,1)=VALF(KK)
  63. endif
  64. 30 CONTINUE
  65. 50 CONTINUE
  66. ** AUX1=0.D0
  67. s1=0.d0
  68. c1=0.d0
  69. kdeb=1
  70. 43 continue
  71. kdebi=kdeb
  72. 44 continue
  73. do 100 im=kdeb/masdim+1,kprem/masdim+1
  74. jm=imasql(im)
  75. if (jm.gt.0) goto 105
  76. if (jm.eq.0) goto 100
  77. jinio=-jm/masdim+1
  78. if (jinio.gt.im+jacc) then
  79. * write (6,*) 'saut kdeb jm ',kdeb,jm
  80. kdeb=-jm
  81. goto 44
  82. endif
  83. 100 continue
  84. 105 continue
  85. ideb=max((im-1)*masdim,kdebi)
  86. kdeb=ideb
  87. 111 continue
  88. do 110 im=kdeb/masdim+1,kprem/masdim+1
  89. jm=imasql(im)
  90. if (jm.le.0) goto 115
  91. if (jm.eq.1) goto 110
  92. jfineo=jm/masdim+1
  93. if(jfineo.gt.im+jacc) then
  94. kdeb=jm
  95. goto 111
  96. endif
  97. 110 continue
  98. 115 continue
  99. im=im-1
  100. ifin=min(im*masdim-1,kprem)
  101. ** write (6,*) ' chole1 kdeb kprem ideb ifin ',kdeb,kprem,ideb,ifin
  102. DO 9 K=ideb,min(ifin,nbnnma-kq)
  103. AUX=VALF(K)
  104. if (aux.eq.0.d0) goto 9
  105. nbo=nbo+1
  106. VALFK=AUX*DAAG(K)
  107. VALF(K)=VALFK
  108. ** AUX1=AUX1+AUX*VALF(K)
  109. ** twoprod
  110. *pv aux2=valf(k)
  111. *pv z=c*aux2
  112. *pv xx=z-(z-aux2)
  113. *pv xy=aux2-xx
  114. *pv z=c*aux
  115. *pv yx=z-(z-aux)
  116. *pv yy=aux-yx
  117. *pv s2=aux2*aux
  118. *pv c2=xy*yy-(((s2-xx*yx)-xx*yy)-xy*yx)
  119. ** twosum
  120. s2=aux*valfk
  121. xx=s1+s2
  122. z=xx-s1
  123. xy=(s1-(xx-z))+(s2-z)
  124. s1=xx
  125. *pv c1=c1+(xy+c2)
  126. c1=c1+xy
  127. 9 CONTINUE
  128. if (ifin.lt.kprem) then
  129. kdeb=ifin+1
  130. goto 43
  131. endif
  132. aux1=s1+c1
  133. ivpof(khg)=kprem+ippkhg
  134. CHOLE1=-AUX1
  135. RETURN
  136. END
  137.  
  138.  
  139.  
  140.  
  141.  
  142.  
  143.  
  144.  
  145.  
  146.  
  147.  
  148.  
  149.  
  150.  
  151.  
  152.  
  153.  
  154.  
  155.  
  156.  
  157.  
  158.  
  159.  
  160.  
  161.  
  162.  
  163.  
  164.  
  165.  
  166.  
  167.  

© Cast3M 2003 - Tous droits réservés.
Mentions légales