Esope source files

   1 : C CCGAHU    SOURCE    GOUNAND   07/07/05    21:15:08     5784           
   2 :       SUBROUTINE CCGAHU(LCOF,NOMLOI,
   3 :      $     FC,
   4 :      $     IMPR,IRET)
   5 :       IMPLICIT REAL*8 (A-H,O-Z)
   6 :       IMPLICIT INTEGER (I-N)
   7 : C***********************************************************************
   8 : C NOM         : CCGAHU
   9 : C DESCRIPTION : Calcul de la loi de comportement aux points de Gauss :
  10 : C               DEDU ADAP avec metrique
  11 : C               On implémente ici la première partie de la fonctionnelle
  12 : C               de Huang
  13 : C     
  14 : C LANGAGE     : ESOPE
  15 : C AUTEUR      : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF)
  16 : C               mél : gounand@semt2.smts.cea.fr
  17 : C***********************************************************************
  18 : C APPELES          :
  19 : C APPELE PAR       :
  20 : C***********************************************************************
  21 : C ENTREES            :
  22 : C ENTREES/SORTIES    :
  23 : C SORTIES            : -
  24 : C TRAVAIL            :
  25 : C***********************************************************************
  26 : C VERSION    : v1, 30/04/07, version initiale
  27 : C HISTORIQUE : v1, 30/04/07, création
  28 : C HISTORIQUE :
  29 : C***********************************************************************
  30 : C Prière de PRENDRE LE TEMPS de compléter les commentaires
  31 : C en cas de modification de ce sous-programme afin de faciliter
  32 : C la maintenance !
  33 : C***********************************************************************
  34 : -INC CCOPTIO
  35 : -INC CCREEL
  36 : CBEGININCLUDE SMCHAEL
  37 :       SEGMENT MCHAEL
  38 :       POINTEUR IMACHE(N1).MELEME
  39 :       POINTEUR ICHEVA(N1).MCHEVA
  40 :       ENDSEGMENT
  41 :       SEGMENT MCHEVA
  42 :       REAL*8 VELCHE(NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM)
  43 :       ENDSEGMENT
  44 :       SEGMENT LCHEVA
  45 :       POINTEUR LISCHE(NBCHE).MCHEVA
  46 :       ENDSEGMENT
  47 : CENDINCLUDE SMCHAEL
  48 :       INTEGER NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM,N1
  49 :       POINTEUR FC.MCHEVA
  50 :       POINTEUR LCOF.LCHEVA
  51 :       POINTEUR JMAJAC.MCHEVA
  52 :       POINTEUR JMIJAC.MCHEVA
  53 :       POINTEUR JDTJAC.MCHEVA
  54 :       POINTEUR JMAREG.MCHEVA
  55 :       POINTEUR JMET.MCHEVA
  56 :       POINTEUR JTHE.MCHEVA
  57 :       POINTEUR JGAM.MCHEVA
  58 :       CHARACTER*8 NOMLOI
  59 :       INTEGER ICOF
  60 : *     
  61 : -INC TMXMAT
  62 : * Objets temporaires
  63 :       POINTEUR JAC.MXMAT,JT.MXMAT
  64 :       POINTEUR G.MXMAT,IG.MXMAT,H.MXMAT,HIG.MXMAT,HM1.MXMAT
  65 :       POINTEUR MM1.MXMAT,JTM.MXMAT,MJ.MXMAT
  66 : * Objets à préconditionner (éventuellement)
  67 :       POINTEUR A.MXMAT,B.MXMAT,C.MXMAT,D.MXMAT,F.MXMAT,M.MXMAT
  68 : *     
  69 :       SEGMENT MCOF
  70 :       POINTEUR COEF(IDIM,IDIM).MCHEVA
  71 :       ENDSEGMENT
  72 :       POINTEUR MET.MCOF
  73 : *     
  74 :       LOGICAL LBID,LFIRST
  75 :       INTEGER LAXSP
  76 :       REAL*8 DEUPI,XR
  77 :       REAL*8 XL,XM,XH
  78 : *     
  79 :       INTEGER IMPR,IRET
  80 : *     
  81 : * Executable statements
  82 : *     
  83 :       IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans ccgahu'
  84 : C      IF (.NOT.(IDIM.EQ.1)) THEN
  85 : C         WRITE(IOIMP,*) 'IDIM=',IDIM,' ?'
  86 : C         GOTO 9999
  87 : C      ENDIF
  88 :       NLFC=FC.VELCHE(/6)
  89 :       NPFC=FC.VELCHE(/5)
  90 :       ICOF=0 
  91 : *     
  92 : *     Récupération des coefficients de la metrique
  93 : *     
  94 :       SEGINI MET
  95 :       DO IIDIM=1,IDIM
  96 :          ICOF=ICOF+1
  97 :          JMET=LCOF.LISCHE(ICOF)
  98 :          IF (ICOF.EQ.1) THEN 
  99 :             NLJM=JMET.VELCHE(/6)
 100 :             NPJM=JMET.VELCHE(/5)
 101 :          ELSE
 102 :             NLJM2=JMET.VELCHE(/6)
 103 :             NPJM2=JMET.VELCHE(/5)
 104 :             IF (NLJM2.NE.NLJM.OR.NPJM2.NE.NPJM) THEN 
 105 :                WRITE(IOIMP,*) 'Erreur grave dims JMET'
 106 :                GOTO 9999
 107 :             ENDIF
 108 :          ENDIF
 109 :          MET.COEF(IIDIM,IIDIM)=JMET
 110 :       ENDDO
 111 :       DO IIDIM=1,IDIM
 112 :          NJ=IDIM-IIDIM
 113 :          IF (NJ.GE.1) THEN
 114 :             DO JIDIM=IIDIM+1,IDIM
 115 :                ICOF=ICOF+1
 116 :                JMET=LCOF.LISCHE(ICOF)
 117 :                NLJM2=JMET.VELCHE(/6)
 118 :                NPJM2=JMET.VELCHE(/5)
 119 :                IF (NLJM2.NE.NLJM.OR.NPJM2.NE.NPJM) THEN 
 120 :                   WRITE(IOIMP,*) 'Erreur grave dims JMET2'
 121 :                   GOTO 9999
 122 :                ENDIF
 123 :                MET.COEF(IIDIM,JIDIM)=JMET
 124 :             ENDDO
 125 :          ENDIF
 126 :       ENDDO
 127 : *     
 128 :       ICOF=ICOF+1
 129 :       JTHE=LCOF.LISCHE(ICOF)
 130 :       NLJT=JTHE.VELCHE(/6)
 131 :       NPJT=JTHE.VELCHE(/5)
 132 :       I1  =JTHE.VELCHE(/4)
 133 :       I2  =JTHE.VELCHE(/3)
 134 :       IF ((I1.NE.1).OR.(I2.NE.1))
 135 :      $     THEN
 136 :          WRITE(IOIMP,*) 'Erreur dims JTHE'
 137 :          GOTO 9999
 138 :       ENDIF
 139 : *     
 140 :       ICOF=ICOF+1
 141 :       JGAM=LCOF.LISCHE(ICOF)
 142 :       NLJG=JTHE.VELCHE(/6)
 143 :       NPJG=JTHE.VELCHE(/5)
 144 :       I1  =JTHE.VELCHE(/4)
 145 :       I2  =JTHE.VELCHE(/3)
 146 :       IF ((I1.NE.1).OR.(I2.NE.1))
 147 :      $     THEN
 148 :          WRITE(IOIMP,*) 'Erreur dims JGAM'
 149 :          GOTO 9999
 150 :       ENDIF
 151 : *     
 152 :       ICOF=ICOF+1
 153 :       JMAJAC=LCOF.LISCHE(ICOF)
 154 :       NLJA=JMAJAC.VELCHE(/6)
 155 :       NPJA=JMAJAC.VELCHE(/5)
 156 :       IREF=JMAJAC.VELCHE(/4)
 157 :       IREL=JMAJAC.VELCHE(/3)
 158 : *     
 159 :       IF (IREL.NE.IDIM) THEN
 160 :          WRITE(IOIMP,*) 'Erreur dims JMAJAC'
 161 :          GOTO 9999
 162 :       ENDIF
 163 : *     
 164 :       ICOF=ICOF+1
 165 :       ICOF=ICOF+1
 166 :       ICOF=ICOF+1
 167 :       JMAREG=LCOF.LISCHE(ICOF)
 168 :       NLJR=JMAREG.VELCHE(/6)
 169 :       NPJR=JMAREG.VELCHE(/5)
 170 :       I1  =JMAREG.VELCHE(/4)
 171 :       I2  =JMAREG.VELCHE(/3)
 172 :       IF ((NLJR.NE.1).OR.(NPJR.NE.1).OR.(I1.NE.IREF).OR.(I2.NE.IREF))
 173 :      $     THEN
 174 :          WRITE(IOIMP,*) 'Erreur dims JMAREG'
 175 :          GOTO 9999
 176 :       ENDIF
 177 : *     
 178 : * Objets temporaires et à préconditionner
 179 : *     
 180 :       LDIM1=IREL
 181 :       LDIM2=IREF
 182 :       SEGINI,JAC
 183 :       SEGINI,MJ
 184 :       SEGINI,A
 185 :       SEGINI,B
 186 :       LDIM1=IREF
 187 :       LDIM2=IREL
 188 :       SEGINI,JT
 189 :       SEGINI,JTM
 190 :       LDIM1=IREF
 191 :       LDIM2=IREF
 192 :       SEGINI,G
 193 :       SEGINI,IG
 194 :       SEGINI,H
 195 :       SEGINI,HM1
 196 :       SEGINI,HIG
 197 :       SEGINI,C
 198 :       LDIM1=IREL
 199 :       LDIM2=IREL
 200 :       SEGINI,M
 201 :       SEGINI,MM1
 202 :       SEGINI,D
 203 :       SEGINI,F
 204 : *     
 205 : * Calcul de la métrique des éléments réguliers
 206 : *     
 206 :       CALL MAMA(JMAREG.VELCHE,IREF,IREF,
 207 :      $     'JTJ     ',H.XMAT,IREF,IREF,
 208 :      $     IMPR,IRET)
 209 :       IF (IRET.NE.0) GOTO 9999
 210 : *      SEGPRT,H
 211 :       LBID=.FALSE.
 211 :       CALL GEOLI2(IREF,1,1,H.XMAT,HM1.XMAT,DETH,LBID,IMPR,IRET)
 212 :       IF (IRET.NE.0) GOTO 9999
 213 : *!!!      XH=1.0D0
 214 :       XH=1.D0/(SQRT(DETH))
 215 : *!!      XH=1.D0/(DETH)
 216 : *     
 217 :       LFIRST = .FALSE.
 218 :       DO ILFC=1,NLFC
 219 :          IF (NLJM.EQ.1) THEN
 220 :             ILJM=1
 221 :          ELSE
 222 :             ILJM=ILFC
 223 :          ENDIF
 224 :          IF (NLJA.EQ.1) THEN
 225 :             ILJA=1
 226 :          ELSE
 227 :             ILJA=ILFC
 228 :          ENDIF
 229 :          IF (NLJT.EQ.1) THEN
 230 :             ILJT=1
 231 :          ELSE
 232 :             ILJT=ILFC
 233 :          ENDIF
 234 :          IF (NLJG.EQ.1) THEN
 235 :             ILJG=1
 236 :          ELSE
 237 :             ILJG=ILFC
 238 :          ENDIF
 239 :          DO IPFC=1,NPFC
 240 :             IF (NPJM.EQ.1) THEN
 241 :                IPJM=1
 242 :             ELSE
 243 :                IPJM=IPFC
 244 :             ENDIF
 245 :             IF (NPJA.EQ.1) THEN
 246 :                IPJA=1
 247 :             ELSE
 248 :                IPJA=IPFC
 249 :             ENDIF
 250 :             IF (NPJT.EQ.1) THEN
 251 :                IPJT=1
 252 :             ELSE
 253 :                IPJT=IPFC
 254 :             ENDIF
 255 :             IF (NPJG.EQ.1) THEN
 256 :                IPJG=1
 257 :             ELSE
 258 :                IPJG=IPFC
 259 :             ENDIF
 260 : *
 261 : * Récupération de l'exposant multiplié par 2
 262 : * On suppose que c'est un entier
 263 : *
 264 :             XTHE=JTHE.VELCHE(1,1,1,1,IPJT,ILJT)
 265 :             XGAM=JGAM.VELCHE(1,1,1,1,IPJG,ILJG)
 266 :             XCO1=XTHE*0.5D0
 267 :             IXP1=NINT(IREF*XGAM)
 268 :             XCO2=(1.D0-XTHE)*((SQRT(DBLE(IREF)))**IXP1)
 269 : *!!!
 270 : *!!!         XCO2=1.D0 
 271 :             IXP2=NINT(1-XGAM)
 272 :             IF (LFIRST) THEN
 273 :                WRITE(IOIMP,*) 'XTHE=',XTHE
 274 :                WRITE(IOIMP,*) 'XGAM=',XGAM
 275 :                WRITE(IOIMP,*) 'XCO1=',XCO1
 276 :                WRITE(IOIMP,*) 'XCO2=',XCO2
 277 :                WRITE(IOIMP,*) 'IXP1=',IXP1,'/2'
 278 :                WRITE(IOIMP,*) 'IXP2=',IXP2
 279 :                LFIRST=.FALSE.
 280 :             ENDIF
 281 : *     
 282 : * Copie des coefficients de la métrique
 283 : *     
 284 :             DO IIDIM=1,IDIM
 285 :                JMET=MET.COEF(IIDIM,IIDIM)
 286 :                M.XMAT(IIDIM,IIDIM)=JMET.VELCHE(1,1,1,1,IPJM,ILJM)
 287 :             ENDDO
 288 :             DO IIDIM=1,IIDIM
 289 :                NJ=IDIM-IIDIM
 290 :                IF (NJ.GE.1) THEN
 291 :                   DO JIDIM=IIDIM+1,IDIM
 292 :                      JMET=MET.COEF(IIDIM,JIDIM)
 293 :                      M.XMAT(IIDIM,JIDIM)=JMET.VELCHE(1,1,1,1,IPJM,ILJM)
 294 :                      M.XMAT(JIDIM,IIDIM)=JMET.VELCHE(1,1,1,1,IPJM,ILJM)
 295 :                   ENDDO
 296 :                ENDIF
 297 :             ENDDO
 298 : *            SEGPRT,M
 299 : * Calcul de l'inverse
 300 : *     
 301 :            LBID=.FALSE.
 301 :            CALL GEOLI2(IREL,1,1,M.XMAT,MM1.XMAT,DETM,LBID,IMPR,IRET)
 302 :            IF (IRET.NE.0) GOTO 9999
 303 :            XDM=1.D0/(SQRT(DETM))
 304 : *1           XDM=1.D0
 305 : ******           XDM=SQRT(DETM)
 306 : *     
 307 : * Copie du jacobien
 308 : *     
 308 :             CALL MAMA(JMAJAC.VELCHE(1,1,1,1,IPJA,ILJA),IREL,IREF,
 309 :      $           'COPIE   ',
 310 :      $           JAC.XMAT,IREL,IREF,
 311 :      $           IMPR,IRET)
 312 :             IF (IRET.NE.0) GOTO 9999
 313 : *            SEGPRT,JAC
 314 : *     
 315 : * Calcul de la métrique G
 316 : *     
 317 : *     Calcul de Jt
 317 :             CALL MAMA(JAC.XMAT,IREL,IREF,
 318 :      $           'TRANSPOS',JT.XMAT,IREF,IREL,
 319 :      $           IMPR,IRET)
 320 :             IF (IRET.NE.0) GOTO 9999
 321 : *     Calcul de MJ
 321 :             CALL MAMAMA(M.XMAT,IREL,IREL,JAC.XMAT,IREL,IREF,
 322 :      $           'FOIS    ',MJ.XMAT,IREL,IREF,
 323 :      $           IMPR,IRET)
 324 :             IF (IRET.NE.0) GOTO 9999
 325 : *     Calcul de JTM
 325 :             CALL MAMA(MJ.XMAT,IREL,IREF,
 326 :      $           'TRANSPOS',JTM.XMAT,IREF,IREL,
 327 :      $           IMPR,IRET)
 328 :             IF (IRET.NE.0) GOTO 9999
 329 : *     Calcul de G=JtMJ
 329 :             CALL MAMAMA(JT.XMAT,IREF,IREL,MJ.XMAT,IREL,IREF,
 330 :      $           'FOIS    ',G.XMAT,IREF,IREF,
 331 :      $           IMPR,IRET)
 332 :             IF (IRET.NE.0) GOTO 9999
 333 : * Calcul de l'inverse, du déterminant et trace de l'inverse de g
 334 :             LBID=.FALSE.
 334 :             CALL GEOLI2(IREF,1,1,G.XMAT,IG.XMAT,DETG,LBID,IMPR,IRET)
 335 :             IF (IRET.NE.0) GOTO 9999
 336 :             XM=SQRT(DETG)
 337 : *!            WRITE(IOIMP,*) 'XM=',XM
 338 : * Calcul de hg-1 et de sa trace
 338 :             CALL MAMAMA(H.XMAT,IREF,IREF,IG.XMAT,IREF,IREF,
 339 :      $           'FOIS    ',HIG.XMAT,IREF,IREF,IMPR,IRET)
 340 :             IF (IRET.NE.0) GOTO 9999
 340 :             CALL MARE(HIG.XMAT,IREF,IREF,'TRACE   ',
 341 :      $           XL,IMPR,IRET)
 342 :             IF (IRET.NE.0) GOTO 9999
 343 : * Calcul de A = (J+)t = MJg-1 (pseudo-inverse)
 343 :             CALL MAMAMA(MJ.XMAT,IREL,IREF,IG.XMAT,IREF,IREF,
 344 :      $           'FOIS    ',A.XMAT,IREL,IREF,IMPR,IRET)
 345 :             IF (IRET.NE.0) GOTO 9999
 346 : * Calcul de B = (J+)t hg-1 = MJg-1hg-1
 346 :             CALL MAMAMA(A.XMAT,IREL,IREF,HIG.XMAT,IREF,IREF,
 347 :      $           'FOIS    ',B.XMAT,IREL,IREF,IMPR,IRET)
 348 :             IF (IRET.NE.0) GOTO 9999
 349 : * Calcul de C = g-1hg-1
 349 :             CALL MAMAMA(IG.XMAT,IREF,IREF,HIG.XMAT,IREF,IREF,
 350 :      $           'FOIS    ',C.XMAT,IREF,IREF,IMPR,IRET)
 351 :             IF (IRET.NE.0) GOTO 9999
 352 : * Calcul de D = MJ g-1hg-1 JtM = B JtM
 352 :             CALL MAMAMA(B.XMAT,IREL,IREF,JTM.XMAT,IREF,IREL,
 353 :      $           'FOIS    ',D.XMAT,IREL,IREL,IMPR,IRET)
 354 :             IF (IRET.NE.0) GOTO 9999
 355 : * Calcul de F = MJ g-1 JtM = A JtM
 355 :             CALL MAMAMA(A.XMAT,IREL,IREF,JTM.XMAT,IREF,IREL,
 356 :      $           'FOIS    ',F.XMAT,IREL,IREL,IMPR,IRET)
 357 :             IF (IRET.NE.0) GOTO 9999
 358 : *     
 359 : * Calcul de la fonctionnelle
 360 : *     
 361 :             IF (NOMLOI.EQ.'AHUF    ') THEN
 362 :                XL1=SQRT(XL)**IXP1
 363 :                FONC1=XL1*XM
 364 : *
 365 :                XM1=XM**IXP2
 366 :                FONC2=XM1
 367 : *
 368 :                CONTRI=(XCO1*FONC1)+(XCO2*FONC2*(XH**(IXP2-1)))
 369 : *!!!!
 370 :                CONTRI=CONTRI*XDM
 371 : *     
 372 : * Calcul du résidu
 373 : *     
 374 :             ELSEIF (NOMLOI(1:4).EQ.'AHUR') THEN
 374 :                CALL CH2INT(NOMLOI(5:5),IDIM1,IMPR,IRET)
 375 :                IF (IRET.NE.0) GOTO 9999
 375 :                CALL CH2INT(NOMLOI(6:6),IDIM2,IMPR,IRET)
 376 :                IF (IRET.NE.0) GOTO 9999
 377 :                K=IDIM1
 378 :                L=IDIM2
 379 :                XL1=SQRT(XL)**IXP1
 380 :                XL2=IXP1*(SQRT(XL)**(IXP1-2))
 381 :                AKL1=XL1*XM*A.XMAT(K,L)
 382 :                AKL2=-1.D0*XL2*XM*B.XMAT(K,L)
 383 :                AKL=AKL1+AKL2
 384 : *
 385 :                XM2=IXP2*(XM**IXP2)
 386 :                CKL=XM2*A.XMAT(K,L)
 387 : *
 388 :                CONTRI=(XCO1*AKL)+(XCO2*CKL*(XH**(IXP2-1)))
 389 : *!!!!
 390 :                CONTRI=CONTRI*XDM
 391 : *     
 392 : * Calcul du jacobien
 393 : *     
 394 :             ELSEIF (NOMLOI(1:4).EQ.'AHUJ') THEN
 394 :                CALL CH2INT(NOMLOI(5:5),IDIM1,IMPR,IRET)
 395 :                IF (IRET.NE.0) GOTO 9999
 395 :                CALL CH2INT(NOMLOI(6:6),IDIM2,IMPR,IRET)
 396 :                IF (IRET.NE.0) GOTO 9999
 396 :                CALL CH2INT(NOMLOI(7:7),IDIM3,IMPR,IRET)
 397 :                IF (IRET.NE.0) GOTO 9999
 397 :                CALL CH2INT(NOMLOI(8:8),IDIM4,IMPR,IRET)
 398 :                IF (IRET.NE.0) GOTO 9999
 399 :                I=IDIM1
 400 :                J=IDIM2
 401 :                K=IDIM3
 402 :                L=IDIM4
 403 :                XL1=SQRT(XL)**IXP1
 404 :                XL2=IXP1*(SQRT(XL)**(IXP1-2))
 405 :                XL3=IXP1*(IXP1-2)*(SQRT(XL)**(IXP1-4))
 406 :                B11=XM*XL1*A.XMAT(I,J)*A.XMAT(K,L)
 407 :                B12=-1.D0*XL2*XM*B.XMAT(I,J)*A.XMAT(K,L)
 408 :                B13=XM*XL1*M.XMAT(K,I)*IG.XMAT(J,L)
 409 :                B14=-1.D0*XM*XL1*F.XMAT(K,I)*IG.XMAT(J,L)
 410 :                B15=-1.D0*XM*XL1*A.XMAT(K,J)*A.XMAT(I,L)
 411 :                BIJKL1=B11+B12+B13+B14+B15
 412 :                B20=+1.D0*XL3*XM*B.XMAT(I,J)*B.XMAT(K,L)
 413 :                B21=-1.D0*XL2*XM*A.XMAT(I,J)*B.XMAT(K,L)
 414 :                B22=-1.D0*XL2*XM*M.XMAT(K,I)*C.XMAT(J,L)
 415 :                B23=+1.D0*XL2*XM*F.XMAT(K,I)*C.XMAT(J,L)
 416 :                B24=+1.D0*XL2*XM*A.XMAT(K,J)*B.XMAT(I,L)
 417 :                B25=+1.D0*XL2*XM*D.XMAT(K,I)*IG.XMAT(J,L)
 418 :                B26=+1.D0*XL2*XM*B.XMAT(K,J)*A.XMAT(I,L)
 419 :                BIJKL2=B20+B21+B22+B23+B24+B25+B26
 420 :                BIJKL=BIJKL1+BIJKL2
 421 : *
 422 :                XM2=IXP2*(XM**IXP2)
 423 :                XM3=IXP2*IXP2*(XM**IXP2)
 424 :                D1=XM3*A.XMAT(I,J)*A.XMAT(K,L)
 425 :                D2=XM2*M.XMAT(K,I)*IG.XMAT(J,L)
 426 :                D3=-1.D0*XM2*F.XMAT(K,I)*IG.XMAT(J,L)
 427 :                D4=-1.D0*XM2*A.XMAT(K,J)*A.XMAT(I,L)
 428 :                DIJKL=D1+D2+D3+D4
 429 : *
 430 :                CONTRI=(XCO1*BIJKL)+(XCO2*DIJKL*(XH**(IXP2-1)))
 431 : *!!!!
 432 :                CONTRI=CONTRI*XDM
 433 : *!               WRITE(IOIMP,*) 'CONTRI=',CONTRI
 434 :             ELSE
 435 :                WRITE(IOIMP,*) 'Erreur grave'
 436 :                GOTO 9999
 437 :             ENDIF
 438 :             FC.VELCHE(1,1,1,1,IPFC,ILFC)=
 439 :      $           FC.VELCHE(1,1,1,1,IPFC,ILFC)+
 440 :      $           CONTRI
 441 :          ENDDO
 442 :       ENDDO
 443 :       SEGSUP,F
 444 :       SEGSUP,D
 445 :       SEGSUP,MM1
 446 :       SEGSUP,M
 447 :       SEGSUP,C
 448 :       SEGSUP,HIG
 449 :       SEGSUP,HM1
 450 :       SEGSUP,H
 451 :       SEGSUP,IG
 452 :       SEGSUP,G
 453 :       SEGSUP,JTM
 454 :       SEGSUP,JT
 455 :       SEGSUP,B
 456 :       SEGSUP,A
 457 :       SEGSUP,MJ
 458 :       SEGSUP,JAC
 459 :       SEGSUP,MET
 460 : *     
 461 : * Normal termination
 462 : *     
 463 :       IRET=0
 464 :       RETURN
 465 : *     
 466 : * Format handling
 467 : *     
 468 : *     
 469 : * Error handling
 470 : *     
 471 :  9999 CONTINUE
 472 :       IRET=1
 473 :       WRITE(IOIMP,*) 'An error was detected in subroutine ccgahu'
 474 :       RETURN
 475 : *     
 476 : * End of subroutine CCGAHU
 477 : *     
 478 :       END
 479 : 
 480 :  
 481 :  
 482 :  
 483 :