Télécharger fpma3d.eso

Retour à la liste

Numérotation des lignes :

fpma3d
  1. C FPMA3D SOURCE OF166741 25/02/06 21:15:05 12146
  2.  
  3. C____________________________________________________________________
  4. C CALCULE LES FORCES DE PRESSIONS SUR LES FACES D ELEMENTS
  5. C MASSIFS TRIDIMENSIONNELS
  6. C
  7. C ENTREES :
  8. C ---------
  9. C
  10. C IPTVPR POINTEUR SUR UN MELVAL CONTENANT LES PRESSIONS APPLIQUEES
  11. C 0 SI ON A DONNE UNE PRESSION CONSTANTE
  12. C IPMAIL POINTEUR SUR UN OBJET GEOMETRIQUE
  13. C IPTINT POINTEUR SUR UN MINTE CONTENANT LES POINTS D INTEGRATION
  14. C ACTIF EN ENTREE ET EN SORTIE SANS MODIFICATION
  15. C IVAFOR POINTEUR SUR UN MPTVAL ET LES MELVAL CONTENANT LES FORCES
  16. C NODALES RESUL
  17. C
  18. C JACQUELINE BROCHARD AVRIL 85
  19. C
  20. C PASSAGE AUX NOUVEAUX CHAMELEM PAR JM CAMPENON LE 17 09 90
  21. C
  22. C______________________________________________________________________
  23.  
  24. SUBROUTINE FPMA3D(IPTVPR,IPMAIL,ipmaim,IPTINT,IVAFOR,XP
  25. & ,netn1,ietn1)
  26.  
  27. IMPLICIT INTEGER(I-N)
  28. IMPLICIT REAL*8 (A-H,O-Z)
  29.  
  30. -INC PPARAM
  31. -INC CCOPTIO
  32.  
  33. -INC SMCHAML
  34. -INC SMELEME
  35. -INC SMINTE
  36. -INC SMCOORD
  37.  
  38. segment netn(notn)
  39. segment ietn(letn)
  40.  
  41. SEGMENT WORK
  42. REAL*8 XE(3,NBNN)
  43. ENDSEGMENT
  44.  
  45. SEGMENT MPTVAL
  46. INTEGER IPOS(NS),NSOF(NS)
  47. INTEGER IVAL(NCOSOU)
  48. CHARACTER*16 TYVAL(NCOSOU)
  49. ENDSEGMENT
  50.  
  51. idimp1 = IDIM+1
  52. * prob optimiseur il faut initialiser melva1
  53. MELVA1=IVAFOR
  54. IF (IPTVPR.NE.0) THEN
  55. MELVA1=IPTVPR
  56. c* SEGACT,MELVA1 <- ACTIF EN E/S
  57. IG11 = MELVA1.VELCHE(/1)
  58. IB12 = MELVA1.VELCHE(/2)
  59. ENDIF
  60.  
  61. MINTE=IPTINT
  62. C* SEGACT,MINTE <- ACTIF EN E/S
  63. NBPGAU=POIGAU(/1)
  64.  
  65. MELEME = IPMAIL
  66. c* SEGACT,MELEME <- ACTIF EN E/S
  67. NBNN = meleme.NUM(/1)
  68. NBELEM = meleme.NUM(/2)
  69.  
  70. SEGINI,WORK
  71.  
  72. netn = netn1
  73. ietn = ietn1
  74. ipt1 = ipmaim
  75. IF (IPT1.GT.0) THEN
  76. if (netn.eq.0 .or. ietn.eq.0) then
  77. write(ioimp,*) 'FPMA3D : incompatibilite netn, ietn & IPMAIM'
  78. endif
  79. c* SEGACT,IPT1 <- ACTIF en E/S
  80. nbnn1 = ipt1.num(/1)
  81. nbel1 = ipt1.num(/2)
  82. ELSE
  83. if (netn.gt.0 .or. ietn.gt.0) then
  84. write(ioimp,*) 'FPMA2D : incompatibilite netn, ietn & IPMAIM'
  85. endif
  86. ENDIF
  87. C
  88. C BOUCLE SUR LES ELEMENTS
  89. C
  90. DO IB = 1, NBELEM
  91.  
  92. CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE)
  93.  
  94. XFLOT = +1.D0
  95. IF (netn.GT.0) THEN
  96. DO inf = 1, nbnn
  97. ip = meleme.num(inf,ib)
  98. ideb = netn(ip)+1
  99. ifin = netn(ip+1)
  100. do itn = ideb, ifin
  101. IEM = ietn(itn)
  102. jne = 0
  103. do i = 1, nbnn
  104. ino = num(i,ib)
  105. do i1 = 1, nbnn1
  106. if (ino.eq.ipt1.num(i1,IEM)) jne=jne+1
  107. enddo
  108. enddo
  109. if (jne.eq.nbnn) goto 170
  110. enddo
  111. ENDDO
  112. CALL ERREUR(26)
  113. GOTO 9900
  114. 170 continue
  115. XG = 0.D0
  116. YG = 0.D0
  117. ZG = 0.D0
  118. DO I = 1, NBNN1
  119. ino = (IPT1.NUM(I,IEM)-1)*idimp1
  120. XG=XG+XCOOR(ino+1)
  121. YG=YG+XCOOR(ino+2)
  122. ZG=ZG+XCOOR(ino+3)
  123. ENDDO
  124. XG=XG / NBNN1
  125. YG=YG / NBNN1
  126. ZG=ZG / NBNN1
  127.  
  128. XK=0.D0
  129. YK=0.D0
  130. ZK=0.D0
  131. DO i = 1,NBNN
  132. XK=XK+XE(1,I)
  133. YK=YK+XE(2,I)
  134. ZK=ZK+XE(3,I)
  135. ENDDO
  136. XK=XK/NBNN
  137. YK=YK/NBNN
  138. ZK=ZK/NBNN
  139.  
  140. V_1 = XG - XK
  141. V_2 = YG - YK
  142. V_3 = ZG - ZK
  143. r_z = 1.D0 / SQRT(V_1*V_1+V_2*V_2+V_3*V_3)
  144. V_1 = V_1 * r_z
  145. V_2 = V_2 * r_z
  146. V_3 = V_3 * r_z
  147. ENDIF
  148. C
  149. C BOUCLE SUR LES POINTS DE GAUSS
  150. C
  151. DO IGAU = 1, NBPGAU
  152.  
  153. VNQSI1 = 0.D0
  154. VNQSI2 = 0.D0
  155. VNQSI3 = 0.D0
  156. VNETA1 = 0.D0
  157. VNETA2 = 0.D0
  158. VNETA3 = 0.D0
  159. DO I = 1, NBNN
  160. VNQSI1 = VNQSI1+SHPTOT(2,I,IGAU)*XE(1,I)
  161. VNQSI2 = VNQSI2+SHPTOT(2,I,IGAU)*XE(2,I)
  162. VNQSI3 = VNQSI3+SHPTOT(2,I,IGAU)*XE(3,I)
  163. VNETA1 = VNETA1+SHPTOT(3,I,IGAU)*XE(1,I)
  164. VNETA2 = VNETA2+SHPTOT(3,I,IGAU)*XE(2,I)
  165. VNETA3 = VNETA3+SHPTOT(3,I,IGAU)*XE(3,I)
  166. ENDDO
  167.  
  168. VNN1 = VNQSI2*VNETA3-VNQSI3*VNETA2
  169. VNN2 = VNQSI3*VNETA1-VNQSI1*VNETA3
  170. VNN3 = VNQSI1*VNETA2-VNQSI2*VNETA1
  171.  
  172. if (igau.eq.1.and.netn.gt.0) then
  173. vnnn = 1.D0 / SQRT(vnn1*vnn1+vnn2*vnn2+vnn3*vnn3)
  174. test = v_1*(vnn1*vnnn) + v_2*(vnn2*vnnn) + v_3*(vnn3*vnnn)
  175. if (test.lt.0d0) xflot = -1.d0
  176. endif
  177.  
  178. r_z = POIGAU(IGAU) * XFLOT
  179. IF (IPTVPR.NE.0) THEN
  180. IGMN=MIN(IGAU,IG11)
  181. IBMN=MIN(IB ,IB12)
  182. r_z = r_z * MELVA1.VELCHE(IGMN,IBMN)
  183. ELSE
  184. r_z = r_z * XP
  185. ENDIF
  186.  
  187. T1 = r_z * VNN1
  188. T2 = r_z * VNN2
  189. T3 = r_z * VNN3
  190.  
  191. MPTVAL=IVAFOR
  192. MELVAL=IVAL(1)
  193. DO i=1,NBNN
  194. VELCHE(i,IB) = VELCHE(i,IB) + SHPTOT(1,i,IGAU)*T1
  195. ENDDO
  196. MELVAL=IVAL(2)
  197. DO i=1,NBNN
  198. VELCHE(i,IB) = VELCHE(i,IB) + SHPTOT(1,i,IGAU)*T2
  199. ENDDO
  200. MELVAL=IVAL(3)
  201. DO i=1,NBNN
  202. VELCHE(i,IB) = VELCHE(i,IB) + SHPTOT(1,i,IGAU)*T3
  203. ENDDO
  204.  
  205. ENDDO
  206.  
  207. ENDDO
  208.  
  209. 9900 CONTINUE
  210. SEGSUP,WORK
  211.  
  212. c* SEGDES,MINTE <- ACTIF en E/S
  213. c* SEGDES,MELEME <- ACTIF en E/S
  214. c* IF (IPTVPR.NE.0) SEGDES,MELVA1 <- ACTIF en E/S
  215.  
  216. RETURN
  217. END
  218.  
  219.  
  220.  

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