1 : C EXDIAG SOURCE GOUNAND 12/02/01 21:15:02 7265
2 : SUBROUTINE EXDIAG(IOPT)
3 : IMPLICIT REAL*8 (A-H,O-Z)
4 : IMPLICIT INTEGER (I-N)
5 : C*************************************************************************
6 : C Operateur EXDIAG
7 : C
8 : C OBJET : Extrait une diagonale d'une matrice au format
9 : C MATRIK
10 : C Pour cela, on doit l'assembler.
11 : C L'assemblage se fait comme dans KRES2.
12 : C Si IOPT=1, on renvoie la diagonale de la matrice.
13 : C Si IOPT=2, on renvoie un inverse approché (SPAI)
14 : C diagonal de la matrice entrée.
15 : C Son expression est :
16 : C m_ii = a_ii / \sum_{j=1,n} aij^2
17 : C
18 : C Voir aussi la notice de KOPS
19 : C
20 : C***********************************************************************
21 : C HISTORIQUE : 21/03/08 : Première version
22 : C
23 : C***********************************************************************
24 : C Prière de PRENDRE LE TEMPS de compléter les commentaires
25 : C en cas de modification de ce sous-programme afin de faciliter
26 : C la maintenance !
27 : C***********************************************************************
28 : -INC CCOPTIO
29 : -INC CCREEL
30 : -INC SMCHPOI
31 : POINTEUR MCHDIA.MCHPOI
32 : -INC SMMATRIK
33 : POINTEUR MATRI2.MATRIK
34 : POINTEUR AMORS.PMORS
35 : POINTEUR AISA.IZA
36 : POINTEUR IDIAGO.IZA
37 : POINTEUR ISCAR.IZA
38 : INTEGER IMPR,IRET
39 : C
40 : CHARACTER*4 MRENU,MMULAG
41 : CHARACTER*8 TYPE
42 : CHARACTER*8 TYMATK,TYRIGI,BLAN
43 : DATA TYMATK,TYRIGI,BLAN/'MATRIK ','RIGIDITE',' '/
44 : *
45 : IMPR=0
46 : *
47 : * Lecture de la matrice
48 : *
49 : TYPE=BLAN
49 : CALL QUETYP(TYPE,1,IRET)
50 : IF (IRET.EQ.0) GOTO 9999
51 : IF (TYPE.EQ.TYRIGI) THEN
52 : IF (IOPT.EQ.2) THEN
53 : * Transformation en matrik en changement de noms d'inconnues
53 : CALL RIMA
54 : IF (IERR.NE.0) RETURN
54 : CALL MACHI2(1)
55 : IF (IERR.NE.0) RETURN
56 : ELSE
56 : CALL LIROBJ('RIGIDITE',IRIG,1,IRET)
57 : IF (IERR.NE.0) RETURN
58 : *
58 : CALL EXDIAR(IRIG,ICHP)
59 : IF (IERR.NE.0) RETURN
60 : *
60 : CALL ECROBJ('CHPOINT ',ICHP)
61 : RETURN
62 : ENDIF
63 : ENDIF
64 : TYPE=TYMATK
64 : CALL LIROBJ(TYPE,MATRIK,1,IRET)
65 : IF(IRET.EQ.0) GOTO 9999
66 : C
67 : C Assemblage proprement dit
68 : C
69 : C Attention, on effectue une recopie du chapeau pour ne garder
70 : C aucune information de préconditionnement (assemblage, numérotation)
71 : C dans la matrice originale sinon une résolution subséquente poserait
72 : C problème !!!!!!
73 : C
74 : SEGINI,MATRI2=MATRIK
75 : MATRIK=MATRI2
76 : MATASS=MATRIK
77 : MRENU='RIEN'
78 : MMULAG='RIEN'
79 : METASS=4
79 : CALL KRES3(MATRIK,MATASS,MRENU,MMULAG,METASS,
80 : $ 0,.FALSE.,
81 : $ IMPR,IRET)
82 : IF (IRET.NE.0) GOTO 9999
83 : C
84 : C Extraction de la diagonale
85 : C
86 : SEGACT MATRIK
87 : AMORS=MATRIK.KIDMAT(4)
88 : AISA=MATRIK.KIDMAT(5)
89 : *
90 : SEGACT AMORS
91 : SEGACT AISA
92 : NTTDDL=AMORS.IA(/1)-1
93 : NBVA=NTTDDL
94 : SEGINI IDIAGO
95 : IF (IOPT.EQ.2) THEN
96 : SEGINI ISCAR
97 : ENDIF
98 : DO ITTDDL=1,NTTDDL
99 : JSTRT=AMORS.IA(ITTDDL)
100 : JSTOP=AMORS.IA(ITTDDL+1)-1
101 : DO J=JSTRT,JSTOP
102 : JTTDDL=AMORS.JA(J)
103 : IF (JTTDDL.EQ.ITTDDL) THEN
104 : IDIAGO.A(ITTDDL)=AISA.A(J)
105 : ENDIF
106 : IF (IOPT.EQ.2) THEN
107 : VAL=AISA.A(J)
108 : ISCAR.A(ITTDDL)=ISCAR.A(ITTDDL)+(VAL*VAL)
109 : ENDIF
110 : ENDDO
111 : ENDDO
112 : SEGSUP AISA
113 : SEGSUP AMORS
114 : *
115 : IF (IOPT.EQ.2) THEN
116 : DO ITTDDL=1,NTTDDL
117 : VAL=ISCAR.A(ITTDDL)
118 : IF (VAL.LE.SQRT(XPETIT)) THEN
119 : WRITE(IOIMP,*) 'La ligne ',ITTDDL,
120 : $ ' de la matrice est nulle : ', VAL
121 : GOTO 9999
122 : ENDIF
123 : IDIAGO.A(ITTDDL)=IDIAGO.A(ITTDDL)/VAL
124 : ENDDO
125 : ENDIF
126 : IF (IOPT.EQ.2) THEN
127 : SEGSUP ISCAR
128 : ENDIF
129 : C
130 : C Transformation en chpoint
131 : C
131 : CALL XDISP(MATRIK,IDIAGO,MCHDIA,IMPR,IRET)
132 : IF (IRET.NE.0) GOTO 9999
133 : SEGSUP IDIAGO
134 : SEGSUP MATRIK
134 : CALL ECROBJ('CHPOINT ',MCHDIA)
135 : *
136 : * Normal termination
137 : *
138 : RETURN
139 : *
140 : * Format handling
141 : *
142 : *
143 : * Error handling
144 : *
145 : 9999 CONTINUE
146 : WRITE(IOIMP,*) 'An error was detected in exdiag.eso'
147 : * 153 2
148 : * Opération illicite dans ce contexte
148 : CALL ERREUR(153)
149 : RETURN
150 : *
151 : * End of EXDIAG
152 : *
153 : END
154 :
155 :
156 :
157 :
158 :
159 :
160 :
161 :
162 :
163 :
164 :
165 :
166 :
167 :
168 :
169 :
170 :
171 :
172 :
173 :
174 :
© Cast3M 2003 - All rights reserved.
Disclaimer