1 : C PSILU0 SOURCE GOUNAND 11/07/21 21:15:47 7046
2 : SUBROUTINE PSILU0(TRANS,N,NNZ,JU,JLU,ALU,Z,R)
3 : IMPLICIT INTEGER(I-N)
4 : IMPLICIT REAL*8 (A-H,O-Z)
5 : C***********************************************************************
6 : C NOM : PSILU0
7 : C DESCRIPTION :
8 : C PSILU0 solves the linear system Mz = r
9 : C M is a (preconditionning) matrix written as :
10 : C M=(D+L)D-1(D+U)
11 : C It can be the ILU(0) (i.e. Crout or Choleski)
12 : C incomplete factorization of a matrix A (computed in meilu0.eso)
13 : C or any other (preconditionning) matrix.
14 : C
15 : C WE DO NOT STORE M but WE STORE L, D-1 and U in the same
16 : C Modified Sparse Row (MSR) matrix.
17 : C (see reference Sparskit version 2 for a description)
18 : C
19 : C L and U correspond to the strictly lower and strictly
20 : C upper part of the factorization.
21 : C
22 : C
23 : C LANGAGE : FORTRAN 77 + chouilla ESOPE (pour les E/S)
24 : C AUTEUR : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/TTMF)
25 : C mél : gounand@semt2.smts.cea.fr
26 : C REFERENCE (bibtex-like) :
27 : C @BOOK{templates,
28 : C AUTHOR={R.Barrett, M.Berry, T.F.Chan, J.Demmel, J.Donato,
29 : C J.Dongarra, V.Eijkhout, R.Pozo, C.Romine,
30 : C H. Van der Vorst},
31 : C TITLE={Templates for the Solution of Linear Systems :
32 : C Building Blocks for Iterative Methods},
33 : C PUBLISHER={SIAM}, YEAR={1994}, ADDRESS={Philadelphia,PA} }
34 : C -> URL : http://www.netlib.org/templates/Templates.html
35 : C Sparskit : a basic tool kit for sparse matrix computations
36 : C Version 2 (Youcef Saad)
37 : C -> URL : http://www.cs.umn.edu/Research/arpa/SPARSKIT/sparskit.html
38 : C
39 : C***********************************************************************
40 : C APPELES : -
41 : C APPELE PAR : CGILU0, BCGSI0, GMRI0
42 : C***********************************************************************
43 : C ENTREES : N, NNZ, JU, JLU, ALU
44 : C R
45 : C ENTREES/SORTIES : Z
46 : C SORTIES : -
47 : C CODE RETOUR (IRET) : -
48 : C N : nombre de degrés de liberté
49 : C NNZ : nombre de valeurs non nulles de la factorisation Morse
50 : C JU,JLU,ALU : factorisation incomplète au format MSR
51 : C R : vecteur second membre du système à résoudre.
52 : C Z : vecteur des inconnues du système à résoudre.
53 : C initialisé dans la subroutine appelante.
54 : C***********************************************************************
55 : C VERSION : 20/12/99
56 : C HISTORIQUE : v1, 01/04/98, création
57 : C HISTORIQUE : 20/12/99
58 : C Modifications de la montée-descente (factorisations incomplètes
59 : C stockées au format MSR (Modified Sparse Row).
60 : C HISTORIQUE :
61 : C***********************************************************************
62 : C Prière de PRENDRE LE TEMPS de compléter les commentaires
63 : C en cas de modification de ce sous-programme afin de faciliter
64 : C la maintenance !
65 : C***********************************************************************
66 : *
67 : * .. Scalar Arguments ..
68 : INTEGER N,NNZ
69 : CHARACTER*1 TRANS
70 : LOGICAL LSAME
71 : EXTERNAL LSAME
72 : * ..
73 : * .. Array Arguments ..
74 : * .. Matrices stockées en Morse
75 : INTEGER JU(N)
76 : INTEGER JLU(NNZ+1)
77 : REAL*8 ALU(NNZ+1)
78 : * .. Vecteurs
79 : REAL*8 Z(N),R(N)
80 : *
81 : *
82 : * .. Variables locales
83 : INTEGER I,K
84 : * .. Executable Statements ..
85 : *
86 : IF (.NOT.LSAME(TRANS,'T')) THEN
87 : *
88 : * Montée
89 : *
90 : DO 1 I=1,N
91 : Z(I)=R(I)
92 : DO 12 K=JLU(I),JU(I)-1
93 : Z(I)=Z(I)-ALU(K)*Z(JLU(K))
94 : 12 CONTINUE
95 : 1 CONTINUE
96 : *
97 : * Descente
98 : *
99 : DO 3 I=N,1,-1
100 : DO 32 K=JU(I),JLU(I+1)-1
101 : Z(I)=Z(I)-ALU(K)*Z(JLU(K))
102 : 32 CONTINUE
103 : Z(I)=ALU(I)*Z(I)
104 : 3 CONTINUE
105 : ELSE
106 : *
107 : * Montée avec U^t
108 : *
109 : DO 5 I=1,N
110 : Z(I)=R(I)
111 : 5 CONTINUE
112 : *
113 : DO 7 I=1,N
114 : Z(I)=Z(I)*ALU(I)
115 : DO 72 K=JU(I),JLU(I+1)-1
116 : Z(JLU(K))=Z(JLU(K))-ALU(K)*Z(I)
117 : 72 CONTINUE
118 : 7 CONTINUE
119 : *
120 : * Descente avec L^t
121 : *
122 : DO 9 I=N,1,-1
123 : DO 92 K=JLU(I),JU(I)-1
124 : Z(JLU(K))=Z(JLU(K))-ALU(K)*Z(I)
125 : 92 CONTINUE
126 : 9 CONTINUE
127 : ENDIF
128 :
129 : *
130 : * Normal termination (if lucky !)
131 : *
132 : RETURN
133 : *
134 : * End of PSILU0.
135 : *
136 : END
137 :
138 :
139 :
140 :
141 :
© Cast3M 2003 - All rights reserved.
Disclaimer