File: sooc_panel_piv.F

package info (click to toggle)
mumps 5.1.2-5
  • links: PTS, VCS
  • area: main
  • in suites: sid
  • size: 15,704 kB
  • sloc: fortran: 310,672; ansic: 12,364; xml: 521; makefile: 469
file content (363 lines) | stat: -rw-r--r-- 12,146 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
C
C  This file is part of MUMPS 5.1.2, released
C  on Mon Oct  2 07:37:01 UTC 2017
C
C
C  Copyright 1991-2017 CERFACS, CNRS, ENS Lyon, INP Toulouse, Inria,
C  University of Bordeaux.
C
C  This version of MUMPS is provided to you free of charge. It is
C  released under the CeCILL-C license:
C  http://www.cecill.info/licences/Licence_CeCILL-C_V1-en.html
C
C     This file contains routines related to OOC,
C     panels, and pivoting. They are used to store
C     permutation information of what is already on
C     disk to be able to permute things back at the
C     solve stage.
C     They do not need to be in the MUMPS_OOC
C     module (most of them do not use any variable
C     from the module, or are called from routines
C     where we do not necessarily want to do a
C     USE SMUMPS_OOC).
      INTEGER FUNCTION SMUMPS_OOC_GET_PANEL_SIZE
     &     ( HBUF_SIZE, NNMAX, K227, K50 )
      IMPLICIT NONE
C
C     Arguments:
C     =========
C
      INTEGER, INTENT(IN) :: NNMAX, K227, K50
      INTEGER(8), INTENT(IN) :: HBUF_SIZE
C
C     Purpose:
C     =======
C
C     - Compute the effective size (maximum number of pivots in a panel)
C     for a front with NNMAX entries in its row (for U) /
C     column (for L).
C     - Be able to adapt the fixed number of columns in panel
C     depending on NNMAX, and size of IO buffer HBUF_SIZE
C
C     Local variables
C     ===============
C
      INTEGER K227_LOC
      INTEGER NBCOL_MAX
      INTEGER EFFECTIVE_SIZE
      NBCOL_MAX=int(HBUF_SIZE / int(NNMAX,8))
C     KEEP(227): Maximum size (nb of col/row) of a panel
      K227_LOC = abs(K227)
      IF (K50.EQ.2) THEN
C        for 2x2 pivots we may end-up having the first part
C        of a 2x2 pivot in the last col of the panel; the
C        adopted solution consists in adding the next column
C        to the panel; therefore we need be able to
C        dynamically increase the panel size by one.
C        note that we also maintain property:
C        KEEP(227): Maximum size (nb of col/row) of a panel
         K227_LOC=max(K227_LOC,2)
         EFFECTIVE_SIZE =  min(NBCOL_MAX-1, K227_LOC-1)
cN       - during bwd the effective size is useless
      ELSE
C        complete buffer space can be used for a panel
         EFFECTIVE_SIZE =  min(NBCOL_MAX, K227_LOC)
      ENDIF
      IF (EFFECTIVE_SIZE.LE.0) THEN
         write(6,*) 'Internal buffers too small to store ',
     &        ' ONE col/row of size', NNMAX
         CALL MUMPS_ABORT()
      ENDIF
      SMUMPS_OOC_GET_PANEL_SIZE = EFFECTIVE_SIZE
      RETURN
      END FUNCTION SMUMPS_OOC_GET_PANEL_SIZE
C
      SUBROUTINE SMUMPS_PERMUTE_PANEL( IPIV, LPIV, ISHIFT,
     &     THE_PANEL, NBROW, NBCOL, KbeforePanel )
      IMPLICIT NONE
