Télécharger la plaquette Cast3M


fr en

Télécharger arpitl.eso

Retour à la liste

Numérotation des lignes : 

arpitl

  1. C ARPITL    SOURCE    PV        22/04/15    13:20:06     11344          
  2.       SUBROUTINE ARPITL (IPRTRA,IPMAUP,SIGMA,INVER,SYM,EPSI,OUT)
  3.  
  4.  
  5. ***********************************************************************
  6. *
  7. *                          A R P I T L
  8. *
  9. * FONCTION:
  10. * ---------
  11. *
  12. *     STEP DE LA FACTORISATION D'ARNOLDI POUR UN PROBLEME LINEAIRE.
  13. *
  14.  
  15. * REMARQUE:
  16. * ---------
  17. *
  18. *     ON NOTE:
  19. *
  20. *     A=IPRIGI
  21. *     B=IPMASS
  22. *
  23. *     IPRIGI : RIGIDITE
  24. *     IPMASS : MASSE OU KSIGMA
  25. *
  26. *
  27. * PARAMETRES:  (E)=ENTREE   (S)=SORTIE
  28. * -----------
  29. *
  30. *
  31. *     IPRTRA  ENTIER     (E)    OPERATEURS DE TRAVAIL
  32. *
  33. *     IPMAUP  ENTIER     (E/S)  POINTEUR VARIABLES ARPACK
  34. *
  35. *     SIGMA   COMPLEX DP (E)    VALEUR SU SHIFT (PEUT ETRE NUL)
  36. *
  37. *     INVER   LOGIQUE    (E)   .TRUE. -> PRODUIT SCALAIRE X'KX
  38. *                              .FALSE. -> PRODUIT SCALAIRE X'MX
  39. *
  40. *     SYM     LOGIQUE    (E)    PROBLEME SYMETRIQUE OU NON
  41. *
  42. *     EPSI    REEL DP    (E)    TOLERANCE EIGENPAIRS
  43. *
  44. *     OUT     LOGIQUE    (S)    FLAG DE CONVERGENCE
  45. *
  46. *
  47. * SOUS-PROGRAMMES APPELES:
  48. * ------------------------
  49. *
  50. *     DSAUPD,DNAUPD,ARPCH1,MUCPRI,RESOU1,LDMT,DECALE,DTCHPO,ARPERR
  51. *
  52. * AUTEUR,DATE DE CREATION:
  53. * -------------------------
  54. *
  55. *     PASCAL BOUDA     29 JUIN 2015
  56. *
  57. * LANGAGE:
  58. * --------
  59. *
  60. *     FORTRAN 77 & 90
  61. *
  62. ***********************************************************************
  63.  
  64.       IMPLICIT INTEGER(I-N)
  65.       IMPLICIT REAL*8 (A-H,O-Z)
  66.  
  67.  
  68. -INC PPARAM
  69. -INC CCOPTIO
  70. -INC CCHAMP
  71. -INC SMRIGID
  72. -INC TARWORK
  73. -INC CCREEL
  74. c -INC TARTRAK
  75.  
  76.       SEGMENT idemen(0)
  77.  
  78.       INTEGER IPRTRA
  79.       INTEGER IPMAUP
  80.       COMPLEX*16 SIGMA
  81.       LOGICAL INVER
  82.       LOGICAL SYM
  83.       REAL*8 EPSI
  84.       LOGICAL OUT
  85.  
  86.  
  87.       INTEGER IPRIGI,IPMASS,IPKSIM
  88.       INTEGER TEST
  89.       CHARACTER*1 SCAL
  90.       INTEGER OPT
  91.       INTEGER IPCTRA(4)
  92.       INTEGER NOID,NOEN
  93.       INTEGER ndim,ncv,lworkl
  94.  
  95.       xspetl = xspeti
  96.       SEGINI IDEMEN
  97.       IDEMEN(**)=0
  98.       NOID=0
  99.       NOEN=1
  100.  
  101.       OUT=.FALSE.
  102.  
  103.       MAUP=IPMAUP
  104.       SEGACT MAUP*MOD
  105.  
  106.       MRITRA=IPRTRA
  107.       SEGACT MRITRA
  108.  
  109. *Recuperation des operateurs de travail
  110.       IPRIGI=RIGI(1)
  111.       IPMASS=RIGI(2)
  112.       IPKSIM=RIGI(4)
  113.  
  114.  
  115. *Récupération de la dimension des tableaux
  116.       ndim=resid(/1)
  117.       ncv=v(/2)
  118.       lworkl=workl(/1)
  119.  
  120. *Si le probleme est symétrique, on appelle la routine spécifique aux
  121. *problemes symetriques, sinon on appelle celle pour les problemes
  122. *quelconques
  123.  
  124.       IF (SYM) THEN
  125.         CALL DSAUPD (ido,bmat,ndim,which,nev,EPSI,resid,
  126.      &      ncv,v,ldv,iparam,ipntr,workd,workl,lworkl,ITRAK,info)
  127.  
  128.       ELSE
  129.  
  130.         CALL DNAUPD (ido,bmat,ndim,which,nev,EPSI,resid,
  131.      &      ncv,v,ldv,iparam,ipntr,workd,workl,lworkl,ITRAK,info)
  132.  
  133.       ENDIF
  134.  
  135. *Reverse communication: On récupère les paramètres de sortie et on
  136. *effectue des actions en fonction de leurs valeurs
  137.       TEST=ido
  138.       SCAL=bmat
  139.       OPT=iparam(7)
  140.  
  141.       IPMAUP=MAUP
  142. c       SEGDES MAUP
  143.  
  144. *On verifie si on doit stopper le processus
  145.       CALL ARPERR (IPMAUP,SYM,OUT)
  146.       IF (OUT) RETURN
  147.  
  148.  
  149. *Initialisation des chpoints de travail
  150.       DO i=1,4
  151.         IPCTRA(i)=0
  152.       ENDDO
  153.  
  154.  
  155. *SCAL: type de probleme
  156. *'I' si standard
  157. *'G' si generalise
  158.  
  159.       IF (SCAL .EQ. 'I') THEN
  160.  
  161.           IF (IIMPI.GE.10) WRITE(IOIMP,*) '*   PB AUX V.P. STANDARD  *'
  162.  
  163.           IF (TEST .EQ. -1 .OR. TEST .EQ. 1) THEN
  164.  
  165. *        &---------------------------------------------------&
  166. *        |          Calcul du produit matrice vecteur        |
  167. *        |           Y <---- inv(inv(B)*A-SIGMA*I)*X         |
  168. *        |                                                   |
  169. *        |               X : workd(ipntr(1))                 |
  170. *        |               Y : workd(ipntr(2))                 |
  171. *        &---------------------------------------------------&
  172.  
  173. ************************************************************************
  174. * 28/08/2015 : Dans ce cas, le shift est obligatoirement nul
  175. * decalage spectral avec une matrice identite non implemente
  176. ************************************************************************
  177.  
  178.             CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(3),1,3)
  179.  
  180.             CALL MUCPRI (IPCTRA(3),IPMASS,IPCTRA(2))
  181.  
  182.             CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(2),1,2)
  183.  
  184. *Mise a sero des inconnues en FLX
  185.             CALL ARCORC (IPCTRA(2),10)
  186. *Mise a zero des inconnues en PI  inutile ?
  187.             CALL ARCORC (IPCTRA(2),15)
  188.  
  189.             IDEMEN(1)=IPCTRA(2)
  190.  
  191.             IF (SYME(1) .EQ. 0) THEN
  192.               CALL RESOU1 (IPRIGI,IDEMEN,NOID,NOEN,xspetl,0,0)
  193.             ELSEIF (SYME(1) .EQ. 1) THEN
  194.               CALL LDMT (IPRIGI,IDEMEN,NOID,NOEN,xspetl,0)
  195.             ENDIF
  196.             IF(IERR.NE.0) RETURN
  197.  
  198.             IPCTRA(1)=IDEMEN(1)
  199.  
  200.             CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(1),2,1)
  201.  
  202.  
  203.           ENDIF
  204.  
  205.       ELSEIF (SCAL .EQ. 'G') THEN
  206.  
  207.         IF (IIMPI.GE.10) 
  208.      &  WRITE(IOIMP,*) '*   PB AUX V.P. GENERALISE  *',TEST,INVER
  209.  
  210.  
  211.         IF (TEST .EQ. -1)  THEN
  212.  
  213. *        &--------------------------------------------------&
  214. *        |          Calcul du produit matrice vecteur       |
  215. *        |                                                  |
  216. *        |        Y <---- inv(A-SIGMA*B)*B*X                |
  217. *        |                                                  |
  218. *        |              X : workd(ipntr(1))                 |
  219. *        |              Y : workd(ipntr(2))                 |
  220. *        &--------------------------------------------------&
  221.  
  222. c             WRITE(*,*) 'X1 :'
  223.             CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(1),1,3)
  224.  
  225.             CALL MUCPRI (IPCTRA(1),IPMASS,IPCTRA(3))
  226. *Mise a sero des inconnues en FLX
  227.             CALL ARCORC (IPCTRA(3),10)
  228. *Mise a zero des inconnues en PI inutile ?
  229.             CALL ARCORC (IPCTRA(3),15)
  230.  
  231. c             WRITE(*,*) '{B*X1} :'
  232.             CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(3),1,2)        
  233.  
  234. c             WRITE(*,*) 'avant RESOU :',SYME(4)
  235.             IDEMEN(1)=IPCTRA(3)
  236.             IF (SYME(4) .EQ. 0) THEN
  237.               CALL RESOU1 (IPKSIM,IDEMEN,NOID,NOEN,xspetl,0,0)
  238.             ELSEIF (SYME(4) .EQ. 1) THEN
  239.               CALL LDMT (IPKSIM,IDEMEN,NOID,NOEN,xspetl,0)
  240.             ENDIF
  241.             IF (IERR.NE.0) RETURN
  242.  
  243. c             WRITE(*,*) 'Y1=[OP^-1]*{B*X1} :'
  244.             IPCTRA(2)=IDEMEN(1)
  245. cbp            CALL ARPCH1 (IPKSIM,IPRIGI,IPMAUP,IPCTRA(2),2,1)
  246.             CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(2),2,1)
  247.  
  248.  
  249.         ELSEIF (TEST .EQ. 1) THEN
  250.  
  251. *        &--------------------------------------------------&
  252. *        |        Calcul du produit matrice vecteur         |
  253. *        |                                                  |
  254. *        |        si INVER :                                |
  255. *        |        Y <---- inv(A-SIGMA*B)*B*X                |
  256. *        |                                                  |
  257. *        |              X : workd(ipntr(1))                 |
  258. *        |              Y : workd(ipntr(2))                 |
  259. *        |                                                  |
  260. *        |        sinon :                                   |
  261. *        |        Y <---- inv(A-SIGMA*B)*X                  |
  262. *        |                                                  |
  263. *        |              X : workd(ipntr(3))                 |
  264. *        |              Y : workd(ipntr(2))                 |
  265. *        &--------------------------------------------------&
  266.  
  267.           IF (INVER) THEN
  268.  
  269.               CALL ARPCH1(IPRIGI,IPRIGI,IPMAUP,IPCTRA(1),1,3)
  270.  
  271.               CALL MUCPRI (IPCTRA(1),IPMASS,IPCTRA(3))
  272. *Mise a sero des inconnues en FLX
  273.               CALL ARCORC (IPCTRA(3),10)
  274. *Mise a zero des inconnues en PI
  275.               CALL ARCORC (IPCTRA(3),15)
  276.  
  277.               CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(3),1,2)
  278.  
  279.           ELSE
  280.  
  281. c               WRITE(*,*) 'X2 :'
  282. cbp              CALL ARPCH1 (IPKSIM,IPRIGI,IPMAUP,IPCTRA(3),3,4)
  283.               CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(3),3,4)
  284.  
  285. *Mise a sero des inconnues en FLX
  286.               CALL ARCORC (IPCTRA(3),10)
  287. *Mise a zero des inconnues en PI
  288.               CALL ARCORC (IPCTRA(3),15)
  289. c               WRITE(*,*) '{X2} : uniquement le chpoint '
  290.  
  291.           ENDIF
  292.  
  293.           IDEMEN(1)=IPCTRA(3)
  294.  
  295.           IF (SYME(4) .EQ. 0) THEN
  296.             CALL RESOU1 (IPKSIM,IDEMEN,NOID,NOEN,xspetl,0,0)
  297.           ELSEIF (SYME(4) .EQ. 1) THEN
  298.             CALL LDMT (IPKSIM,IDEMEN,NOID,NOEN,xspetl,0)
  299.           ENDIF
  300.           IF (IERR.NE.0) RETURN
  301.  
  302. c           WRITE(*,*) 'Y2=[OP^-1]*{X2} :'
  303.           IPCTRA(2)=IDEMEN(1)
  304. cbp          CALL ARPCH1 (IPKSIM,IPRIGI,IPMAUP,IPCTRA(2),2,1)
  305.           CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(2),2,1)
  306.  
  307.  
  308.         ELSEIF (TEST .EQ. 2) THEN
  309.  
  310. *        &-------------------------------------&
  311. *        |  Calcul du produit matrice vecteur  |
  312. *        |                                     |
  313. *        |    Si INVER                         |
  314. *        |        Y <---- A*X                  |
  315. *        |                                     |
  316. *        |    Sinon                            |
  317. *        |        Y <---- B*X                  |
  318. *        |                                     |
  319. *        |         X : workd(ipntr(1))         |
  320. *        |         Y : workd(ipntr(2))         |
  321. *        &-------------------------------------&
  322.  
  323.             CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(1),1,3)
  324.  
  325.             IF (.NOT. INVER) THEN
  326.  
  327.               CALL MUCPRI (IPCTRA(1),IPMASS,IPCTRA(2))
  328. c               WRITE(*,*) 'Y3=B*X3 :'
  329.  
  330.             ELSE
  331.  
  332.               CALL MUCPRI (IPCTRA(1),IPRIGI,IPCTRA(2))
  333.  
  334. *Mise a sero des inconnues en FLX
  335.               CALL ARCORC (IPCTRA(2),10)
  336. *Mise a zero des inconnues en PI
  337.               CALL ARCORC (IPCTRA(2),15)
  338. c               WRITE(*,*) 'Y3={B*X3} :'
  339.             ENDIF
  340.  
  341.             CALL ARPCH1 (IPRIGI,IPRIGI,IPMAUP,IPCTRA(2),2,2)
  342.  
  343.  
  344.           ENDIF
  345.  
  346.       ENDIF
  347.  
  348. *Destruction des chpoints de travail
  349.         DO i=1,4
  350.           IF (IPCTRA(i) .NE. 0) THEN
  351.             CALL DTCHPO (IPCTRA(i))
  352.           ENDIF
  353.         ENDDO
  354.  
  355.       SEGDES MRITRA
  356.  
  357.       END
  358.  
  359.  
  360.  
  361.  
  362.  
  363.  
  364.  
  365.  
  366.  
  367.  
  368.  

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