Télécharger bantou.eso

Retour à la liste

Numérotation des lignes :

  1. C BANTOU SOURCE PV 09/03/12 21:15:33 6325
  2. SUBROUTINE BANTOU
  3. # (TABCOR,TABVAL,LADIM,ZSIG,VALMOY,DELZET,OMMAX)
  4. C
  5. C-----------------------------------------------------------------------
  6. C Cet opérateur génère un champ aléatoire par la méthode des
  7. C bandes tournantes.
  8. C Renvoi un pointeur actif TABVAL sur la liste des valeurs
  9. C Les coordonnées en entrée ont été adimunsiofn©es par Lambda,
  10. C c-à-d qu'on s'est ramené au cas isotrope avec lambda = 1
  11. C-----------------------------------------------------------------------
  12. C
  13. C---------------------------
  14. C Paramètres Entrée/Sortie :
  15. C---------------------------
  16. C
  17. C E/ TABCOR : Pointeur sur la TABLE contenant les coordonnées des
  18. C points où l'on va générer les valeurs (il peut y avoir
  19. C répétition)
  20. C /S TABVAL : Pointeur actif sur la TABLE contenant les valeurs en ces
  21. C points
  22. C E/ LADIM : Indicateur de la dimension relative aux propriétés
  23. C aléatoires du champ à générer
  24. C E/ ZSIG : Ecart-type du champ aléatoire à générer
  25. C E/ DELZET : Incrément d'espace sur la bande (Pas de discrétisation)
  26. C E/ OMMAX : Vecteur d'onde de coupure
  27. C
  28. C--------------------------
  29. C Variables internes :
  30. C--------------------------
  31. C
  32. C NLINE : Nombre de bandes
  33. C XLDOM : Dimension maximale du domaine
  34. C DOMEGA : Incrément de vecteur d'onde
  35. C
  36. C---------------------------
  37. C Ajouts intéressants
  38. C---------------------------
  39. C
  40. C=======================================================================
  41. C
  42. IMPLICIT INTEGER(I-N)
  43. IMPLICIT REAL*8(A-H,O-Z)
  44. C
  45. -INC CCOPTIO
  46. -INC SMCOORD
  47. -INC SMCHAML
  48. -INC SMELEME
  49. -INC CCREEL
  50. C
  51. C TABLE DES COORDONNEES DES POINTS SUPPORTS
  52. SEGMENT TABCOR
  53. REAL*8 COR(NBL,3)
  54. ENDSEGMENT
  55. C TABLE DES VALEURS EN CES POINTS
  56. SEGMENT TABVAL
  57. REAL*8 VAL(NBL)
  58. ENDSEGMENT
  59.  
  60. SEGMENT TH(NLINE)
  61. SEGMENT PHI(NLINE)
  62. SEGMENT ZETMIN(NLINE)
  63. SEGMENT ZETMAX(NLINE)
  64. SEGMENT XLINE(NZET)
  65. C
  66. C EPSILON SERVANT A DIFFERENTS TESTS
  67. C
  68. EPS=1.D-12
  69. C
  70. C INITIALISATION DU CHAMP ALEATOIRE A ZERO
  71. C
  72. SEGACT TABCOR
  73. NBL = COR(/1)
  74. SEGINI TABVAL
  75. DO 20 I=1,NBL
  76. VAL(I) = 0.D0
  77. 20 CONTINUE
  78. C
  79. C ===========================================================
  80. C = Génération des angles décrivant chaque point des lignes =
  81. C ===========================================================
  82. C
  83. C 1D
  84. C ----
  85. C
  86. IF (LADIM.EQ.1) THEN
  87. C Le champ aléatoire a des propriétés unidirectionnelles
  88. C On génère une seule bande (NLINE = 1) de direction (1 0 0)
  89. C On n'a pas besoin de TH ni de PHI
  90. NLINE=1
  91. C
  92. C 2D
  93. C ----
  94. C
  95. ELSEIF (LADIM.EQ.2) THEN
  96. C Le champ aléatoire a des propriétés bidimensionnelles
  97. C Générer des lignes de direction répartie régulièrement entre 0
  98. C et 2*PI
  99. C NLINE est figé à une valeur fixe : 100
  100. C
  101. NLINE=100
  102. SEGINI TH
  103. C
  104. C THETA VARIE ENTRE 0 ET 2*PI
  105. C
  106. DO 1 I=1,NLINE
  107. TH(I) = 2.D0*XPI * FLOAT(I-1) / FLOAT(NLINE)
  108. 1 CONTINUE
  109. C
  110. C 3D
  111. C ----
  112. C
  113. ELSEIF (LADIM.EQ.3) THEN
  114. C Le champ aléatoire a des propriétés tridimensionnelles
  115. C Générer des lignes de direction aléatoire répartie entre 0 et
  116. C 2*PI pour le site, et 0 et PI pour l'azimuth
  117. C NLINE est figé à une valeur fixe : 1000
  118. C
  119. NLINE=1000
  120. SEGINI TH
  121. SEGINI PHI
  122. C
  123. C THETA varie aléatoirement entre 0 et 2*PI, PHI entre 0 et PI
  124. C THETA et COS(PHI) suivent des lois uniformes
  125. C
  126. DO 2 I=1,NLINE
  127. CALL TDRAND(XRAN)
  128. TH(I)=XRAN*XPI*2.D0
  129. CALL TDRAND(XRAN)
  130. PHI(I)=ACOS(1.D0-2.D0*XRAN)
  131. 2 CONTINUE
  132. ENDIF
  133. C
  134. C ===========================================================
  135. C = Détermination de la dimension maximale du domaine XLDOM =
  136. C = et des pas de discrétisation de l'espace dual DOMEGA =
  137. C ===========================================================
  138. C utile uniquement en dimension > 1
  139. C XLDOM est le maximum sur toutes les bandes de l'intervalle des
  140. C valeurs algébriques des positions de la projection de chaque point
  141. C du maillage sur ces droites.
  142.  
  143. SEGINI ZETMIN
  144. SEGINI ZETMAX
  145. XLDOM=-1.D0
  146. C
  147. IF (LADIM.EQ.1) THEN
  148. L = 1
  149. ZETMAX(L) = -1 * XGRAND
  150. ZETMIN(L) = XGRAND
  151. DO 9 I=1,NBL
  152. ZET = COR(I,1)
  153. ZETMIN(L) = MIN (ZETMIN(L), ZET)
  154. ZETMAX(L) = MAX (ZETMAX(L), ZET)
  155. 9 CONTINUE
  156. XLDOM = MAX (XLDOM, ZETMAX(L)-ZETMIN(L))
  157. C
  158. ELSEIF (LADIM.EQ.2) THEN
  159. C
  160. DO 4 L=1,NLINE
  161. ZETMAX(L) = -1 * XGRAND
  162. ZETMIN(L) = XGRAND
  163. THETA = TH(L)
  164. XLI = COS(THETA)
  165. YLI = SIN(THETA)
  166. DO 5 I=1,NBL
  167. X = COR(I,1)
  168. Y = COR(I,2)
  169. ZET = X*XLI + Y*YLI
  170. ZETMIN(L) = MIN (ZETMIN(L), ZET)
  171. ZETMAX(L) = MAX (ZETMAX(L), ZET)
  172. 5 CONTINUE
  173. XLDOM = MAX (XLDOM, ZETMAX(L)-ZETMIN(L))
  174. 4 CONTINUE
  175. C
  176. ELSE
  177. C
  178. DO 6 L=1,NLINE
  179. ZETMAX(L) = -1 * XGRAND
  180. ZETMIN(L) = XGRAND
  181. THETA = TH(L)
  182. FI = PHI(L)
  183. XLI = COS(THETA) * SIN(FI)
  184. YLI = SIN(THETA) * SIN(FI)
  185. ZLI = COS(FI)
  186. DO 7 I=1,NBL
  187. X = COR(I,1)
  188. Y = COR(I,2)
  189. Z = COR(I,3)
  190. ZET = X*XLI + Y*YLI + Z*ZLI
  191. ZETMIN(L) = MIN (ZETMIN(L), ZET)
  192. ZETMAX(L) = MAX (ZETMAX(L), ZET)
  193. 7 CONTINUE
  194. XLDOM = MAX (XLDOM, ZETMAX(L)-ZETMIN(L))
  195. 6 CONTINUE
  196. C
  197. ENDIF
  198. C
  199. C DOMEGA doit être inférieur à pi / L
  200. DOMEGA = .9 * (XPI / XLDOM)
  201. C
  202. C ===========================================================
  203. C = Génération du champ aléatoire =
  204. C ===========================================================
  205. C Sommation sur toutes les bandes de la valeur trouvée sur la bande
  206. C en cours au point projeté du point considéré sur cette bande,
  207. C et ce, pour chaque point champ à créer.
  208. C Les valeurs en deux points confondus sont égales.
  209. C
  210. C 1D
  211. C ----
  212. C
  213. IF (LADIM.EQ.1) THEN
  214.  
  215. C Valeurs aléatoires sur la bande
  216. NZET = INT ((ZETMAX(1) - ZETMIN(1)) / DELZET) + 2
  217. SEGINI XLINE
  218. ZETMIL = ZETMIN(1)
  219. CALL BANT1D(DOMEGA,OMMAX,DELZET,ZETMIL,XLINE)
  220. C
  221. C PROJECT LINE PROCESS TO GENERATE 2-D ARRAY
  222. C
  223. DO 8 I=1,NBL
  224. ZET = COR(I,1)
  225. INDEX = INT ((ZET - ZETMIN(1)) / DELZET + 0.5D0) + 1
  226. VAL(I) = VAL(I) + XLINE(INDEX)
  227. 8 CONTINUE
  228.  
  229. SEGSUP XLINE
  230. C
  231. C 2D
  232. C ----
  233. C
  234. ELSEIF (LADIM.EQ.2) THEN
  235. C
  236. DO 12 L=1,NLINE
  237.  
  238. C Valeurs aléatoires sur la bande
  239. NZET = INT ((ZETMAX(L) - ZETMIN(L)) / DELZET) + 2
  240. SEGINI XLINE
  241. C
  242. ZETMIL = ZETMIN(L)
  243. CALL BANT2D(DOMEGA,OMMAX,DELZET,ZETMIL,XLINE)
  244. C
  245. C PROJECT LINE PROCESS TO GENERATE 2-D ARRAY
  246. C
  247. THETA = TH(L)
  248. XLI = COS(THETA)
  249. YLI = SIN(THETA)
  250. DO 13 I=1,NBL
  251. X = COR(I,1)
  252. Y = COR(I,2)
  253. ZET = X*XLI + Y*YLI
  254. INDEX = INT ((ZET - ZETMIN(L)) / DELZET + 0.5D0) + 1
  255. VAL(I) = VAL(I) + XLINE(INDEX)
  256. 13 CONTINUE
  257. C
  258. SEGSUP XLINE
  259. C
  260. 12 CONTINUE
  261.  
  262. SEGSUP TH
  263. C
  264. C 3D
  265. C ----
  266. C
  267. ELSEIF (LADIM.EQ.3) THEN
  268. C
  269. DO 15 L=1,NLINE
  270.  
  271. C Valeurs aléatoires sur la bande
  272. NZET = INT ((ZETMAX(L) - ZETMIN(L)) / DELZET) + 2
  273. SEGINI XLINE
  274. ZETMIL = ZETMIN(L)
  275. CALL BANT3D(DOMEGA,OMMAX,DELZET,ZETMIL,XLINE)
  276. C
  277. C PROJECT LINE PROCESS TO GENERATE 3-D ARRAY
  278. C
  279. THETA = TH(L)
  280. FI = PHI(L)
  281. XLI = COS(THETA) * SIN(FI)
  282. YLI = SIN(THETA) * SIN(FI)
  283. ZLI = COS(FI)
  284. DO 16 I=1,NBL
  285. X = COR(I,1)
  286. Y = COR(I,2)
  287. Z = COR(I,3)
  288. ZET = X*XLI + Y*YLI + Z*ZLI
  289. INDEX = INT ((ZET - ZETMIN(L)) / DELZET + 0.5D0) + 1
  290. VAL(I) = VAL(I) + XLINE(INDEX)
  291. 16 CONTINUE
  292. C
  293. SEGSUP XLINE
  294. C
  295. 15 CONTINUE
  296.  
  297. SEGSUP TH
  298. SEGSUP PHI
  299. ENDIF
  300. C
  301. C.....DIVIDE BY SQRT(FLOAT(NLINE))
  302. C et mise à la moyenne
  303. C
  304. XX1 = ZSIG / SQRT(FLOAT(NLINE))
  305.  
  306. DO 18 I=1,NBL
  307. VAL(I) = (VAL(I)*XX1) + VALMOY
  308. 18 CONTINUE
  309. C
  310. SEGSUP ZETMIN
  311. SEGSUP ZETMAX
  312. C
  313. RETURN
  314. END
  315.  
  316.  
  317.  
  318.  
  319.  
  320.  
  321.  

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