C
C     Purpose:
C     =======
C
C     Permute rows of a panel, stored by columns, according
C     to permutation array IPIV.
C     IPIV is such that, for I = 1 to LPIV, row ISHIFT + I
C     in the front must be permuted with row IPIV( I )
C
C     Since the panel is not necessary at the beginning of
C     the front, let KbeforePanel be the number of pivots in the
C     front before the first pivot of the panel.
C
C     In the panel, row ISHIFT+I-KbeforePanel is permuted with
C     row IPIV(I)-KbeforePanel
C
C     Note:
C     ====
C
C     This routine can also be used to permute the columns of
C     a matrix (U) stored by rows. In that case, the argument
C     NBROW represents the number of columns, and NBCOL represents
C     the number of rows.
C
C
C     Arguments:
C     =========
C
      INTEGER LPIV, ISHIFT, NBROW, NBCOL, KbeforePanel
      INTEGER IPIV(LPIV)
      REAL THE_PANEL(NBROW, NBCOL)
C
C     Local variables:
C     ===============
C
      INTEGER I, IPERM
C
C     Executable statements
C     =====================
C
      DO I = 1, LPIV
C        Swap rows ISHIFT + I and PIV(I)
         IPERM=IPIV(I)
         IF ( I+ISHIFT.NE.IPERM) THEN
            CALL sswap(NBCOL,
     &           THE_PANEL(I+ISHIFT-KbeforePanel,1), NBROW,
     &           THE_PANEL(IPERM-KbeforePanel,1), NBROW)
         ENDIF
      END DO
      RETURN
      END SUBROUTINE SMUMPS_PERMUTE_PANEL
      SUBROUTINE SMUMPS_GET_OOC_PERM_PTR(TYPEF,
     &     NBPANELS,
     &     I_PIVPTR, I_PIV, IPOS, IW, LIW)
      USE MUMPS_OOC_COMMON ! To access TYPEF_L and TYPEF_U
      IMPLICIT NONE
      INCLUDE 'mumps_headers.h'
C
C     Purpose:
C     =======
C
C     Get the pointers in IW on pivoting information to be stored
C     during factorization and used during the solve phase. This
C     routine is both for the symmetric (TYPEF=TYPEF_L) and unsymmetric
C     cases (TYPEF=TYPEF_L or TYPEF_U).
C     The total size of this space is estimated during
C     fac_ass.F / fac_ass_ELT.F and must be:
C     * Symmetric case: 1 for NASS + 1 for NBPANELS_L + NBPANELS_L + NASS
C     * Unsymmetric case: 1 + (1+NBPANELS_L+NASS) + (1+NBPANELS_U+NASS)
C     Size computation is in routine SMUMPS_OOC_GET_PP_SIZES.
C
C     At the end of the standard description of the structure of a node
C     (header, nb slaves, <slaves_list>, row indices, col indices), we
C     add, when panel version with pivoting is used:
C
C     NASS (nb of fully summed variables)
C     NBPANELS_L
C     PIVRPTR(1:NBPANELS_L)
C     PIV_L     (1:NASS)             NASS (=IW(IPOS)(or NASS-PIVRPTR(1) in
C     the future, after compression)
C     NBPANELS_U
C     PIVRPTR(1:NBPANELS_U)
C     PIV_U     (1:NASS)             NASS (=IW(IPOS)(or NASS-PIVRPTR(1) in
C     the future, after compression)
C
C
C     Output parameters:
C     =================
C     NBPANELS : nb of panels as estimated during assembly
C     I_PIVPTR : position in  IW of the starting of the pointer list
C     (of size NBPANELS) of the pointers to the list of pivots
C     I_PIV    : position in  IW of the starting of the pivot permutation list
C
      INTEGER, intent(out) :: NBPANELS, I_PIVPTR, I_PIV
      INTEGER, intent(in) :: TYPEF ! TYPEF_L or TYPEF_U
      INTEGER, intent(in) :: LIW, IPOS
      INTEGER IW(LIW)
C     Locals
      INTEGER I_NBPANELS, I_NASS
C
      I_NASS       = IPOS
      I_NBPANELS   = I_NASS + 1 ! L
      NBPANELS     = IW(I_NBPANELS) ! L
      I_PIVPTR     = I_NBPANELS + 1 ! L
      I_PIV        = I_PIVPTR + NBPANELS ! L
C     ... of size NASS = IW(I_NASS)
      IF (TYPEF==TYPEF_U) THEN
         I_NBPANELS   = I_PIV+IW(I_NASS) ! U
         NBPANELS     = IW(I_NBPANELS) ! U
         I_PIVPTR     = I_NBPANELS + 1 ! U
         I_PIV        = I_PIVPTR + NBPANELS ! U
      ENDIF
      RETURN
      END SUBROUTINE SMUMPS_GET_OOC_PERM_PTR
      SUBROUTINE SMUMPS_OOC_PP_SET_PTR(K50,NBPANELS_L,NBPANELS_U,
     &     NASS, IPOS, IW, LIW )
      IMPLICIT NONE
C
C     Purpose:
C     =======
C
C     Initialize the contents of PIV/PIVPTR/etc. that will store
C     pivoting information during the factorization.
C     NASS and NBPANELS are recorded. PIVPTR(1:NBPANELS)
C     is initialized to NASS+1. This will be modified during
C     the factorization in cases where permutations have to
C     be performed during the solve phase.
C
C     Arguments:
C     =========
C
      INTEGER K50
      INTEGER IPOS, NASS, NBPANELS_L, NBPANELS_U, LIW
      INTEGER IW(LIW)
C
C     Local variables:
C     ===============
C
      INTEGER IPOS_U
C     Executable statements
      IF (K50.EQ.1) THEN
         WRITE(*,*) "Internal error: SMUMPS_OOC_PP_SET_PTR called"
      ENDIF
      IW(IPOS)=NASS
      IW(IPOS+1)=NBPANELS_L
      IW(IPOS+2:IPOS+1+NBPANELS_L)=NASS+1
      IF (K50 == 0) THEN
         IPOS_U=IPOS+2+NASS+NBPANELS_L
         IW(IPOS_U)=NBPANELS_U
         IW(IPOS_U+1:IPOS_U+NBPANELS_U)=NASS+1
      ENDIF
      RETURN
      END SUBROUTINE SMUMPS_OOC_PP_SET_PTR
      SUBROUTINE SMUMPS_OOC_PP_TRYRELEASE_SPACE (
     &     IWPOS, IOLDPS, IW, LIW, MonBloc, NFRONT, KEEP
     &     )
      USE SMUMPS_OOC
      IMPLICIT NONE
      INCLUDE 'mumps_headers.h'
C
C     Purpose:
C     =======
C     If space used was at the top of the stack then
C     try to free space by detecting that
C     no permutation needs to be applied during
C     solve on panels.
C     One position is left (I_NASS) and set to -1
C     to indicate that permutation not needed at solve.
C
C     Arguments:
C     =========
C
      INTEGER, INTENT(IN)    :: IOLDPS, LIW, NFRONT,
     &     KEEP(500)
      INTEGER, INTENT(INOUT) :: IWPOS, IW(LIW)
      TYPE(IO_BLOCK), INTENT(IN):: MonBloc
C
C     Local variables:
C     ===============
C
      INTEGER :: NBPANELS_L,I_PIVRPTR_L, I_PIVR_L, NBPANELS_U,
     &     I_PIVRPTR_U, I_PIVR_U, XSIZE, IBEGOOC
      LOGICAL FREESPACE    ! set to true when permutation not needed
C     Executable statements
      IF (KEEP(50).EQ.1) RETURN ! no pivoting
C     --------------------------------
C     quick return if record is not at
C     the top of stack of L factors
      IF ((IOLDPS+IW(IOLDPS+XXI)).NE.IWPOS) RETURN
C     ---------------------------------------------
C     Panel+pivoting: get pointers on each subarray
C     ---------------------------------------------
      XSIZE   = KEEP(IXSZ)
      IBEGOOC = IOLDPS+2*NFRONT+6+IW(IOLDPS+5+XSIZE)+XSIZE
C     -- get L related data
      CALL SMUMPS_GET_OOC_PERM_PTR(TYPEF_L, NBPANELS_L,
     &     I_PIVRPTR_L, I_PIVR_L,
     &     IBEGOOC, IW, LIW)
      FREESPACE =
     &     (MonBloc%LastPiv.EQ.(IW(I_PIVRPTR_L)-1))
      IF (KEEP(50).EQ.0) THEN
C     -- get U related dataA
         CALL SMUMPS_GET_OOC_PERM_PTR(TYPEF_U, NBPANELS_U,
     &        I_PIVRPTR_U, I_PIVR_U,
     &        IBEGOOC, IW, LIW)
         FREESPACE =  FREESPACE .AND.
     &        (MonBloc%LastPiv.EQ.(IW(I_PIVRPTR_U)-1))
      ENDIF
C     ---------------------------------
C     Check if permutations eed be
C     performed on panels during solve
C     --------------------------------
      IF (FREESPACE) THEN
C     -- compress memory for that node: keep one entry set to -7777
         IW(IBEGOOC) = -7777    ! will be tested during solve
         IW(IOLDPS+XXI) = IBEGOOC
     &        - IOLDPS + 1      ! new size of inode's record
         IWPOS = IBEGOOC+1      ! move back to top of stack
      ENDIF
      RETURN
      END SUBROUTINE SMUMPS_OOC_PP_TRYRELEASE_SPACE
C
      SUBROUTINE SMUMPS_OOC_GET_PP_SIZES(K50, NBROW_L, NBCOL_U, NASS,
     &     NBPANELS_L, NBPANELS_U, LREQ)
      USE SMUMPS_OOC       ! To call SMUMPS_OOC_PANEL_SIZE
      IMPLICIT NONE
C
C     Purpose
C     =======
C
C     Compute the size of the workspace required to store the permutation
C     information during factorization, so that solve can permute back
C     what has to be permuted (this could not be done during factorization
C     because it was already on disk).
C
C     Arguments
C     =========
C
      INTEGER, intent(IN)  :: K50, NBROW_L, NBCOL_U, NASS
      INTEGER, intent(OUT) :: NBPANELS_L, NBPANELS_U, LREQ
      NBPANELS_L=-99999
      NBPANELS_U=-99999
C
C     Quick return in SPD case (no pivoting)
C
      IF (K50.EQ.1) THEN
         LREQ = 0
         RETURN
      ENDIF
C
C     L information is always computed
C
      NBPANELS_L = (NASS / SMUMPS_OOC_PANEL_SIZE(NBROW_L))+1
      LREQ =    1               ! Store NASS
     &     + 1                  ! Store NBPANELS_L
     &     + NASS               ! Store permutations
     &     + NBPANELS_L         ! Store pointers on permutations
      IF (K50.eq.0) THEN
C
C     Also take U information into account
C
         NBPANELS_U = (NASS / SMUMPS_OOC_PANEL_SIZE(NBCOL_U) ) +1
         LREQ = LREQ + 1        ! Store NBPANELS_U
     &        + NASS            ! Store permutations
     &        + NBPANELS_U      ! Store pointers on permutations
      ENDIF
      RETURN
      END SUBROUTINE SMUMPS_OOC_GET_PP_SIZES
      SUBROUTINE SMUMPS_OOC_PP_CHECK_PERM_FREED
     &           (IW_LOCATION, MUST_BE_PERMUTED)
      IMPLICIT NONE
      INTEGER, INTENT(IN) :: IW_LOCATION
      LOGICAL, INTENT(INOUT) :: MUST_BE_PERMUTED
C
C     Purpose
C     =======
C
C     Reset MUST_BE_PERMUTED to .FALSE. when we detect
C     that the SMUMPS_OOC_PP_TRY_RELEASE_SPACE has freed
C     the permutation information (see that routine).
C
      IF (IW_LOCATION .EQ. -7777) THEN
        MUST_BE_PERMUTED = .FALSE.
      ENDIF
      RETURN
      END SUBROUTINE SMUMPS_OOC_PP_CHECK_PERM_FREED