OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
cfac_front_aux.F
Go to the documentation of this file.
1C
2C This file is part of MUMPS 5.5.1, released
3C on Tue Jul 12 13:17:24 UTC 2022
4C
5C
6C Copyright 1991-2022 CERFACS, CNRS, ENS Lyon, INP Toulouse, Inria,
7C Mumps Technologies, University of Bordeaux.
8C
9C This version of MUMPS is provided to you free of charge. It is
10C released under the CeCILL-C license
11C (see doc/CeCILL-C_V1-en.txt, doc/CeCILL-C_V1-fr.txt, and
12C https://cecill.info/licences/Licence_CeCILL-C_V1-en.html)
13C
14 MODULE cmumps_fac_front_aux_m
15 CONTAINS
16 SUBROUTINE cmumps_fac_h(NFRONT,NASS,IW,LIW,A,LA,
17 & INOPV,NOFFW,
18 & DET_EXPW, DET_MANTW, DET_SIGNW,
19 & IOLDPS,POSELT,UU,SEUIL,KEEP, KEEP8, DKEEP,
20 & PP_FIRST2SWAP_L, PP_LastPanelonDisk_L,
21 & PP_LastPIVRPTRFilled_L,
22 & PP_FIRST2SWAP_U, PP_LastPanelonDisk_U,
23 & PP_LastPIVRPTRFilled_U,MAXFROMN,IS_MAXFROMN_AVAIL,
24 & Inextpiv, OOC_EFFECTIVE_ON_FRONT, NVSCHUR
25 &)
26!$ USE omp_lib
27 USE mumps_ooc_common
28 IMPLICIT NONE
29 INTEGER NFRONT,NASS,LIW,INOPV
30 INTEGER(8) :: LA
31 INTEGER :: KEEP(500)
32 INTEGER(8) :: KEEP8(150)
33 REAL :: DKEEP(230)
34 REAL UU, SEUIL
35 COMPLEX A(LA)
36 INTEGER IW(LIW)
37 REAL, intent(in) :: MAXFROMN
38 LOGICAL, intent(inout) :: IS_MAXFROMN_AVAIL
39 INTEGER, intent(inout) :: Inextpiv
40 LOGICAL, intent(in) :: OOC_EFFECTIVE_ON_FRONT
41 INTEGER, intent(in) :: NVSCHUR
42 REAL AMROW
43 REAL RMAX, SEUIL_LOC
44 COMPLEX SWOP
45 INTEGER(8) :: APOS, POSELT
46 INTEGER(8) :: J1, J2, J3_8, JJ, IDIAG
47 INTEGER(8) :: J1_ini
48 INTEGER(8) :: NFRONT8
49 INTEGER IOLDPS
50 INTEGER NPIV,IPIV,IPIV_SHIFT
51 INTEGER, intent(inout) :: NOFFW
52 INTEGER, intent(inout) :: DET_EXPW, DET_SIGNW
53 COMPLEX, intent(inout) :: DET_MANTW
54 INTEGER J, J3
55 INTEGER NPIVP1,JMAX,ISW,ISWPS1
56 INTEGER ISWPS2,KSW,XSIZE
57 INTEGER I_PIVRPTR_L, I_PIVR_L, NBPANELS_L
58 INTEGER I_PIVRPTR_U, I_PIVR_U, NBPANELS_U
59 INTEGER PP_FIRST2SWAP_L, PP_LastPanelonDisk_L,
60 & pp_lastpivrptrfilled_l,
61 & pp_first2swap_u, pp_lastpanelondisk_u,
62 & pp_lastpivrptrfilled_u
63 INTEGER ISHIFT, K206
64 INTEGER CMUMPS_IXAMAX
65 include 'mumps_headers.h'
66 INTRINSIC max
67 REAL, PARAMETER :: RZERO = 0.0e0
68#if defined(_OPENMP)
69 INTEGER :: NOMP, CHUNK, K360
70 K360 = keep(360)
71 nomp = omp_get_max_threads()
72#endif
73 seuil_loc = max(dkeep(1), seuil)
74 nfront8 = int(nfront,8)
75 inopv = 0
76 xsize = keep(ixsz)
77 npiv = iw(ioldps+1+xsize)
78 npivp1 = npiv + 1
79 k206 = keep(206)
80 IF ((keep(50).NE.1).AND.ooc_effective_on_front) THEN
81 CALL cmumps_get_ooc_perm_ptr(typef_l, nbpanels_l,
82 & i_pivrptr_l, i_pivr_l,
83 & ioldps+2*nfront+6+iw(ioldps+5+xsize)
84 & +keep(ixsz),
85 & iw, liw)
86 CALL cmumps_get_ooc_perm_ptr(typef_u, nbpanels_u,
87 & i_pivrptr_u, i_pivr_u,
88 & ioldps+2*nfront+6+iw(ioldps+5+xsize)+xsize,
89 & iw, liw)
90 ENDIF
91 ishift = 0
92 IF (k206.GE.1) THEN
93 IF (inextpiv.GT.npivp1.AND.inextpiv.LE.nass) THEN
94 ishift = inextpiv - npivp1
95 ENDIF
96 IF (ishift.GT.0.AND.is_maxfromn_avail) THEN
97 ipiv = npivp1
98 apos = poselt + nfront8*int(npiv,8) + int(ipiv-1,8)
99 idiag = apos + int(ipiv - npivp1,8)*nfront8
100 IF (abs(a(idiag)) .GE. uu*maxfromn .AND.
101 & abs(a(idiag)) .GT. max(seuil_loc,tiny(rmax))
102 & ) THEN
103 ishift = 0
104 ENDIF
105 ENDIF
106 IF ( ishift .GT. 0) THEN
107 is_maxfromn_avail = .false.
108 ENDIF
109 ENDIF
110 DO 460 ipiv_shift=npivp1+ishift,nass+ishift
111 IF (ipiv_shift .LE. nass) THEN
112 ipiv=ipiv_shift
113 ELSE
114 ipiv=ipiv_shift-nass-1+npivp1
115 ENDIF
116 apos = poselt + nfront8*int(npiv,8) + int(ipiv-1,8)
117 jmax = 1
118 amrow = rzero
119 j1 = apos
120 j3 = nass -npiv
121 jmax = cmumps_ixamax(j3,a(j1),nfront,keep(360))
122 jj = j1 + int(jmax-1,8)*nfront8
123 amrow = abs(a(jj))
124 rmax = amrow
125 j1 = apos + int(nass-npiv,8) * nfront8
126 j3 = nfront - nass - keep(253)-nvschur
127 IF (is_maxfromn_avail) THEN
128 rmax = max(maxfromn,rmax)
129 is_maxfromn_avail = .false.
130 ELSE
131 IF (j3.EQ.0) GOTO 370
132 IF (keep(351).EQ.1) THEN
133 j1_ini = j1
134!$ CHUNK = max(K360/2,(J3+NOMP-1)/NOMP)
135!$OMP PARALLEL DO schedule(static, CHUNK)
136!$OMP& FIRSTPRIVATE(J1_ini,NFRONT8,J3)
137!$OMP& REDUCTION(max:RMAX) IF (J3.GE.K360)
138 DO j=1,j3
139 rmax = max(abs(a(j1_ini + int(j-1,8) * nfront8)),
140 & rmax)
141 END DO
142!$OMP END PARALLEL DO
143 ELSE
144 DO j=1,j3
145 rmax = max(abs(a(j1)), rmax)
146 j1 = j1 + nfront8
147 END DO
148 ENDIF
149 END IF
150 370 IF (rmax.LE.tiny(rmax)) GO TO 460
151 idiag = apos + int(ipiv - npivp1,8)*nfront8
152 IF (abs(a(idiag)) .GE. uu*rmax .AND.
153 & abs(a(idiag)) .GT. max(seuil_loc,tiny(rmax)) ) THEN
154 jmax = ipiv - npiv
155 GO TO 380
156 ENDIF
157 IF ( .NOT. ( amrow .GE. uu*rmax .AND.
158 & amrow .GT. max(seuil_loc,tiny(rmax))
159 & )
160 & ) GO TO 460
161 noffw = noffw + 1
162 380 CONTINUE
163 IF (k206.GE.1) THEN
164 inextpiv = ipiv + 1
165 ENDIF
166 CALL cmumps_update_minmax_pivot
167 & ( abs(a(apos + int(jmax - 1,8) * nfront8 )),
168 & dkeep, keep, .false.)
169 IF (keep(258) .NE. 0) THEN
170 CALL cmumps_updatedeter(
171 & a(apos + int(jmax - 1,8) * nfront8 ),
172 & det_mantw, det_expw )
173 ENDIF
174 IF (ipiv.EQ.npivp1) GO TO 400
175 IF (keep(405) .EQ.0) THEN
176 keep8(80) = keep8(80)+1
177 ELSE
178!$OMP ATOMIC UPDATE
179 keep8(80) = keep8(80)+1
180!$OMP END ATOMIC
181 ENDIF
182 det_signw = - det_signw
183 j1 = poselt + int(npiv,8)
184 j3_8 = poselt + int(ipiv-1,8)
185 DO j= 1,nfront
186 swop = a(j1)
187 a(j1) = a(j3_8)
188 a(j3_8) = swop
189 j1 = j1 + nfront8
190 j3_8 = j3_8 + nfront8
191 END DO
192 iswps1 = ioldps + 5 + npivp1 + nfront + xsize
193 iswps2 = ioldps + 5 + ipiv + nfront + xsize
194 isw = iw(iswps1)
195 iw(iswps1) = iw(iswps2)
196 iw(iswps2) = isw
197 400 IF (jmax.EQ.1) GO TO 420
198 det_signw = -det_signw
199 j1 = poselt + int(npiv,8) * nfront8
200 j2 = poselt + int(npiv + jmax - 1,8) * nfront8
201 DO ksw=1,nfront
202 swop = a(j1)
203 a(j1) = a(j2)
204 a(j2) = swop
205 j1 = j1 + 1_8
206 j2 = j2 + 1_8
207 END DO
208 iswps1 = ioldps + 5 + npiv + 1 + xsize
209 iswps2 = ioldps + 5 + npiv + jmax + xsize
210 isw = iw(iswps1)
211 iw(iswps1) = iw(iswps2)
212 iw(iswps2) = isw
213 GO TO 420
214 460 CONTINUE
215 inopv = 1
216 GOTO 430
217 420 CONTINUE
218 IF (ooc_effective_on_front) THEN
219 IF (keep(251).EQ.0) THEN
220 CALL cmumps_store_perminfo( iw(i_pivrptr_l),
221 & nbpanels_l,
222 & iw(i_pivr_l), nass, npivp1, npiv+jmax,
223 & pp_lastpanelondisk_l,
224 & pp_lastpivrptrfilled_l)
225 ENDIF
226 CALL cmumps_store_perminfo( iw(i_pivrptr_u),
227 & nbpanels_u,
228 & iw(i_pivr_u), nass, npivp1, ipiv,
229 & pp_lastpanelondisk_u,
230 & pp_lastpivrptrfilled_u)
231 ENDIF
232 430 CONTINUE
233 is_maxfromn_avail = .false.
234 RETURN
235 END SUBROUTINE cmumps_fac_h
236 SUBROUTINE cmumps_fac_n(NFRONT,NASS,IW,LIW,A,LA,
237 & IOLDPS,POSELT,IFINB,XSIZE,
238 & KEEP,MAXFROMN,IS_MAXFROMN_AVAIL,NVSCHUR)
239!$ USE OMP_LIB
240 IMPLICIT NONE
241 include 'mumps_headers.h'
242 INTEGER NFRONT,NASS,LIW,IFINB
243 INTEGER(8) :: LA
244 COMPLEX A(LA)
245 INTEGER IW(LIW)
246 COMPLEX ALPHA,VALPIV
247 INTEGER(8) :: APOS, POSELT, UUPOS, LPOS, IRWPOS
248 INTEGER(8) :: NFRONT8
249 INTEGER IOLDPS,NPIV,XSIZE
250 INTEGER, intent(in) :: KEEP(500)
251 REAL, intent(inout) :: MAXFROMN
252 LOGICAL, intent(inout) :: IS_MAXFROMN_AVAIL
253 INTEGER, intent(in) :: NVSCHUR
254 INTEGER NEL,IROW,NEL2,JCOL,NELMAXM
255 INTEGER NPIVP1
256 COMPLEX, PARAMETER :: ONE=(1.0e0,0.0e0)
257#if defined(_OPENMP)
258 LOGICAL:: OMP_FLAG
259 INTEGER:: NOMP, K360, CHUNK
260 nomp = omp_get_max_threads()
261 k360 = keep(360)
262#endif
263 nfront8=int(nfront,8)
264 npiv = iw(ioldps+1+xsize)
265 npivp1 = npiv + 1
266 nel = nfront - npivp1
267 nelmaxm= nel -keep(253)-nvschur
268 nel2 = nass - npivp1
269 ifinb = 0
270 IF (npivp1.EQ.nass) ifinb = 1
271 apos = poselt + int(npiv,8)*(nfront8 + 1_8)
272 valpiv = one/a(apos)
273#if defined(_OPENMP)
274 omp_flag = .false.
275 chunk = max(nel,1)
276 IF (nomp.GT.1) THEN
277 IF (nel.LT.k360) THEN
278 IF (nel*nel2.GE.keep(361)) THEN
279 omp_flag = .true.
280 chunk = max(20, (nel+nomp-1)/nomp)
281 ENDIF
282 ELSE
283 omp_flag = .true.
284 chunk = max(k360/2, (nel+nomp-1)/nomp)
285 ENDIF
286 ENDIF
287#endif
288 IF (keep(351).EQ.2) THEN
289 maxfromn = 0.0e0
290 IF (nel2 > 0) THEN
291 is_maxfromn_avail = .true.
292 ENDIF
293!$OMP PARALLEL DO schedule(static, CHUNK)
294!$OMP& PRIVATE(LPOS, UUPOS, IRWPOS, ALPHA, JCOL)
295!$OMP& FIRSTPRIVATE(APOS,NFRONT8,VALPIV,NEL,NEL2)
296!$OMP& REDUCTION(max:MAXFROMN) IF(OMP_FLAG)
297 DO irow = 1, nel
298 lpos = apos + nfront8*int(irow,8)
299 a(lpos) = a(lpos)*valpiv
300 alpha = -a(lpos)
301 irwpos = lpos + 1_8
302 uupos = apos + 1_8
303 IF (nel2 > 0) THEN
304 a(irwpos) = a(irwpos) + alpha*a(uupos)
305 IF (irow.LE.nelmaxm)
306 & maxfromn=max(maxfromn, abs(a(irwpos)))
307 irwpos = irwpos+1_8
308 uupos = uupos+1_8
309 DO jcol = 2, nel2
310 a(irwpos) = a(irwpos) + alpha*a(uupos)
311 irwpos = irwpos+1_8
312 uupos = uupos+1_8
313 ENDDO
314 ENDIF
315 END DO
316!$OMP END PARALLEL DO
317 ELSE
318!$OMP PARALLEL DO schedule(static, CHUNK)
319!$OMP& FIRSTPRIVATE(APOS,NFRONT8,VALPIV,NEL,NEL2)
320!$OMP& PRIVATE(LPOS, UUPOS, IRWPOS, ALPHA, JCOL) IF(OMP_FLAG)
321 DO irow = 1, nel
322 lpos = apos + nfront8*int(irow,8)
323 a(lpos) = a(lpos)*valpiv
324 alpha = -a(lpos)
325 irwpos = lpos + 1_8
326 uupos = apos + 1_8
327 DO jcol = 1, nel2
328 a(irwpos) = a(irwpos) + alpha*a(uupos)
329 irwpos = irwpos+1_8
330 uupos = uupos+1_8
331 ENDDO
332 ENDDO
333!$OMP END PARALLEL DO
334 ENDIF
335 RETURN
336 END SUBROUTINE cmumps_fac_n
337 SUBROUTINE cmumps_fac_pt_setlock427( K427_OUT, K427,
338 & K405, K222, NEL1, NASS )
339 INTEGER, INTENT(IN) :: K427, K405, K222, NEL1, NASS
340 INTEGER, INTENT(OUT) :: K427_OUT
341 k427_out = k427
342 IF ( k427_out .GT. 0 ) k427_out = 0
343 IF ( k427_out .LT. 0 ) k427_out = -1
344 IF ( k427_out .GT. 99 ) k427_out = 0
345 IF ( k427_out .LT. -100 ) k427_out = -1
346 RETURN
347 END SUBROUTINE cmumps_fac_pt_setlock427
348 SUBROUTINE cmumps_fac_p(A,LA,NFRONT,
349 & NPIV,NASS,POSELT,CALL_UTRSM, KEEP, INODE,
350 & CALL_OOC, IWFAC, LIWFAC, LAFAC, MonBloc, MYID, KEEP8,
351 & LNextPiv2beWritten, UNextPiv2beWritten,
352 & IFLAG )
353 USE cmumps_ooc, ONLY : io_block, typef_both_lu,
354 & cmumps_ooc_io_lu_panel
355 USE mumps_ooc_common, ONLY : strat_try_write
356 IMPLICIT NONE
357 INTEGER(8) :: LA,POSELT,LAFAC
358 COMPLEX A(LA)
359 INTEGER NFRONT, NPIV, NASS
360 LOGICAL, INTENT(IN) :: CALL_UTRSM
361 INTEGER, INTENT(INOUT) :: IFLAG
362 LOGICAL, INTENT(IN) :: CALL_OOC
363 INTEGER LIWFAC, MYID,
364 & lnextpiv2bewritten, unextpiv2bewritten
365 INTEGER IWFAC(LIWFAC)
366 TYPE(IO_BLOCK) :: MonBloc
367 INTEGER :: KEEP(500)
368 INTEGER(8) :: KEEP8(150)
369 INTEGER(8) :: LPOS, LPOS1, LPOS2, UPOS
370 INTEGER NEL1, NEL11, IFLAG_OOC
371 INTEGER :: INODE
372 COMPLEX ALPHA, ONE
373 parameter(one=(1.0e0,0.0e0), alpha=(-1.0e0,0.0e0))
374 include 'mumps_headers.h'
375 nel1 = nfront - nass
376 nel11 = nfront - npiv
377 lpos2 = poselt + int(nass,8)*int(nfront,8)
378 lpos = lpos2 + int(npiv,8)
379 lpos1 = poselt + int(npiv,8)
380 upos = poselt + int(nass,8)
381 IF ( call_utrsm ) THEN
382 CALL ctrsm('R', 'U', 'N', 'U', nel1, npiv, one,
383 & a(poselt), nfront, a(upos), nfront)
384 ENDIF
385 CALL ctrsm('L','L','N','N',npiv,nel1,one,a(poselt),nfront,
386 & a(lpos2),nfront)
387 IF (call_ooc) THEN
388 CALL cmumps_ooc_io_lu_panel
389 & ( strat_try_write, typef_both_lu,
390 & a(poselt), lafac, monbloc,
391 & lnextpiv2bewritten, unextpiv2bewritten,
392 & iwfac, liwfac,
393 & myid, keep8(31), iflag_ooc,
394 & .false. )
395 IF (iflag_ooc .LT. 0) THEN
396 iflag = iflag_ooc
397 GOTO 500
398 ENDIF
399 ENDIF
400 CALL cgemm('N','N',nel11,nel1,npiv,alpha,a(lpos1),
401 & nfront,a(lpos2),nfront,one,a(lpos),nfront)
402 IF ((call_utrsm).AND.(nass-npiv.GT.0)) THEN
403 lpos2 = poselt + int(npiv,8)*int(nfront,8)
404 lpos = lpos2 + int(nass,8)
405 CALL cgemm('N','N',nel1,nass-npiv,npiv,alpha,a(upos),
406 & nfront,a(lpos2),nfront,one,a(lpos),nfront)
407 ENDIF
408 500 CONTINUE
409 RETURN
410 END SUBROUTINE cmumps_fac_p
411 SUBROUTINE cmumps_fac_t(A,LA,NPIVB,NFRONT,
412 & NPIV,NASS,POSELT)
413 IMPLICIT NONE
414 INTEGER NPIVB,NASS
415 INTEGER(8) :: LA
416 COMPLEX A(LA)
417 INTEGER(8) :: APOS, POSELT
418 INTEGER NFRONT, NPIV, NASSL
419 INTEGER(8) :: LPOS, LPOS1, LPOS2
420 INTEGER NEL1, NEL11, NPIVE
421 COMPLEX ALPHA, ONE
422 parameter(one=(1.0e0,0.0e0), alpha=(-1.0e0,0.0e0))
423 nel1 = nfront - nass
424 nel11 = nfront - npiv
425 npive = npiv - npivb
426 nassl = nass - npivb
427 apos = poselt + int(npivb,8)*int(nfront,8)
428 & + int(npivb,8)
429 lpos2 = apos + int(nassl,8)
430 CALL ctrsm('r','u','n','u',NEL1,NPIVE,ONE,A(APOS),NFRONT,
431 & A(LPOS2),NFRONT)
432 LPOS = LPOS2 + int(NFRONT,8)*int(NPIVE,8)
433 LPOS1 = APOS + int(NFRONT,8)*int(NPIVE,8)
434 CALL cgemm('n','n',NEL1,NEL11,NPIVE,ALPHA,A(LPOS2),
435 & NFRONT,A(LPOS1),NFRONT,ONE,A(LPOS),NFRONT)
436 RETURN
437 END SUBROUTINE CMUMPS_FAC_T
438 SUBROUTINE CMUMPS_FAC_SQ(IBEG_BLOCK, IEND_BLOCK, NPIV,
439 & NFRONT, LAST_ROW, LAST_COL, A, LA, POSELT,
440 & FIRST_COL, CALL_LTRSM, CALL_UTRSM, CALL_GEMM,
441 & WITH_COMM_THREAD, LR_ACTIVATED
442 & )
443!$ USE OMP_LIB
444#if defined(_OPENMP)
445 USE CMUMPS_BUF
446#endif
447 IMPLICIT NONE
448 INTEGER, intent(in) :: IBEG_BLOCK, IEND_BLOCK
449 INTEGER, intent(in) :: NPIV, NFRONT, LAST_ROW, LAST_COL
450 INTEGER, intent(in) :: FIRST_COL
451 INTEGER(8), intent(in) :: LA
452 COMPLEX, intent(inout) :: A(LA)
453 INTEGER(8), intent(in) :: POSELT
454 LOGICAL, intent(in) :: CALL_LTRSM, CALL_UTRSM, CALL_GEMM
455 LOGICAL, intent(in) :: WITH_COMM_THREAD, LR_ACTIVATED
456 INTEGER(8) :: NFRONT8, LPOSN, LPOS2N
457 INTEGER(8) :: LPOS, LPOS1, LPOS2, UPOS, POSELT_LOCAL
458 INTEGER :: NELIM, LKJIW, NEL1, NEL11, UTRSM_NCOLS
459 COMPLEX ALPHA, ONE
460 PARAMETER (ONE=(1.0E0,0.0E0), ALPHA=(-1.0E0,0.0E0))
461!$ INTEGER :: NOMP
462!$ LOGICAL :: TRSM_GEMM_FINISHED
463!$ LOGICAL :: SAVE_NESTED, SAVE_DYNAMIC
464 NFRONT8= int(NFRONT,8)
465 NELIM = IEND_BLOCK - NPIV
466 NEL1 = LAST_ROW - IEND_BLOCK
467 IF ( NEL1 < 0 ) THEN
468 WRITE(*,*)
469 & "Internal error 1 in CMUMPS_FAC_SQ,IEND_BLOCK>LAST_ROW",
470 & IEND_BLOCK, LAST_ROW
471 CALL MUMPS_ABORT()
472 ENDIF
473 LKJIW = NPIV - IBEG_BLOCK + 1
474 NEL11 = LAST_COL - NPIV
475 LPOS2 = POSELT + int(IEND_BLOCK,8)*NFRONT8 + int(IBEG_BLOCK-1,8)
476 UTRSM_NCOLS = LAST_COL - FIRST_COL
477 UPOS = POSELT + int(IBEG_BLOCK-1,8)*NFRONT8 + int(FIRST_COL,8)
478 POSELT_LOCAL = POSELT + int(IBEG_BLOCK-1,8)*NFRONT8
479 & + int(IBEG_BLOCK-1,8)
480.NE..AND..NE. IF ((NEL10)(LKJIW0)) THEN
481.EQV. IF (WITH_COMM_THREAD .FALSE.) THEN
482 IF (CALL_LTRSM) THEN
483 CALL ctrsm('l','l','n','n',LKJIW,NEL1,ONE,
484 & A(POSELT_LOCAL),NFRONT,
485 & A(LPOS2),NFRONT)
486 ENDIF
487 IF (CALL_UTRSM) THEN
488 CALL ctrsm('r','u','n','u',UTRSM_NCOLS,LKJIW,ONE,
489 & A(POSELT_LOCAL),NFRONT,
490 & A(UPOS),NFRONT)
491 LPOS2N = POSELT + int(NPIV,8)*NFRONT8 + int(IBEG_BLOCK-1,8)
492 LPOSN = POSELT + int(NPIV,8)*NFRONT8 + int(FIRST_COL,8)
493 CALL cgemm('n','n',UTRSM_NCOLS,NELIM,
494 & LKJIW,ALPHA,A(UPOS),NFRONT,A(LPOS2N),
495 & NFRONT,ONE,A(LPOSN),NFRONT)
496 ENDIF
497 IF (CALL_GEMM) THEN
498 LPOS = LPOS2 + int(LKJIW,8)
499 LPOS1 = POSELT_LOCAL + int(LKJIW,8)
500 CALL cgemm('n','n',NEL11,NEL1,LKJIW,ALPHA,A(LPOS1),
501 & NFRONT,A(LPOS2),NFRONT,ONE,A(LPOS),NFRONT)
502 ENDIF
503 ELSE
504!$ NOMP = OMP_GET_MAX_THREADS()
505!$ CALL OMP_SET_NUM_THREADS(2)
506!$ SAVE_NESTED = OMP_GET_NESTED()
507!$ SAVE_DYNAMIC = OMP_GET_DYNAMIC()
508!$ CALL OMP_SET_NESTED(.TRUE.)
509!$ CALL OMP_SET_DYNAMIC(.FALSE.)
510!$ TRSM_GEMM_FINISHED = .FALSE.
511!$OMP PARALLEL SHARED(TRSM_GEMM_FINISHED)
512.EQ.!$ IF (OMP_GET_THREAD_NUM() 1) THEN
513#if defined(WORKAROUNDINTELILP64OPENMPLIMITATION)
514!$ CALL OMP_SET_NUM_THREADS(int(NOMP,4))
515#else
516!$ CALL OMP_SET_NUM_THREADS(NOMP)
517#endif
518 IF (CALL_LTRSM) THEN
519 CALL ctrsm('l','l','n','n',LKJIW,NEL1,ONE,
520 & A(POSELT_LOCAL),NFRONT,
521 & A(LPOS2),NFRONT)
522 ENDIF
523 IF (CALL_UTRSM) THEN
524 CALL ctrsm('r','u','n','u',UTRSM_NCOLS,LKJIW,ONE,
525 & A(POSELT_LOCAL),NFRONT,
526 & A(UPOS),NFRONT)
527 LPOS2N = POSELT + int(NPIV,8)*NFRONT8 + int(IBEG_BLOCK-1,8)
528 LPOSN = POSELT + int(NPIV,8)*NFRONT8 + int(FIRST_COL,8)
529 CALL cgemm('n','n',UTRSM_NCOLS,NELIM,
530 & LKJIW,ALPHA,A(UPOS),NFRONT,A(LPOS2N),
531 & NFRONT,ONE,A(LPOSN),NFRONT)
532 ENDIF
533 IF (CALL_GEMM) THEN
534 LPOS = LPOS2 + int(LKJIW,8)
535 LPOS1 = POSELT_LOCAL + int(LKJIW,8)
536 CALL cgemm('n','n',NEL11,NEL1,LKJIW,ALPHA,A(LPOS1),
537 & NFRONT,A(LPOS2),NFRONT,ONE,A(LPOS),NFRONT)
538 END IF
539!$ TRSM_GEMM_FINISHED = .TRUE.
540!$ ELSE
541.NOT.!$ DO WHILE ( TRSM_GEMM_FINISHED)
542!$ CALL CMUMPS_BUF_TEST()
543!$ CALL MUMPS_USLEEP(10000)
544!$ END DO
545!$ END IF
546!$OMP END PARALLEL
547!$ CALL OMP_SET_NESTED(SAVE_NESTED)
548!$ CALL OMP_SET_DYNAMIC(SAVE_DYNAMIC)
549#if defined(WORKAROUNDINTELILP64OPENMPLIMITATION)
550!$ CALL OMP_SET_NUM_THREADS(int(NOMP,4))
551#else
552!$ CALL OMP_SET_NUM_THREADS(NOMP)
553#endif
554 ENDIF
555 ELSE
556.AND..NE. IF (CALL_UTRSMUTRSM_NCOLS0) THEN
557 CALL ctrsm('r','u','n','u',utrsm_ncols,lkjiw,one,
558 & a(poselt_local),nfront,
559 & a(upos),nfront)
560 lpos2n = poselt + int(npiv,8)*nfront8 + int(ibeg_block-1,8)
561 lposn = poselt + int(npiv,8)*nfront8 + int(first_col,8)
562 CALL cgemm('N','N',utrsm_ncols,nelim,
563 & lkjiw,alpha,a(upos),nfront,a(lpos2n),
564 & nfront,one,a(lposn),nfront)
565 ENDIF
566 ENDIF
567 RETURN
568 END SUBROUTINE cmumps_fac_sq
569 SUBROUTINE cmumps_fac_mq(IBEG_BLOCK,IEND_BLOCK,
570 & NFRONT, NASS, NPIV, LAST_COL, A, LA, POSELT, IFINB,
571 & LR_ACTIVATED
572 & )
573 IMPLICIT NONE
574 INTEGER, intent(in) :: IBEG_BLOCK, IEND_BLOCK, NFRONT,
575 & nass, npiv, last_col
576 INTEGER, intent(out) :: IFINB
577 INTEGER(8), intent(in) :: LA, POSELT
578 COMPLEX, intent(inout) :: A(LA)
579 LOGICAL, intent(in) :: LR_ACTIVATED
580 COMPLEX :: VALPIV
581 INTEGER(8) :: APOS, UUPOS, LPOS
582 INTEGER(8) :: NFRONT8
583 COMPLEX :: ONE, ALPHA
584 INTEGER :: NEL2,NPIVP1,KROW,NEL
585 parameter(one=(1.0e0,0.0e0), alpha=(-1.0e0,0.0e0))
586 nfront8= int(nfront,8)
587 npivp1 = npiv + 1
588 nel = last_col - npivp1
589 ifinb = 0
590 nel2 = iend_block - npivp1
591 IF (nel2.EQ.0) THEN
592 IF (iend_block.EQ.nass) THEN
593 ifinb = -1
594 ELSE
595 ifinb = 1
596 ENDIF
597 ELSE
598 apos = poselt + int(npiv,8)*(nfront8 + 1_8)
599 valpiv = one/a(apos)
600 lpos = apos + nfront8
601 DO 541 krow = 1,nel2
602 a(lpos) = a(lpos)*valpiv
603 lpos = lpos + nfront8
604 541 CONTINUE
605 lpos = apos + nfront8
606 uupos = apos + 1_8
607#if defined(MUMPS_USE_BLAS2)
608 CALL cgeru(nel,nel2,alpha,a(uupos),1,a(lpos),nfront,
609 & a(lpos+1_8),nfront)
610#else
611 CALL cgemm('N','N',nel,nel2,1,alpha,a(uupos),nel,
612 & a(lpos),nfront,one,a(lpos+1_8),nfront)
613#endif
614 ENDIF
615 RETURN
616 END SUBROUTINE cmumps_fac_mq
617 SUBROUTINE cmumps_fac_fr_update_cbrows( INODE, NFRONT, NASS,
618 & CALL_UTRSM, A, LA, LAFAC, POSELT, IW, LIW, IOLDPS,
619 & MonBloc, MYID, NOFFW,
620 & DET_EXPW, DET_MANTW, DET_SIGNW,
621 & LIWFAC,
622 & PP_FIRST2SWAP_L, PP_FIRST2SWAP_U,
623 & LNextPiv2beWritten, UNextPiv2beWritten,
624 & PP_LastPIVRPTRFilled_L, PP_LastPIVRPTRFilled_U,
625 &
626 & XSIZE, SEUIL, UU, DKEEP, KEEP8, KEEP, IFLAG,
627 & OOC_EFFECTIVE_ON_FRONT, NVSCHUR)
628 USE cmumps_ooc, ONLY: io_block
629 IMPLICIT NONE
630 INTEGER, intent(in) :: INODE, NFRONT, NASS,
631 & liw, myid, xsize, ioldps, liwfac
632 INTEGER(8), intent(in) :: LA, POSELT
633 INTEGER, intent(inout) :: NOFFW
634 INTEGER, intent(inout) :: DET_EXPW, DET_SIGNW
635 COMPLEX, intent(inout) :: DET_MANTW
636 INTEGER, intent(inout) :: PP_FIRST2SWAP_L, PP_FIRST2SWAP_U,
637 & LNextPiv2beWritten, UNextPiv2beWritten,
638 & pp_lastpivrptrfilled_l, pp_lastpivrptrfilled_u,
639 & iflag
640 LOGICAL, intent(in) :: CALL_UTRSM
641 INTEGER, intent(inout) :: IW(LIW)
642 COMPLEX, intent(inout) :: A(LA)
643 REAL, intent(in) :: SEUIL, UU, DKEEP(230)
644 INTEGER, intent(in) :: KEEP( 500 )
645 INTEGER(8), intent(inout) :: LAFAC
646 INTEGER(8) :: KEEP8(150)
647 INTEGER, intent(in) :: NVSCHUR
648 TYPE(IO_BLOCK), intent(inout) :: MonBloc
649 LOGICAL, intent(in) :: OOC_EFFECTIVE_ON_FRONT
650 INTEGER :: NPIV, NEL1, IBEG_BLOCK, IFINB, INOPV
651 INTEGER Inextpiv
652 REAL :: MAXFROMN
653 LOGICAL :: IS_MAXFROMN_AVAIL
654 NPIV = iw(ioldps+1+xsize)
655 nel1 = nfront - nass
656 IF (keep(206).GE.1) THEN
657 inextpiv = 1
658 ELSE
659 inextpiv = 0
660 ENDIF
661 IF ((npiv.GT.0).AND.(nel1.GT.0)) THEN
662 IF (ooc_effective_on_front) THEN
663 monbloc%LastPiv = npiv
664 ENDIF
665 CALL cmumps_fac_p(a,la,nfront, npiv, nass, poselt,
666 & call_utrsm, keep, inode,
667 & ooc_effective_on_front, iw(ioldps),
668 & liwfac, lafac,
669 & monbloc, myid, keep8,
670 & lnextpiv2bewritten, unextpiv2bewritten,
671 & iflag)
672 ENDIF
673 npiv = iw(ioldps+1+xsize)
674 ibeg_block = npiv
675 IF (nass.EQ.npiv) GOTO 500
676 is_maxfromn_avail = .false.
677 120 CALL cmumps_fac_h(nfront,nass,iw,liw,a,la,
678 & inopv, noffw,
679 & det_expw, det_mantw, det_signw,
680 & ioldps,poselt,uu,seuil,
681 & keep, keep8, dkeep,
682 & pp_first2swap_l, monbloc%LastPanelWritten_L,
683 & pp_lastpivrptrfilled_l,
684 & pp_first2swap_u, monbloc%LastPanelWritten_U,
685 & pp_lastpivrptrfilled_u, maxfromn, is_maxfromn_avail,
686 & inextpiv, ooc_effective_on_front, nvschur
687 & )
688 IF (inopv.NE.1) THEN
689 CALL cmumps_fac_n(nfront,nass,iw,liw,a,la,
690 & ioldps,poselt,ifinb,xsize,
691 & keep, maxfromn, is_maxfromn_avail,
692 & nvschur)
693 iw(ioldps+1+xsize) = iw(ioldps+1+xsize) + 1
694 IF (ifinb.EQ.0) GOTO 120
695 ENDIF
696 npiv = iw(ioldps+1+xsize)
697 nel1 = nfront - nass
698 IF ((npiv.LE.ibeg_block).OR.(nel1.EQ.0)) GO TO 500
699 CALL cmumps_fac_t(a,la,ibeg_block,
700 & nfront,npiv,nass,poselt)
701 500 CONTINUE
702 RETURN
703 END SUBROUTINE cmumps_fac_fr_update_cbrows
704 SUBROUTINE cmumps_fac_i(NFRONT,NASS,LAST_ROW,
705 & IBEG_BLOCK, IEND_BLOCK,
706 & N,INODE,IW,LIW,A,LA,
707 & INOPV,NOFFW,NBTINYW,
708 & DET_EXPW, DET_MANTW, DET_SIGNW,
709 & IFLAG,IOLDPS,POSELT,UU,SEUIL,KEEP,KEEP8,
710 & DKEEP,PIVNUL_LIST,LPN_LIST,
711 &
712 & PP_FIRST2SWAP_L, PP_LastPanelonDisk_L,
713 & PP_LastPIVRPTRFilled_L,
714 & PP_FIRST2SWAP_U, PP_LastPanelonDisk_U,
715 & PP_LastPIVRPTRFilled_U,
716 & PIVOT_OPTION, LR_ACTIVATED, IEND_BLR, Inextpiv,
717 & OOC_EFFECTIVE_ON_FRONT, NVSCHUR, PARPIV_T1,
718 & TIPIV
719 & )
720!$ USE OMP_LIB
721 USE mumps_ooc_common
722 IMPLICIT NONE
723 INTEGER, intent(in) :: IBEG_BLOCK, IEND_BLOCK
724 INTEGER, intent(inout), OPTIONAL :: TIPIV(:)
725 INTEGER(8), intent(in) :: LA
726 COMPLEX, intent(inout) :: A(LA)
727 INTEGER, intent(in) :: NFRONT,NASS,N,LIW,INODE,LAST_ROW
728 INTEGER, intent(inout) :: IFLAG,INOPV,NOFFW, NBTINYW
729 INTEGER, intent(inout) :: DET_EXPW, DET_SIGNW
730 COMPLEX, intent(inout) :: DET_MANTW
731 REAL, intent(in) :: UU, SEUIL
732 INTEGER, intent(inout) :: IW(LIW)
733 INTEGER, intent(in) :: IOLDPS
734 INTEGER(8), intent(in) :: POSELT
735 INTEGER KEEP(500)
736 INTEGER(8) KEEP8(150)
737 INTEGER, intent(in) :: LPN_LIST
738 INTEGER, intent(inout) :: PIVNUL_LIST(LPN_LIST)
739 REAL DKEEP(230)
740 INTEGER PP_FIRST2SWAP_L, PP_LastPanelonDisk_L,
741 & pp_lastpivrptrfilled_l,
742 & pp_first2swap_u, pp_lastpanelondisk_u,
743 & pp_lastpivrptrfilled_u
744 INTEGER, intent(in) :: PIVOT_OPTION, IEND_BLR
745 LOGICAL, intent(in) :: LR_ACTIVATED
746 INTEGER, intent(inout) :: Inextpiv
747 LOGICAL, intent(in) :: OOC_EFFECTIVE_ON_FRONT
748 INTEGER, intent(in) :: NVSCHUR
749 INTEGER, intent(in) :: PARPIV_T1
750 INCLUDE 'mumps_headers.h'
751 COMPLEX SWOP
752 INTEGER XSIZE
753 INTEGER(8) :: APOS, IDIAG
754 INTEGER(8) :: J1, J2, JJ, J3
755 INTEGER(8) :: NFRONT8
756 INTEGER ILOC
757 COMPLEX ZERO
758 PARAMETER( ZERO = (0.0e0,0.0e0) )
759 REAL RZERO, RMAX, AMROW, MAX_PREV_in_PARPIV
760 INTEGER(8) :: APOSMAX, APOSROW
761 REAL :: RMAX_NORELAX
762 REAL PIVNUL, ABS_PIVOT
763 COMPLEX FIXA, CSEUIL, PIVOT
764 INTEGER NPIV,IPIV, LRLOC
765 INTEGER NPIVP1,JMAX,J,ISW,ISWPS1
766 INTEGER ISWPS2,KSW, HF, IPIVNUL
767 INTEGER CMUMPS_IXAMAX
768 INTEGER :: ISHIFT, K206
769 INTEGER :: IPIV_SHIFT,IPIV_END
770 INTRINSIC max
771 DATA rzero /0.0e0/
772#if defined(_OPENMP)
773 INTEGER :: NOMP,CHUNK,K361
774#endif
775 INTEGER I_PIVRPTR_L, I_PIVR_L, NBPANELS_L
776 INTEGER I_PIVRPTR_U, I_PIVR_U, NBPANELS_U
777#if defined(_OPENMP)
778 nomp = omp_get_max_threads()
779 k361 = keep(361)
780#endif
781 pivnul = dkeep(1)
782 fixa = cmplx(dkeep(2),kind=kind(fixa))
783 cseuil = cmplx(seuil,kind=kind(cseuil))
784 nfront8 = int(nfront,8)
785 k206 = keep(206)
786 xsize = keep(ixsz)
787 npiv = iw(ioldps+1+xsize)
788 hf = 6 + iw(ioldps+5+xsize)+xsize
789 npivp1 = npiv + 1
790 aposmax = poselt+nfront8*nfront8-1_8
791 IF (ooc_effective_on_front) THEN
792 CALL cmumps_get_ooc_perm_ptr(typef_l, nbpanels_l,
793 & i_pivrptr_l, i_pivr_l,
794 & ioldps+2*nfront+6+iw(ioldps+5+xsize)+xsize,
795 & iw, liw)
796 CALL cmumps_get_ooc_perm_ptr(typef_u, nbpanels_u,
797 & i_pivrptr_u, i_pivr_u,
798 & ioldps+2*nfront+6+iw(ioldps+5+xsize)+xsize,
799 & iw, liw)
800 ENDIF
801 IF ( present(tipiv) ) THEN
802 iloc = npivp1 - ibeg_block + 1
803 tipiv(iloc) = iloc
804 ENDIF
805 IF (inopv .EQ. -1) THEN
806 apos = poselt + nfront8*int(npivp1-1,8) + int(npiv,8)
807 pivot = a(apos)
808 abs_pivot = abs(pivot)
809 idiag = apos
810 CALL cmumps_update_minmax_pivot
811 & ( abs_pivot, dkeep, keep, .true.)
812 IF(abs_pivot.LT.seuil) THEN
813 IF (real(pivot) .GE. rzero) THEN
814 a(apos) = cseuil
815 ELSE
816 a(apos) = -cseuil
817 ENDIF
818 nbtinyw = nbtinyw + 1
819 ELSE IF (keep(258) .NE. 0) THEN
820 CALL cmumps_updatedeter(pivot, det_mantw, det_expw )
821 ENDIF
822 IF (ooc_effective_on_front) THEN
823 IF (keep(251).EQ.0) THEN
824 CALL cmumps_store_perminfo( iw(i_pivrptr_l),
825 & nbpanels_l,
826 & iw(i_pivr_l), nass, npivp1, npivp1,
827 & pp_lastpanelondisk_l,
828 & pp_lastpivrptrfilled_l)
829 ENDIF
830 CALL cmumps_store_perminfo( iw(i_pivrptr_u),
831 & nbpanels_u,
832 & iw(i_pivr_u), nass, npivp1, npivp1,
833 & pp_lastpanelondisk_u,
834 & pp_lastpivrptrfilled_u)
835 ENDIF
836 GO TO 420
837 ENDIF
838 inopv = 0
839 ishift = 0
840 ipiv_end = iend_block
841 IF (k206.GE.1) THEN
842 IF (inextpiv.GT.npivp1.AND.inextpiv.LE.iend_block) THEN
843 ishift = inextpiv - npivp1
844 ENDIF
845 IF ( k206.EQ.1
846 & .OR. (k206 .GT.1 .AND. iend_block.EQ.iend_blr) ) THEN
847 ipiv_end = iend_block + ishift
848 ENDIF
849 ENDIF
850 DO 460 ipiv_shift = npivp1+ishift, ipiv_end
851 IF (ipiv_shift .LE. iend_block) THEN
852 ipiv=ipiv_shift
853 ELSE
854 ipiv = ipiv_shift-iend_block-1+npivp1
855 IF (ibeg_block.EQ.npivp1) THEN
856 EXIT
857 ENDIF
858 ENDIF
859 apos = poselt + nfront8*int(ipiv-1,8) + int(npiv,8)
860 jmax = 1
861 IF ((pivot_option.EQ.0).OR.(uu.EQ.rzero)) THEN
862 abs_pivot = abs(a(apos))
863 IF(abs_pivot.LT.seuil) THEN
864 CALL cmumps_update_minmax_pivot
865 & ( abs_pivot, dkeep, keep, .true.)
866 IF (real(a(apos)) .GE. rzero) THEN
867 a(apos) = cseuil
868 ELSE
869 a(apos) = -cseuil
870 ENDIF
871 nbtinyw = nbtinyw + 1
872 GO TO 420
873 ELSE IF (abs_pivot.EQ.rzero) THEN
874 GO TO 630
875 ENDIF
876 GO TO 380
877 ENDIF
878 amrow = rzero
879 j1 = apos
880 IF (pivot_option.EQ.1 .OR. (lr_activated .AND.
881 & (keep(480).GE.2
882 & ))) THEN
883 j = iend_blr - npiv
884 ELSE
885 j = nass - npiv
886 ENDIF
887 j2 = j1 + j - 1_8
888 jmax = cmumps_ixamax(j,a(j1),1,keep(361))
889 jj = j1 + int(jmax - 1,8)
890 amrow = abs(a(jj))
891 rmax = amrow
892 IF (pivot_option.GE.2) THEN
893 j1 = j2 + 1_8
894 IF (pivot_option.GE.3
895 & ) THEN
896 j2 = apos +
897 & int(- npiv + nfront - 1 - keep(253) - nvschur,8)
898 ELSE
899 j2 = apos +int(- npiv + nass - 1 ,8)
900 ENDIF
901 IF (j2.LT.j1) GO TO 370
902 IF (keep(351).EQ.1) THEN
903!$ CHUNK = max(K361/2,(int(J2-J1)+NOMP-1)/NOMP)
904!$OMP PARALLEL DO schedule(static, CHUNK) PRIVATE(JJ)
905!$OMP& FIRSTPRIVATE(J1,J2)
906!$OMP& REDUCTION(max:RMAX) IF ((J2-J1).GE.K361)
907 DO jj=j1,j2
908 rmax = max(abs(a(jj)),rmax)
909 ENDDO
910!$OMP END PARALLEL DO
911 ELSE
912 DO 360 jj=j1,j2
913 rmax = max(abs(a(jj)),rmax)
914 360 CONTINUE
915 ENDIF
916 370 CONTINUE
917 ENDIF
918 idiag = apos + int(ipiv - npivp1,8)
919 abs_pivot = abs(a(idiag))
920 IF (parpiv_t1.NE.0) THEN
921 rmax_norelax = real(a(aposmax+int(ipiv,8)))
922 ELSE
923 rmax_norelax = rzero
924 ENDIF
925 rmax = max(rmax,rmax_norelax)
926 IF ( rmax .LE. pivnul ) THEN
927 IF (last_row.EQ.nfront) THEN
928 lrloc = last_row -keep(253)-nvschur
929 ELSE
930 lrloc = last_row
931 ENDIF
932 IF (nfront - keep(253) .EQ. nass) THEN
933 IF (iend_block.NE.nass ) THEN
934 GOTO 460
935 ENDIF
936 j1=poselt+int(ipiv-1,8)+int(npiv,8)*nfront8
937 j2=poselt+int(ipiv-1,8)+int(lrloc-1,8)*nfront8
938 ELSE
939 j1=poselt+int(ipiv-1,8)
940 j2=poselt+int(ipiv-1,8)+int(lrloc-1,8)*nfront8
941 ENDIF
942 DO jj=j1, j2, nfront8
943 IF ( abs(a(jj)) .GT. pivnul ) THEN
944 GOTO 460
945 END IF
946 ENDDO
947 IF ((parpiv_t1.NE.0)
948 & .AND.(parpiv_t1.NE.-1)
949 & .AND.(rmax_norelax.LT.0)
950 & .AND.(ipiv.GT.1)) THEN
951 max_prev_in_parpiv = rzero
952 DO jj=1,ipiv-1
953 max_prev_in_parpiv= max( max_prev_in_parpiv,
954 & real(a(aposmax+int(jj,8))) )
955 ENDDO
956 IF (max_prev_in_parpiv.GT.pivnul) THEN
957 aposrow = poselt + nfront8*int(ipiv-1,8)
958 DO jj=1,ipiv-1
959 IF (abs(a(aposrow+jj-1)).GT.pivnul) GOTO 460
960 ENDDO
961 ENDIF
962 ENDIF
963 CALL cmumps_update_minmax_pivot
964 & ( abs_pivot, dkeep, keep, .true.)
965!$OMP ATOMIC CAPTURE
966 keep(109) = keep(109)+1
967 ipivnul = keep(109)
968!$omp END atomic
969 pivnul_list(ipivnul) = iw( ioldps+hf+npiv+ipiv-npivp1 )
970 IF(real(fixa).GT.rzero) THEN
971 IF(real(a(idiag)) .GE. rzero) THEN
972 a(idiag) = fixa
973 ELSE
974 a(idiag) = -fixa
975 ENDIF
976 ELSE
977 j1 = apos
978 j2 = apos +
979 & int(- npiv + nfront - 1 - keep(253),8)
980 DO jj=j1,j2
981 a(jj) = zero
982 ENDDO
983 a(idiag) = -fixa
984 ENDIF
985 jmax = ipiv - npiv
986 GOTO 385
987 ENDIF
988 rmax = max(rmax,abs(rmax_norelax))
989 IF (abs_pivot .GE. uu*rmax .AND.
990 & abs_pivot .GT. max(seuil,tiny(rmax))) THEN
991 jmax = ipiv - npiv
992 GO TO 380
993 ENDIF
994 IF ( .NOT. (amrow .GE. uu*rmax .AND.
995 & amrow .GT. max(seuil,tiny(rmax))) ) GO TO 460
996 noffw = noffw + 1
997 380 CONTINUE
998 IF (k206.GE.1) THEN
999 inextpiv = ipiv + 1
1000 ENDIF
1001 CALL cmumps_update_minmax_pivot
1002 & ( abs(a(apos+int(jmax-1,8))),
1003 & dkeep, keep, .false.)
1004 IF (keep(258) .NE. 0) THEN
1005 CALL cmumps_updatedeter( a(apos+int(jmax-1,8)),
1006 & det_mantw,
1007 & det_expw )
1008 ENDIF
1009 385 CONTINUE
1010 IF (ipiv.EQ.npivp1) GO TO 400
1011 IF (keep(405) .EQ. 0) THEN
1012 keep8(80) = keep8(80)+1
1013 ELSE
1014!$OMP ATOMIC UPDATE
1015 keep8(80) = keep8(80)+1
1016!$OMP END ATOMIC
1017 ENDIF
1018 IF (parpiv_t1.NE.0) THEN
1019 swop = a(aposmax+int(npivp1,8))
1020 a(aposmax+int(npivp1,8)) = a(aposmax+int(ipiv,8))
1021 a(aposmax+int(ipiv,8)) = swop
1022 ENDIF
1023 det_signw = - det_signw
1024 j1 = poselt + int(npiv,8)*nfront8
1025 j2 = j1 + nfront8 - 1_8
1026 j3 = poselt + int(ipiv-1,8)*nfront8
1027 DO 390 jj=j1,j2
1028 swop = a(jj)
1029 a(jj) = a(j3)
1030 a(j3) = swop
1031 j3 = j3 + 1_8
1032 390 CONTINUE
1033 iswps1 = ioldps + hf - 1 + npivp1
1034 iswps2 = ioldps + hf - 1 + ipiv
1035 isw = iw(iswps1)
1036 iw(iswps1) = iw(iswps2)
1037 iw(iswps2) = isw
1038 400 IF (jmax.EQ.1) GO TO 420
1039 det_signw = - det_signw
1040 IF ( present(tipiv) ) THEN
1041 tipiv(iloc) = iloc + jmax - 1
1042 ENDIF
1043 j1 = poselt + int(npiv,8)
1044 j2 = poselt + int(npiv + jmax - 1,8)
1045 DO 410 ksw=1,last_row
1046 swop = a(j1)
1047 a(j1) = a(j2)
1048 a(j2) = swop
1049 j1 = j1 + nfront8
1050 j2 = j2 + nfront8
1051 410 CONTINUE
1052 iswps1 = ioldps + hf - 1 + nfront + npiv + 1
1053 iswps2 = ioldps + hf - 1 + nfront + npiv + jmax
1054 isw = iw(iswps1)
1055 iw(iswps1) = iw(iswps2)
1056 iw(iswps2) = isw
1057 GO TO 420
1058 460 CONTINUE
1059 IF (k206 .GE. 1) THEN
1060 inextpiv=iend_block+1
1061 ENDIF
1062 IF (iend_block.EQ.nass) THEN
1063 inopv = 1
1064 ELSE
1065 inopv = 2
1066 ENDIF
1067 GO TO 430
1068 630 CONTINUE
1069 iflag = -10
1070 GOTO 430
1071 420 CONTINUE
1072 IF (ooc_effective_on_front) THEN
1073 IF (keep(251).EQ.0) THEN
1074 CALL cmumps_store_perminfo( iw(i_pivrptr_l),
1075 & nbpanels_l,
1076 & iw(i_pivr_l), nass, npivp1, ipiv,
1077 & pp_lastpanelondisk_l,
1078 & pp_lastpivrptrfilled_l)
1079 ENDIF
1080 CALL cmumps_store_perminfo( iw(i_pivrptr_u),
1081 & nbpanels_u,
1082 & iw(i_pivr_u), nass, npivp1, npiv+jmax,
1083 & pp_lastpanelondisk_u,
1084 & pp_lastpivrptrfilled_u)
1085 ENDIF
1086 430 CONTINUE
1087 RETURN
1088 END SUBROUTINE cmumps_fac_i
1089 SUBROUTINE cmumps_fac_i_ldlt
1090 & ( nfront,nass,inode,ibeg_block,iend_block,
1091 & iw,liw, a,la, inopv,
1092 & nnegw, nb22t1w, nbtinyw,
1093 & det_expw, det_mantw, det_signw,
1094 & iflag,ioldps,poselt,uu, seuil,keep,keep8,pivsiz,
1095 & dkeep,pivnul_list,lpn_list, xsize,
1096 & pp_first2swap_l, pp_lastpanelondisk,
1097 & pp_lastpivrptrindexfilled,maxfromm,is_maxfromm_avail,
1098 & pivot_option, iend_blr, inextpiv,
1099 & ooc_effective_on_front,
1100 & nvschur, parpiv_t1, lr_activated
1101 & )
1102!$ USE OMP_LIB
1103 USE mumps_ooc_common
1104 IMPLICIT NONE
1105 INTEGER(8) :: POSELT, LA
1106 INTEGER NFRONT,NASS,LIW,INODE,IFLAG,INOPV,
1107 & ioldps
1108 INTEGER, intent(inout) :: NNEGW, NB22T1W, NBTINYW
1109 INTEGER, intent(inout) :: DET_EXPW, DET_SIGNW
1110 COMPLEX, intent(inout) :: DET_MANTW
1111 INTEGER, intent(in) :: IBEG_BLOCK, IEND_BLOCK
1112 INTEGER, intent(in) :: PIVOT_OPTION,IEND_BLR
1113 INTEGER, intent(inout) :: Inextpiv
1114 LOGICAL, intent(in) :: OOC_EFFECTIVE_ON_FRONT
1115 INTEGER PIVSIZ,LPIV, XSIZE
1116 COMPLEX A(LA)
1117 REAL UU, UULOC, SEUIL
1118 INTEGER IW(LIW)
1119 INTEGER KEEP(500)
1120 INTEGER(8) KEEP8(150)
1121 INTEGER LPN_LIST
1122 INTEGER PIVNUL_LIST(LPN_LIST)
1123 REAL DKEEP(230)
1124 INTEGER PP_FIRST2SWAP_L, PP_LastPanelonDisk
1125 INTEGER PP_LastPIVRPTRIndexFilled
1126 REAL, intent(in) :: MAXFROMM
1127 LOGICAL, intent(inout) :: IS_MAXFROMM_AVAIL
1128 INTEGER, intent(in) :: NVSCHUR
1129 INTEGER, intent(in) :: PARPIV_T1
1130 LOGICAL, intent(in) :: LR_ACTIVATED
1131 include 'mpif.h'
1132 INTEGER (8) :: POSPV1,POSPV2,OFFDAG,APOSJ
1133 INTEGER JMAX, LIM, LIM_SWAP
1134 REAL RMAX,AMAX,TMAX, MAX_PREV_in_PARPIV, ABS_PIVOT
1135 REAL RMAX_NORELAX, TMAX_NORELAX, UULOCM1
1136 INTEGER(8) :: APOSMAX, APOSROW
1137 REAL MAXPIV
1138 REAL PIVNUL
1139 REAL MAXFROMM_UPDATED
1140 COMPLEX FIXA, CSEUIL
1141 COMPLEX PIVOT,DETPIV
1142 REAL ABSDETPIV
1143 INCLUDE 'mumps_headers.h'
1144 INTEGER :: HF, IPIVNUL
1145 INTEGER :: J
1146 INTEGER(8) :: APOS, J1, J2, JJ, NFRONT8, KK, J1_ini, JJ_ini
1147 INTEGER :: LDA
1148 INTEGER(8) :: LDA8
1149 INTEGER NPIV,IPIV
1150 INTEGER NPIVP1,K
1151 INTEGER :: ISHIFT, K206, IPIV_SHIFT, IPIV_END
1152 INTRINSIC max
1153 COMPLEX ZERO, ONE
1154 parameter( zero = (0.0e0,0.0e0) )
1155 parameter( one = (1.0e0,1.0e0) )
1156 REAL RZERO,RONE
1157 parameter(rzero=0.0e0, rone=1.0e0)
1158#if defined(_OPENMP)
1159 LOGICAL :: OMP_FLAG
1160 INTEGER :: NOMP, CHUNK, J1_end
1161#endif
1162 INTEGER I_PIVRPTR, I_PIVR, NBPANELS_L
1163!$ NOMP = OMP_GET_MAX_THREADS()
1164 pivnul = dkeep(1)
1165 fixa = cmplx(dkeep(2),kind=kind(fixa))
1166 cseuil = cmplx(seuil,kind=kind(cseuil))
1167 lda = nfront
1168 lda8 = int(lda,8)
1169 nfront8 = int(nfront,8)
1170 k206 = keep(206)
1171 uuloc = uu
1172 IF (uuloc.GT.rzero) THEN
1173 uulocm1 = rone/uuloc
1174 ELSE
1175 uulocm1 = rone
1176 ENDIF
1177 hf = 6 + xsize
1178 IF (keep(50).NE.1 .AND. ooc_effective_on_front) THEN
1179 CALL cmumps_get_ooc_perm_ptr(typef_l, nbpanels_l,
1180 & i_pivrptr, i_pivr, ioldps+2*nfront+6+keep(ixsz),
1181 & iw, liw)
1182 ENDIF
1183 pivsiz = 1
1184 npiv = iw(ioldps+1+xsize)
1185 npivp1 = npiv + 1
1186 aposmax = poselt+lda8*lda8-1_8
1187 IF(inopv .EQ. -1) THEN
1188 apos = poselt + (lda8+1_8) * int(npiv,8)
1189 pospv1 = apos
1190 CALL cmumps_update_minmax_pivot
1191 & ( abs(a(apos)), dkeep, keep, .true.)
1192 IF(abs(a(apos)).LT.seuil) THEN
1193 IF(real(a(apos)) .GE. rzero) THEN
1194 a(apos) = cseuil
1195 ELSE
1196 a(apos) = -cseuil
1197 ENDIF
1198 nbtinyw = nbtinyw + 1
1199 ELSE
1200 IF (keep(258) .NE. 0) THEN
1201 CALL cmumps_updatedeter( a(apos), det_mantw, det_expw )
1202 ENDIF
1203 ENDIF
1204 IF (keep(50).NE.1 .AND. ooc_effective_on_front) THEN
1205 CALL cmumps_store_perminfo( iw(i_pivrptr), nbpanels_l,
1206 & iw(i_pivr), nass, npivp1, npivp1,
1207 & pp_lastpanelondisk,
1208 & pp_lastpivrptrindexfilled)
1209 ENDIF
1210 GO TO 420
1211 ENDIF
1212 inopv = 0
1213 ishift = 0
1214 ipiv_end = iend_block
1215 IF (k206.GE.1) THEN
1216 IF (inextpiv.GT.npivp1.AND.inextpiv.LE.iend_block) THEN
1217 ishift = inextpiv - npivp1
1218 ENDIF
1219 IF ( k206.EQ.1
1220 & .OR. (k206 .GT.1 .AND. iend_block.EQ.iend_blr) ) THEN
1221 ipiv_end = iend_block + ishift
1222 ENDIF
1223 IF (ishift.GT.0.AND.is_maxfromm_avail) THEN
1224 ipiv = npivp1
1225 apos = poselt + lda8*int(ipiv-1,8) + int(npiv,8)
1226 pospv1 = apos + int(ipiv - npivp1,8)
1227 pivot = a(pospv1)
1228 IF ( maxfromm .GT. pivnul ) THEN
1229 IF (parpiv_t1.NE.0) THEN
1230 maxfromm_updated = max
1231 & ( maxfromm,
1232 & abs(real(a(aposmax+int(ipiv,8))))
1233 & )
1234 ELSE
1235 maxfromm_updated = maxfromm
1236 ENDIF
1237 IF ( (abs(pivot) .GE. uuloc*maxfromm_updated).AND.
1238 & abs(pivot) .GT. max(seuil,tiny(maxfromm_updated))
1239 & ) THEN
1240 ishift = 0
1241 ENDIF
1242 ENDIF
1243 ENDIF
1244 IF ( ishift .GT. 0) THEN
1245 is_maxfromm_avail = .false.
1246 ENDIF
1247 ENDIF
1248 DO 460 ipiv_shift = npivp1+ishift, ipiv_end
1249 IF (ipiv_shift .LE. iend_block) THEN
1250 ipiv=ipiv_shift
1251 ELSE
1252 ipiv = ipiv_shift-iend_block-1+npivp1
1253 IF (ibeg_block.EQ.npivp1) THEN
1254 EXIT
1255 ENDIF
1256 ENDIF
1257 apos = poselt + lda8*int(ipiv-1,8) + int(npiv,8)
1258 pospv1 = apos + int(ipiv - npivp1,8)
1259 pivot = a(pospv1)
1260 abs_pivot = abs(pivot)
1261 IF (uuloc.EQ.rzero.OR.pivot_option.EQ.0) THEN
1262 IF(abs_pivot.LT.seuil) THEN
1263 CALL cmumps_update_minmax_pivot
1264 & ( abs_pivot, dkeep, keep, .true.)
1265 IF(real(pivot) .GE. rzero) THEN
1266 a(pospv1) = cseuil
1267 ELSE
1268 a(pospv1) = -cseuil
1269 ENDIF
1270 nbtinyw = nbtinyw + 1
1271 ELSE IF (abs_pivot.EQ.rzero) THEN
1272 GO TO 630
1273 ELSE
1274 CALL cmumps_update_minmax_pivot
1275 & ( abs_pivot, dkeep, keep, .false.)
1276 IF (keep(258) .NE. 0) THEN
1277 CALL cmumps_updatedeter(pivot, det_mantw, det_expw )
1278 ENDIF
1279 ENDIF
1280 GO TO 420
1281 ENDIF
1282 IF ( is_maxfromm_avail ) THEN
1283 IF ( maxfromm .GT. pivnul ) THEN
1284 IF (parpiv_t1.NE.0) THEN
1285 maxfromm_updated = max
1286 & ( maxfromm,
1287 & abs(real(a(aposmax+int(ipiv,8))))
1288 & )
1289 ELSE
1290 maxfromm_updated = maxfromm
1291 ENDIF
1292 IF ( (abs_pivot .GE. uuloc*maxfromm_updated).AND.
1293 & (abs_pivot .GT. max(seuil,tiny(maxfromm_updated)))
1294 & ) THEN
1295 CALL cmumps_update_minmax_pivot
1296 & ( abs_pivot,
1297 & dkeep, keep, .false.)
1298 IF (keep(258) .NE. 0) THEN
1299 CALL cmumps_updatedeter(pivot, det_mantw, det_expw )
1300 ENDIF
1301 GOTO 415
1302 ENDIF
1303 ENDIF
1304 is_maxfromm_avail = .false.
1305 ENDIF
1306 amax = -rone
1307 jmax = 0
1308 IF (pivot_option.EQ.3
1309 & ) THEN
1310 lim = nfront - keep(253)-nvschur
1311 ELSEIF (pivot_option.GE.2
1312 & ) THEN
1313 lim = nass
1314 ELSEIF (pivot_option.GE.1) THEN
1315 lim = iend_blr
1316 ELSE
1317 write(*,*) 'Internal error in FAC_I_LDLT 1x1:',
1318 & pivot_option
1319 CALL mumps_abort()
1320 ENDIF
1321 j1 = apos
1322 j2 = pospv1 - 1_8
1323 DO jj=j1,j2
1324 IF(abs(a(jj)) .GT. amax) THEN
1325 amax = abs(a(jj))
1326 jmax = ipiv - int(pospv1-jj)
1327 ENDIF
1328 ENDDO
1329 j1 = pospv1 + lda8
1330 DO j=1, iend_block - ipiv
1331 IF(abs(a(j1)) .GT. amax) THEN
1332 amax = abs(a(j1))
1333 jmax = ipiv + j
1334 ENDIF
1335 j1 = j1 + lda8
1336 ENDDO
1337 rmax = rzero
1338 j1_ini = j1
1339#if defined(_OPENMP)
1340 j1_end = lim - iend_block
1341 chunk = max(j1_end,1)
1342 IF ( j1_end.GE.keep(360)) THEN
1343 omp_flag = .true.
1344 chunk = max(keep(360)/2,(j1_end+nomp-1)/nomp)
1345 ELSE
1346 omp_flag = .false.
1347 ENDIF
1348#endif
1349!$OMP PARALLEL DO SCHEDULE(static, CHUNK) PRIVATE(J1)
1350!$OMP& REDUCTION(max:RMAX) IF(OMP_FLAG)
1351 DO j=1, lim - iend_block
1352 j1 = j1_ini + int(j-1,8) * lda8
1353 rmax = max(abs(a(j1)),rmax)
1354 ENDDO
1355!$OMP END PARALLEL DO
1356 IF (parpiv_t1.NE.0) THEN
1357 rmax_norelax = real(a(aposmax+int(ipiv,8)))
1358 ELSE
1359 rmax_norelax = rzero
1360 ENDIF
1361 rmax = max(rmax,rmax_norelax)
1362 IF (max(amax,rmax,abs(pivot)).LE.pivnul) THEN
1363 IF ((parpiv_t1.NE.0)
1364 & .AND.(parpiv_t1.NE.-1)
1365 & .AND.(rmax_norelax.LT.0)
1366 & .AND.(ipiv.GT.1)) THEN
1367 max_prev_in_parpiv = rzero
1368 DO jj=1,ipiv-1
1369 max_prev_in_parpiv= max( max_prev_in_parpiv,
1370 & real(a(aposmax+int(jj,8))) )
1371 ENDDO
1372 IF (max_prev_in_parpiv.GT.pivnul) THEN
1373 aposrow = poselt + nfront8*int(ipiv-1,8)
1374 DO jj=1,ipiv-1
1375 IF (abs(a(aposrow+jj-1)).GT.pivnul) THEN
1376 GOTO 460
1377 ENDIF
1378 ENDDO
1379 ENDIF
1380 ENDIF
1381 CALL cmumps_update_minmax_pivot
1382 & ( abs(a(pospv1)), dkeep, keep, .true.)
1383!$OMP ATOMIC CAPTURE
1384 keep(109) = keep(109) + 1
1385 ipivnul = keep(109)
1386!$OMP END ATOMIC
1387 pivnul_list(ipivnul) = iw( ioldps+hf+npiv+ipiv-npivp1 )
1388 IF(real(fixa).GT.rzero) THEN
1389 IF(real(pivot) .GE. rzero) THEN
1390 a(pospv1) = fixa
1391 ELSE
1392 a(pospv1) = -fixa
1393 ENDIF
1394 ELSE
1395 j1 = apos
1396 j2 = pospv1 - 1_8
1397 DO jj=j1,j2
1398 a(jj) = zero
1399 ENDDO
1400 j1 = pospv1 + lda8
1401 DO j=1, iend_block - ipiv
1402 a(j1) = zero
1403 j1 = j1 + lda8
1404 ENDDO
1405 DO j=1,lim - iend_block
1406 a(j1) = zero
1407 j1 = j1 + lda8
1408 ENDDO
1409 a(pospv1) = one
1410 ENDIF
1411 pivot = a(pospv1)
1412 GO TO 415
1413 ENDIF
1414 rmax = max(rmax,abs(rmax_norelax))
1415 IF ( abs(pivot).GE.uuloc*max(rmax,amax)
1416 & .AND. abs(pivot) .GT. max(seuil,tiny(rmax)) ) THEN
1417 CALL cmumps_update_minmax_pivot
1418 & ( abs(pivot),
1419 & dkeep, keep, .false.)
1420 IF (keep(258) .NE.0 ) THEN
1421 CALL cmumps_updatedeter(pivot, det_mantw, det_expw )
1422 ENDIF
1423 GO TO 415
1424 END IF
1425 IF (npivp1.EQ.iend_block) THEN
1426 GOTO 460
1427 ELSE IF (jmax.EQ.0) THEN
1428 GOTO 460
1429 ENDIF
1430 IF (max(abs(pivot),rmax,amax).LE.tiny(rmax)) THEN
1431 GOTO 460
1432 ENDIF
1433 IF (
1434 & (keep(19).NE.0).AND.(max(amax,rmax,abs(pivot)).LE.seuil)
1435 & )
1436 & THEN
1437 GO TO 460
1438 ENDIF
1439 IF (rmax.LT.amax) THEN
1440 j1 = apos
1441 j2 = pospv1 - 1_8
1442 DO jj=j1,j2
1443 IF(int(pospv1-jj) .NE. ipiv-jmax) THEN
1444 rmax = max(rmax,abs(a(jj)))
1445 ENDIF
1446 ENDDO
1447 j1 = pospv1 + lda8
1448 DO j=1,iend_block-ipiv
1449 IF(ipiv+j .NE. jmax) THEN
1450 rmax = max(abs(a(j1)),rmax)
1451 ENDIF
1452 j1 = j1 + lda8
1453 ENDDO
1454 ENDIF
1455 aposj = poselt + int(jmax-1,8)*lda8 + int(npiv,8)
1456 pospv2 = aposj + int(jmax - npivp1,8)
1457 IF (ipiv.LT.jmax) THEN
1458 offdag = aposj + int(ipiv - npivp1,8)
1459 ELSE
1460 offdag = apos + int(jmax - npivp1,8)
1461 END IF
1462 tmax = rzero
1463#if defined(_OPENMP)
1464 j1_end = lim-jmax
1465 chunk = max(j1_end,1)
1466 IF (j1_end.GE.keep(360)) THEN
1467 omp_flag = .true.
1468 chunk = max(keep(360)/2,(j1_end+nomp-1)/nomp)
1469 ELSE
1470 omp_flag = .false.
1471 ENDIF
1472#endif
1473 IF (jmax .LT. ipiv) THEN
1474 jj_ini = pospv2
1475!$OMP PARALLEL DO SCHEDULE(static, CHUNK) IF (OMP_FLAG)
1476!$OMP& PRIVATE(JJ) REDUCTION(max:TMAX)
1477 DO k = 1, lim - jmax
1478 jj = jj_ini+ int(k,8)*nfront8
1479 IF (jmax+k.NE.ipiv) THEN
1480 tmax=max(tmax,abs(a(jj)))
1481 ENDIF
1482 ENDDO
1483!$OMP END PARALLEL DO
1484 DO kk = aposj, pospv2-1_8
1485 tmax = max(tmax,abs(a(kk)))
1486 ENDDO
1487 ELSE
1488 jj_ini = pospv2
1489!$OMP PARALLEL DO SCHEDULE(static, CHUNK) PRIVATE(JJ)
1490!$omp& reduction(max:tmax) IF(omp_flag)
1491 DO k = 1, lim-jmax
1492 jj = jj_ini + int(k,8)*nfront8
1493 tmax=max(tmax,abs(a(jj)))
1494 ENDDO
1495!$OMP END PARALLEL DO
1496 DO kk = aposj, pospv2 - 1_8
1497 IF (kk.NE.offdag) THEN
1498 tmax = max(tmax,abs(a(kk)))
1499 ENDIF
1500 ENDDO
1501 ENDIF
1502 IF (parpiv_t1.NE.0) THEN
1503 tmax_norelax = max(seuil*uulocm1,
1504 & abs(real(a(aposmax+int(jmax,8))))
1505 & )
1506 ELSE
1507 tmax_norelax = seuil*uulocm1
1508 ENDIF
1509 tmax = max(tmax,tmax_norelax)
1510 detpiv = a(pospv1)*a(pospv2) - a(offdag)**2
1511 absdetpiv = abs(detpiv)
1512 IF (seuil.GT.rzero) THEN
1513 IF (sqrt(absdetpiv) .LE. seuil ) THEN
1514 GOTO 460
1515 ENDIF
1516 ENDIF
1517 maxpiv = max(abs(a(pospv1)),abs(a(pospv2)))
1518 IF (maxpiv.EQ.rzero) maxpiv = rone
1519 IF ((abs(a(pospv2))*rmax+amax*tmax)*uuloc.GT.
1520 & absdetpiv .OR. (absdetpiv .EQ. rzero) ) THEN
1521 GO TO 460
1522 ENDIF
1523 IF ((abs(a(pospv1))*tmax+amax*rmax)*uuloc.GT.
1524 & absdetpiv .OR. (absdetpiv.EQ. rzero) ) THEN
1525 GO TO 460
1526 ENDIF
1527 CALL cmumps_update_minmax_pivot
1528 & ( sqrt(absdetpiv),
1529 & dkeep, keep, .false.)
1530 IF (keep(258) .NE.0 ) THEN
1531 CALL cmumps_updatedeter(detpiv, det_mantw, det_expw )
1532 ENDIF
1533 pivsiz = 2
1534 nb22t1w = nb22t1w + 1
1535 415 CONTINUE
1536 IF (k206.GE.1) THEN
1537 inextpiv = max(npivp1+pivsiz, ipiv+1)
1538 ENDIF
1539 DO k=1,pivsiz
1540 IF (pivsiz .EQ. 2) THEN
1541 IF (k==1) THEN
1542 lpiv = min(ipiv,jmax)
1543 ELSE
1544 lpiv = max(ipiv,jmax)
1545 ENDIF
1546 ELSE
1547 lpiv = ipiv
1548 ENDIF
1549 IF (lpiv.EQ.npivp1) GOTO 416
1550 IF (keep(405) .EQ. 0) THEN
1551 keep8(80) = keep8(80)+1
1552 ELSE
1553!$OMP ATOMIC UPDATE
1554 keep8(80) = keep8(80)+1
1555!$OMP END ATOMIC
1556 ENDIF
1557 lim_swap = nfront
1558 CALL cmumps_swap_ldlt( a, la, iw, liw,
1559 & ioldps, npivp1, lpiv, poselt, lim_swap,
1560 & lda, nfront, 1, parpiv_t1, keep(50),
1561 & keep(ixsz), -9999)
1562 416 CONTINUE
1563 IF (keep(50).NE.1 .AND. ooc_effective_on_front) THEN
1564 CALL cmumps_store_perminfo( iw(i_pivrptr), nbpanels_l,
1565 & iw(i_pivr), nass, npivp1, lpiv, pp_lastpanelondisk,
1566 & pp_lastpivrptrindexfilled)
1567 ENDIF
1568 npivp1 = npivp1 + 1
1569 ENDDO
1570 IF(pivsiz .EQ. 2) THEN
1571 a(poselt+(lda8+1_8)*int(npiv,8)+1_8) = detpiv
1572 ENDIF
1573 GOTO 420
1574 460 CONTINUE
1575 IF (k206 .GE. 1) THEN
1576 inextpiv=iend_block+1
1577 ENDIF
1578 IF (iend_block.EQ.nass) THEN
1579 inopv = 1
1580 ELSE
1581 inopv = 2
1582 ENDIF
1583 GO TO 420
1584 630 CONTINUE
1585 pivsiz = 0
1586 iflag = -10
1587 420 CONTINUE
1588 is_maxfromm_avail = .false.
1589 RETURN
1590 END SUBROUTINE cmumps_fac_i_ldlt
1591 SUBROUTINE cmumps_fac_mq_ldlt(IEND_BLOCK,
1592 & NFRONT,NASS,NPIV,INODE,
1593 & A,LA,LDA,
1594 & POSELT,IFINB,PIVSIZ,
1595 & MAXFROMM, IS_MAXFROMM_AVAIL, IS_MAX_USEFUL,
1596 & PARPIV_T1, LAST_ROW, IEND_BLR, NVSCHUR_K253,
1597 & LR_ACTIVATED
1598 & )
1599 IMPLICIT NONE
1600 INTEGER, intent(out):: IFINB
1601 INTEGER, intent(in) :: INODE, NFRONT, NASS, NPIV
1602 INTEGER, intent(in) :: IEND_BLOCK
1603 INTEGER, intent(in) :: LDA
1604 INTEGER(8), intent(in) :: LA
1605 COMPLEX, intent(inout) :: A(LA)
1606 INTEGER, intent(in) :: LAST_ROW
1607 INTEGER, intent(in) :: IEND_BLR
1608 INTEGER(8) :: POSELT
1609 REAL, intent(out) :: MAXFROMM
1610 LOGICAL, intent(out) :: IS_MAXFROMM_AVAIL
1611 LOGICAL, intent(in) :: IS_MAX_USEFUL
1612 INTEGER, intent(in) :: PARPIV_T1
1613 INTEGER, INTENT(in) :: NVSCHUR_K253
1614 LOGICAL, intent(in) :: LR_ACTIVATED
1615 COMPLEX VALPIV
1616 REAL :: MAXFROMMTMP
1617 INTEGER NCB1
1618 INTEGER(8) :: NFRONT8
1619 INTEGER(8) :: LDA8
1620 INTEGER(8) :: K1POS
1621 INTEGER NEL2
1622 COMPLEX ONE, ZERO
1623 COMPLEX A11,A22,A12
1624 INTEGER(8) :: APOS, LPOS, LPOS1, LPOS2
1625 INTEGER(8) :: POSPV1, POSPV2
1626 INTEGER :: PIVSIZ,NPIV_NEW,J2,I
1627 INTEGER(8) :: OFFDAG, OFFDAG_OLD, K1, K2, IROW
1628#if defined(__ve__)
1629 INTEGER(8) :: J2_8, KU1, KU2
1630#else
1631 INTEGER(8) :: IBEG, IEND, JJ_LOC, JJ, ROW_SHIFT
1632 INTEGER(8) :: IBEG_LOC, IEND_LOC
1633#endif
1634 COMPLEX SWOP,DETPIV,MULT1,MULT2
1635 INTEGER(8) :: APOSMAX
1636!$ LOGICAL :: OMP_FLAG
1637 include 'mumps_headers.h'
1638 parameter(one = (1.0e0,0.0e0),
1639 & zero = (0.0e0,0.0e0))
1640 lda8 = int(lda,8)
1641 nfront8 = int(nfront,8)
1642 npiv_new = npiv + pivsiz
1643 ifinb = 0
1644 is_maxfromm_avail = .false.
1645 ncb1 = last_row - iend_block
1646 nel2 = iend_block - npiv_new
1647 IF (nel2.EQ.0) THEN
1648 IF (iend_block.EQ.nass) THEN
1649 ifinb = -1
1650 ELSE
1651 ifinb = 1
1652 ENDIF
1653 ENDIF
1654 maxfromm = 0.0e0
1655 IF(pivsiz .EQ. 1) THEN
1656 apos = poselt + int(npiv,8)*(nfront8 + 1_8)
1657 valpiv = one/a(apos)
1658 lpos = apos + lda8
1659#if defined(__ve__)
1660 IF (nel2+ncb1.GT.0) THEN
1661!$ OMP_FLAG = (NCB1 + NEL2> 300)
1662!$OMP PARALLEL DO PRIVATE(I,K1POS) IF (OMP_FLAG)
1663!NEC$ IVDEP
1664 DO i=1, nel2 + ncb1
1665 k1pos = lpos+ int(i-1,8)*lda8
1666 a(apos+int(i,8))=a(k1pos)
1667 ENDDO
1668!$OMP END PARALLEL DO
1669!$OMP PARALLEL DO PRIVATE(I,K1POS) IF (OMP_FLAG)
1670!NEC$ IVDEP
1671 DO i=1, nel2 + ncb1
1672 k1pos = lpos+ int(i-1,8)*lda8
1673 a(k1pos) = a(k1pos) * valpiv
1674 ENDDO
1675!$OMP END PARALLEL DO
1676 IF (.NOT. is_max_useful) THEN
1677!$ OMP_FLAG = (NCB1 > 300).AND.(NEL2.GE.2)
1678!$omp parallel DO private(i,j2,j2_8,k1pos) IF (omp_flag)
1679!NEC$ IVDEP
1680 DO j2 = 1, nel2
1681 j2_8 = int(j2,8)
1682!NEC$ IVDEP
1683 DO i=j2, nel2 + ncb1
1684 k1pos = lpos+ int(i-1,8)*lda8
1685 a(k1pos+j2_8)=a(k1pos+j2_8)-(a(k1pos)*a(apos+j2_8))
1686 ENDDO
1687 ENDDO
1688!$OMP END PARALLEL DO
1689 ELSE
1690 IF (nel2.GT.0) THEN
1691 maxfrommtmp=0.0e0
1692!$ OMP_FLAG = (NCB1+NEL2 > 300)
1693!$OMP PARALLEL DO PRIVATE(I,K1POS) IF (OMP_FLAG)
1694!$OMP& REDUCTION(max:MAXFROMMTMP)
1695!NEC$ IVDEP
1696 DO i=1, nel2 + ncb1 - nvschur_k253
1697 k1pos = lpos+ int(i-1,8)*lda8
1698 a(k1pos+1_8)=a(k1pos+1_8) - a(k1pos) * a(apos+1_8)
1699 maxfrommtmp=max(maxfrommtmp, abs(a(k1pos+1_8)))
1700 ENDDO
1701 !$OMP END PARALLEL DO
1702 is_maxfromm_avail = .true.
1703 maxfromm=max(maxfromm, maxfrommtmp)
1704 IF (nvschur_k253.GT.0) THEN
1705 DO i= nel2 + ncb1- nvschur_k253 +1, nel2 + ncb1
1706 k1pos = lpos+ int(i-1,8)*lda8
1707 a(k1pos+1_8)=a(k1pos+1_8) - a(k1pos) * a(apos+1_8)
1708 ENDDO
1709 ENDIF
1710 ENDIF
1711 IF (nel2.GT.1) THEN
1712!$ OMP_FLAG = (NCB1+NEL2 > 300).AND.(NEL2.GE.3)
1713!$OMP PARALLEL DO PRIVATE(I,J2,J2_8,K1POS) IF (OMP_FLAG)
1714!NEC$ IVDEP
1715 DO j2 = 2, nel2
1716 j2_8 = int(j2,8)
1717!NEC$ IVDEP
1718 DO i=j2, nel2 + ncb1
1719 k1pos = lpos+ int(i-1,8)*lda8
1720 a(k1pos+j2_8)=a(k1pos+j2_8)-(a(k1pos)*a(apos+j2_8))
1721 ENDDO
1722 ENDDO
1723!$OMP END PARALLEL DO
1724 ENDIF
1725 ENDIF
1726 ENDIF
1727#else
1728 IF (nel2 > 0) THEN
1729 IF (.NOT. is_max_useful) THEN
1730 DO i=1, nel2
1731 k1pos = lpos + int(i-1,8)*lda8
1732 a(apos+int(i,8))=a(k1pos)
1733 a(k1pos) = a(k1pos) * valpiv
1734 DO jj=1_8, int(i,8)
1735 a(k1pos+jj)=a(k1pos+jj) - a(k1pos) * a(apos+jj)
1736 ENDDO
1737 ENDDO
1738 ELSE
1739 is_maxfromm_avail = .true.
1740 DO i=1, nel2
1741 k1pos = lpos + int(i-1,8)*lda8
1742 a(apos+int(i,8))=a(k1pos)
1743 a(k1pos) = a(k1pos) * valpiv
1744 a(k1pos+1_8)=a(k1pos+1_8) - a(k1pos) * a(apos+1_8)
1745 maxfromm=max( maxfromm,abs(a(k1pos+1_8)) )
1746 DO jj = 2_8, int(i,8)
1747 a(k1pos+jj)=a(k1pos+jj) - a(k1pos) * a(apos+jj)
1748 ENDDO
1749 ENDDO
1750 ENDIF
1751 ENDIF
1752 IF (ncb1.GT.0) THEN
1753 IF (.NOT. is_max_useful) THEN
1754!$OMP PARALLEL DO PRIVATE(JJ,K1POS) IF (NCB1 > 300)
1755 DO i=nel2+1, nel2 + ncb1
1756 k1pos = lpos+ int(i-1,8)*lda8
1757 a(apos+int(i,8))=a(k1pos)
1758 a(k1pos) = a(k1pos) * valpiv
1759 DO jj = 1_8, int(nel2,8)
1760 a(k1pos+jj)=a(k1pos+jj) - a(k1pos) * a(apos+jj)
1761 ENDDO
1762 ENDDO
1763!$OMP END PARALLEL DO
1764 ELSE
1765 maxfrommtmp=0.0e0
1766!$ OMP_FLAG = (NCB1-NVSCHUR_K253>300)
1767!$OMP PARALLEL DO PRIVATE(JJ,K1POS)
1768!$OMP& REDUCTION(max:MAXFROMMTMP) IF (OMP_FLAG)
1769 DO i=nel2+1, nel2 + ncb1 - nvschur_k253
1770 k1pos = lpos+ int(i-1,8)*lda8
1771 a(apos+int(i,8))=a(k1pos)
1772 a(k1pos) = a(k1pos) * valpiv
1773 IF (nel2 > 0) THEN
1774 a(k1pos+1_8) = a(k1pos+1_8) - a(k1pos) * a(apos+1_8)
1775 maxfrommtmp=max(maxfrommtmp, abs(a(k1pos+1_8)))
1776 DO jj = 2_8, int(nel2,8)
1777 a(k1pos+jj)=a(k1pos+jj) - a(k1pos) * a(apos+jj)
1778 ENDDO
1779 ENDIF
1780 ENDDO
1781!$OMP END PARALLEL DO
1782 DO i = nel2 + ncb1 - nvschur_k253 + 1, nel2 + ncb1
1783 k1pos = lpos+ int(i-1,8)*lda8
1784 a(apos+int(i,8))=a(k1pos)
1785 a(k1pos) = a(k1pos) * valpiv
1786 DO jj = 1_8, int(nel2,8)
1787 a(k1pos+jj)=a(k1pos+jj) - a(k1pos) * a(apos+jj)
1788 ENDDO
1789 ENDDO
1790 maxfromm=max(maxfromm, maxfrommtmp)
1791 ENDIF
1792 ENDIF
1793#endif
1794 ELSE
1795 pospv1 = poselt + int(npiv,8)*(nfront8 + 1_8)
1796 pospv2 = pospv1 + nfront8 + 1_8
1797 offdag_old = pospv2 - 1_8
1798 offdag = pospv1 + 1_8
1799 swop = a(pospv2)
1800 detpiv = a(offdag)
1801 a22 = a(pospv1)/detpiv
1802 a11 = swop/detpiv
1803 a12 = -a(offdag_old)/detpiv
1804 a(offdag) = a(offdag_old)
1805 a(offdag_old) = zero
1806 lpos1 = pospv2 + lda8 - 1_8
1807 lpos2 = lpos1 + 1_8
1808#if defined(__ve__)
1809#if defined(check)
1810 IF (lda8.NE.nfront8) CALL mumps_abort()
1811 IF (nel2 + ncb1.NE.last_row-npiv_new) CALL mumps_abort()
1812#endif
1813 CALL ccopy(last_row-npiv_new, a(lpos1), lda, a(pospv1+2_8), 1)
1814 CALL ccopy(last_row-npiv_new, a(lpos2), lda, a(pospv2+1_8), 1)
1815!$ OMP_FLAG = (NEL2+NCB1 > 300)
1816!$OMP PARALLEL DO PRIVATE(J2,J2_8,I,K1,K2,KU1,KU2)
1817!$OMP& IF (OMP_FLAG)
1818!NEC$ IVDEP
1819 DO j2=1, nel2 + ncb1
1820 j2_8 = int(j2,8)
1821 ku1 = pospv1 + 2_8 + (j2_8-1_8)
1822 ku2 = pospv2 + 1_8 + (j2_8-1_8)
1823 k1 = lpos1 + (j2_8-1_8)*nfront8
1824 k2 = k1 + 1_8
1825 a(k1) = a11*a(ku1)+a12*a(ku2)
1826 a(k2) = a12*a(ku1)+a22*a(ku2)
1827 ENDDO
1828 IF (nel2.GT.0) THEN
1829!$ OMP_FLAG = (NCB1+NEL2 > 300).AND.(NEL2.GE.2)
1830!$OMP PARALLEL DO PRIVATE(I,J2,J2_8,K1,K2,MULT1,MULT2,IROW)
1831!$OMP& IF (OMP_FLAG)
1832!NEC$ IVDEP
1833 DO j2 = 1,nel2
1834 j2_8 = int(j2,8)
1835 mult1 = -a(pospv1 + 2_8 + j2_8-1_8)
1836 mult2 = -a(pospv2 + 1_8 + j2_8-1_8)
1837!NEC$ IVDEP
1838 DO i= j2, nel2 + ncb1
1839 k1 = lpos1 + (int(i,8)-1_8)*nfront8
1840 k2 = k1 + 1_8
1841 irow = k2 + j2_8
1842 a(irow) = a(irow) + mult1*a(k1) +
1843 & mult2*a(k2)
1844 ENDDO
1845 ENDDO
1846 ENDIF
1847#else
1848 jj = pospv2 + nfront8-1_8
1849 ibeg = jj + 2_8
1850 iend = ibeg
1851 DO j2 = 1,nel2
1852 k1 = jj
1853 k2 = jj+1_8
1854 mult1 = - (a11*a(k1)+a12*a(k2))
1855 mult2 = - (a12*a(k1)+a22*a(k2))
1856 a(pospv1 + 2_8 + (int(j2,8)-1_8)) = a(k1)
1857 a(pospv2 + 1_8 + (int(j2,8)-1_8)) = a(k2)
1858 k1 = pospv1 + 2_8
1859 k2 = pospv2 + 1_8
1860 DO irow = ibeg, iend
1861 a(irow) = a(irow) + mult1*a(k1) + mult2*a(k2)
1862 k1 = k1 + 1_8
1863 k2 = k2 + 1_8
1864 ENDDO
1865 a( jj ) = -mult1
1866 a( jj + 1_8 ) = -mult2
1867 ibeg = ibeg + nfront8
1868 iend = iend + nfront8 + 1_8
1869 jj = jj+nfront8
1870 ENDDO
1871 iend = iend-1_8
1872!$OMP PARALLEL DO PRIVATE(J2, K1, K2, MULT1, MULT2, IROW, JJ_LOC,
1873!$OMP& ROW_SHIFT, IBEG_LOC, IEND_LOC) IF (LAST_ROW-IEND_BLOCK>300)
1874 DO j2 = 1,last_row-iend_block
1875 row_shift = (j2-1_8)*nfront8
1876 jj_loc = jj + row_shift
1877 ibeg_loc = ibeg + row_shift
1878 iend_loc = iend + row_shift
1879 k1 = jj_loc
1880 k2 = jj_loc+1_8
1881 mult1 = - (a11*a(k1)+a12*a(k2))
1882 mult2 = - (a12*a(k1)+a22*a(k2))
1883 a(pospv1 + 2_8 + nel2 + (j2-1_8)) = a(k1)
1884 a(pospv2 + 1_8 + nel2 + (j2-1_8)) = a(k2)
1885 k1 = pospv1 + 2_8
1886 k2 = pospv2 + 1_8
1887 DO irow = ibeg_loc, iend_loc
1888 a(irow) = a(irow) + mult1*a(k1) + mult2*a(k2)
1889 k1 = k1 + 1_8
1890 k2 = k2 + 1_8
1891 ENDDO
1892 a( jj_loc ) = -mult1
1893 a( jj_loc + 1_8 ) = -mult2
1894 ENDDO
1895!$OMP END PARALLEL DO
1896#endif
1897 ENDIF
1898 IF ((is_maxfromm_avail).AND.(nel2.GT.0)) THEN
1899 IF (parpiv_t1.NE.0) THEN
1900 aposmax = poselt+lda8*lda8-1_8 + int(npiv_new+1,8)
1901 maxfromm = max(maxfromm,
1902 & real(a(aposmax))
1903 & )
1904 ENDIF
1905 ENDIF
1906 RETURN
1907 END SUBROUTINE cmumps_fac_mq_ldlt
1908 SUBROUTINE cmumps_fac_sq_ldlt(IBEG_BLOCK,IEND_BLOCK,NPIV,
1909 & NFRONT,NASS,INODE,A,LA,
1910 & LDA,
1911 & POSELT,
1912 & KEEP,KEEP8,
1913 & FIRST_ROW_TRSM, LAST_ROW_TRSM,
1914 & LAST_COL_GEMM, LAST_ROW_GEMM,
1915 & CALL_TRSM, CALL_GEMM, LR_ACTIVATED,
1916 & IW, LIW, OFFSET_IW
1917 & )
1918 IMPLICIT NONE
1919 INTEGER, intent(in) :: NPIV
1920 INTEGER, intent(in) :: NFRONT, NASS, IBEG_BLOCK, IEND_BLOCK
1921 INTEGER(8), intent(in) :: LA
1922 COMPLEX, intent(inout) :: A(LA)
1923 INTEGER, intent(in) :: INODE
1924 INTEGER :: KEEP(500)
1925 INTEGER(8) :: KEEP8(150)
1926 INTEGER(8), intent(in) :: POSELT
1927 INTEGER, intent(in) :: LDA
1928 INTEGER, intent(in) :: LAST_COL_GEMM
1929 INTEGER, intent(in) :: LAST_ROW_GEMM, LAST_ROW_TRSM,
1930 & first_row_trsm
1931 LOGICAL, intent(in) :: CALL_TRSM, CALL_GEMM, LR_ACTIVATED
1932 INTEGER :: OFFSET_IW, LIW
1933 INTEGER :: IW(LIW)
1934 INTEGER(8) :: LDA8
1935 INTEGER NPIV_BLOCK, NEL1
1936 INTEGER NRHS_TRSM
1937 INTEGER(8) :: LPOS, UPOS, APOS
1938 INTEGER IROW
1939 INTEGER Block
1940 INTEGER BLSIZE
1941 COMPLEX ONE, ALPHA
1942 include 'mumps_headers.h'
1943 parameter(one=(1.0e0,0.0e0), alpha=(-1.0e0,0.0e0))
1944 lda8 = int(lda,8)
1945 nel1 = last_col_gemm - iend_block
1946 nrhs_trsm = last_row_trsm-first_row_trsm
1947 npiv_block = npiv - ibeg_block + 1
1948 IF (npiv_block.EQ.0) GO TO 500
1949 IF (nel1.NE.0) THEN
1950 IF (call_trsm) THEN
1951 apos = poselt + lda8*int(ibeg_block-1,8) + int(ibeg_block-1,8)
1952 lpos = poselt + lda8*int(first_row_trsm,8)+int(ibeg_block-1,8)
1953 upos = poselt + lda8*int(ibeg_block-1,8)+int(first_row_trsm,8)
1954 CALL ctrsm('L', 'U', 'T', 'U', npiv_block, nrhs_trsm,
1955 & one, a(apos), lda, a(lpos), lda)
1956 CALL cmumps_fac_ldlt_copy2u_scalel(nrhs_trsm, 1, keep(424),
1957 & nfront, npiv_block, liw, iw, offset_iw, la, a,
1958 & poselt, lpos, upos, apos, .NOT.lr_activated)
1959 ENDIF
1960 IF (call_gemm) THEN
1961#if defined(GEMMT_AVAILABLE)
1962 IF ( keep(421).EQ. -1) THEN
1963 lpos = poselt + lda8*int(iend_block,8) + int(ibeg_block-1,8)
1964 upos = poselt + lda8*int(ibeg_block-1,8) + int(iend_block,8)
1965 apos = poselt + lda8*int(iend_block,8) + int(iend_block,8)
1966 CALL cgemmt( 'U','N','N', nel1,
1967 & npiv_block,
1968 & alpha, a( upos ), lda,
1969 & a( lpos ), lda, one, a( apos ), lda )
1970 ELSE
1971#endif
1972 IF ( last_col_gemm - iend_block > keep(7) ) THEN
1973 blsize = keep(8)
1974 ELSE
1975 blsize = last_col_gemm - iend_block
1976 END IF
1977 IF ( last_col_gemm - iend_block .GT. 0 ) THEN
1978#if defined(SAK_BYROW)
1979 DO irow = iend_block+1, last_col_gemm, blsize
1980 block = min( blsize, last_col_gemm - irow + 1 )
1981 lpos = poselt + int(irow - 1,8) * lda8 +
1982 & int(ibeg_block - 1,8)
1983 upos = poselt + int(ibeg_block - 1,8) * lda8 +
1984 & int(irow - 1,8)
1985 apos = poselt + int(irow - 1,8) * lda8 +
1986 & int(iend_block,8)
1987 CALL cgemm( 'N','N', irow + block - iend_block - 1,
1988 & block, npiv_block,
1989 & alpha, a( upos ), lda,
1990 & a( lpos ), lda, one, a( apos ), lda )
1991 ENDDO
1992#else
1993 DO irow = iend_block+1, last_col_gemm, blsize
1994 block = min( blsize, last_col_gemm - irow + 1 )
1995 lpos = poselt + int( irow - 1,8) * lda8 +
1996 & int(ibeg_block - 1,8)
1997 upos = poselt + int(ibeg_block - 1,8) * lda8 +
1998 & int( irow - 1,8)
1999 apos = poselt + int( irow - 1,8) * lda8 + int( irow - 1,8)
2000 CALL cgemm( 'N','N', block, last_col_gemm - irow + 1,
2001 & npiv_block, alpha, a( upos ), lda,
2002 & a( lpos ), lda, one, a( apos ), lda )
2003 END DO
2004#endif
2005 END IF
2006#if defined(GEMMT_AVAILABLE)
2007 END IF
2008#endif
2009 lpos = poselt + int(last_col_gemm,8)*lda8 + int(ibeg_block-1,8)
2010 upos = poselt + int(ibeg_block-1,8) * lda8 +
2011 & int(iend_block,8)
2012 apos = poselt + int(last_col_gemm,8)*lda8 + int(iend_block,8)
2013 IF (last_row_gemm .GT. last_col_gemm) THEN
2014 CALL cgemm('N', 'N', nel1, last_row_gemm-last_col_gemm,
2015 & npiv_block, alpha, a(upos), lda, a(lpos), lda,
2016 & one, a(apos), lda)
2017 ENDIF
2018 ENDIF
2019 ENDIF
2020 500 CONTINUE
2021 RETURN
2022 END SUBROUTINE cmumps_fac_sq_ldlt
2023 SUBROUTINE cmumps_swap_ldlt( A, LA, IW, LIW,
2024 & IOLDPS, NPIVP1, IPIV, POSELT, LASTROW2SWAP,
2025 & LDA, NFRONT, LEVEL, PARPIV, K50, XSIZE,
2026 & IBEG_BLOCK_TO_SEND )
2027 IMPLICIT NONE
2028 INTEGER(8) :: POSELT, LA
2029 INTEGER LIW, IOLDPS, NPIVP1, IPIV
2030 INTEGER LDA, NFRONT, LEVEL, PARPIV, K50, XSIZE
2031 INTEGER LASTROW2SWAP
2032 COMPLEX A( LA )
2033 INTEGER IW( LIW )
2034 INTEGER, INTENT(IN) :: IBEG_BLOCK_TO_SEND
2035 INCLUDE 'mumps_headers.h'
2036 INTEGER :: IBEG
2037 INTEGER ISW, ISWPS1, ISWPS2, HF
2038 INTEGER(8) :: IDIAG, APOS
2039 INTEGER(8) :: LDA8
2040 COMPLEX SWOP
2041 lda8 = int(lda,8)
2042 apos = poselt + lda8*int(ipiv-1,8) + int(npivp1-1,8)
2043 idiag = apos + int(ipiv - npivp1,8)
2044 hf = 6 + iw( ioldps + 5 + xsize) + xsize
2045 iswps1 = ioldps + hf + npivp1 - 1
2046 iswps2 = ioldps + hf + ipiv - 1
2047 isw = iw(iswps1)
2048 iw(iswps1) = iw(iswps2)
2049 iw(iswps2) = isw
2050 isw = iw(iswps1+nfront)
2051 iw(iswps1+nfront) = iw(iswps2+nfront)
2052 iw(iswps2+nfront) = isw
2053 IF ( level .eq. 2 ) THEN
2054 ibeg = ibeg_block_to_send
2055 CALL cswap( npivp1 - 1 - ibeg + 1,
2056 & a( poselt + int(npivp1-1,8) +
2057 & int(ibeg-1,8) * lda8), lda,
2058 & a( poselt + int(ipiv-1,8) +
2059 & int(ibeg-1,8) * lda8), lda )
2060 END IF
2061 CALL cswap( npivp1-1,
2062 & a( poselt+int(npivp1-1,8) * lda8 ), 1,
2063 & a( poselt + int(ipiv-1,8) * lda8 ), 1 )
2064 CALL cswap( ipiv - npivp1 - 1,
2065 & a( poselt+int(npivp1,8) * lda8 + int(npivp1-1,8) ),
2066 & lda, a( apos + 1_8 ), 1 )
2067 swop = a(idiag)
2068 a(idiag) = a( poselt+int(npivp1-1,8)*lda8+int(npivp1-1,8) )
2069 a( poselt + int(npivp1-1,8)*lda8 + int(npivp1-1,8) ) = swop
2070 IF (lastrow2swap - ipiv.GT.0) THEN
2071 CALL cswap( lastrow2swap - ipiv,
2072 & a( apos + lda8 ), lda,
2073 & a( idiag + lda8 ), lda )
2074 ENDIF
2075 IF (parpiv.NE.0 .AND.k50.EQ.2) THEN
2076 IF ( level .eq. 2 .OR. level.eq.1) THEN
2077 apos = poselt+lda8*lda8-1_8
2078 swop = a(apos+int(npivp1,8))
2079 a(apos+int(npivp1,8))= a(apos+int(ipiv,8))
2080 a(apos+int(ipiv,8)) = swop
2081 ENDIF
2082 ENDIF
2083 RETURN
2084 END SUBROUTINE cmumps_swap_ldlt
2085 SUBROUTINE cmumps_fac_ldlt_copy2u_scalel( IROWMAX, IROWMIN,
2086 & SIZECOPY, LDA, NCOLS, LIW, IW, OFFSET_IW,
2087 & LA, A, POSELT, A_LPOS, A_UPOS, A_DPOS,
2088 & COPY_NEEDED )
2089!$ USE OMP_LIB
2090 INTEGER, INTENT(IN) :: IROWMAX, IROWMIN
2091 INTEGER, INTENT(IN) :: SIZECOPY
2092 INTEGER, INTENT(IN) :: LDA, NCOLS
2093 INTEGER, INTENT(IN) :: LIW
2094 INTEGER, INTENT(IN) :: IW(LIW)
2095 INTEGER, INTENT(IN) :: OFFSET_IW
2096 INTEGER(8), INTENT(IN) :: LA
2097 COMPLEX, INTENT(INOUT) :: A(LA)
2098 INTEGER(8), INTENT(IN) :: POSELT, A_LPOS, A_UPOS, A_DPOS
2099 LOGICAL, INTENT(IN) :: COPY_NEEDED
2100 INTEGER(8) :: LPOS, UPOS
2101 INTEGER(8) :: DPOS, POSPV1, POSPV2, OFFDAG
2102 INTEGER(8) :: LDA8
2103 INTEGER :: IROWEND, IROW, Block2
2104 INTEGER :: I, J
2105 COMPLEX :: MULT1, MULT2, A11, DETPIV, A22, A12
2106 INTEGER :: BLSIZECOPY
2107 COMPLEX :: ONE
2108 parameter(one=(1.0e0,0.0e0))
2109 INTEGER(8) :: LPOSI, UPOSI
2110 LOGICAL :: PIVOT_2X2
2111!$ LOGICAL :: OMP_FLAG
2112!$ INTEGER :: NOMP, CHUNK
2113 lda8 = int(lda,8)
2114 IF (sizecopy.NE.0) THEN
2115 blsizecopy = sizecopy
2116 ELSE
2117 blsizecopy = 250
2118 ENDIF
2119!$ NOMP = OMP_GET_MAX_THREADS()
2120!$ OMP_FLAG = .FALSE.
2121!$ CHUNK = (64/4)
2122!$ IF (NOMP .GT. 1 .AND. NCOLS .GE. 4*CHUNK) THEN
2123!$ OMP_FLAG = .TRUE.
2124!$ CHUNK = max(2*CHUNK, NCOLS/NOMP)
2125!$ ENDIF
2126 DO irowend = irowmax, irowmin, -blsizecopy
2127 block2 = min(blsizecopy, irowend)
2128 irow = irowend - block2 + 1
2129 lpos = a_lpos + int(irow-1,8)*lda8
2130 upos = a_upos + int(irow-1,8)
2131!$OMP PARALLEL DO PRIVATE(PIVOT_2X2, A11, DPOS,
2132!$OMP& POSPV1, POSPV2, OFFDAG, A22, A12, DETPIV, J, MULT1, MULT2
2133!$OMP& , LPOSI, UPOSI
2134!$OMP& ) FIRSTPRIVATE(Block2, LDA, LDA8, LPOS, UPOS, A_DPOS)
2135!$OMP& SCHEDULE(STATIC,CHUNK) IF(OMP_FLAG)
2136 DO i=1, ncols
2137 pivot_2x2 = .false.
2138 IF(iw(offset_iw+i-1) .LE. 0) THEN
2139 pivot_2x2 = .true.
2140 ELSE
2141 IF (i .GT. 1) THEN
2142 IF (iw(offset_iw+i-2) .LE. 0) THEN
2143 cycle
2144 ENDIF
2145 ENDIF
2146 ENDIF
2147 dpos = a_dpos + lda8*int(i-1,8) + int(i-1,8)
2148 IF(.not. pivot_2x2) THEN
2149 a11 = one/a(dpos)
2150 lposi = lpos+int(i-1,8)
2151 IF (copy_needed) THEN
2152 uposi = upos+int(i-1,8)*lda8
2153#if defined(__ve__)
2154!NEC$ IVDEP
2155#endif
2156 DO j = 1, block2
2157 a(uposi+int(j-1,8)) = a(lposi+int(j-1,8)*lda8)
2158 END DO
2159 ENDIF
2160#if defined(__ve__)
2161!NEC$ IVDEP
2162#endif
2163 DO j = 1, block2
2164 a(lposi+int(j-1,8)*lda8) = a(lposi+int(j-1,8)*lda8)*a11
2165 END DO
2166 ELSE
2167 IF (copy_needed) THEN
2168 CALL ccopy(block2, a(lpos+int(i-1,8)),
2169 & lda, a(upos+int(i-1,8)*lda8), 1)
2170 CALL ccopy(block2, a(lpos+int(i,8)),
2171 & lda, a(upos+int(i,8)*lda8), 1)
2172 ENDIF
2173 pospv1 = dpos
2174 pospv2 = dpos + int(lda+1,8)
2175 offdag = pospv1+1_8
2176 a11 = a(pospv1)
2177 a22 = a(pospv2)
2178 a12 = a(offdag)
2179 detpiv = a11*a22 - a12**2
2180 a22 = a11/detpiv
2181 a11 = a(pospv2)/detpiv
2182 a12 = -a12/detpiv
2183#if defined(__ve__)
2184!NEC$ IVDEP
2185#endif
2186 DO j = 1,block2
2187 mult1 = a11*a(lpos+int(j-1,8)*lda8+int(i-1,8))
2188 & + a12*a(lpos+int(j-1,8)*lda8+int(i,8))
2189 mult2 = a12*a(lpos+int(j-1,8)*lda8+int(i-1,8))
2190 & + a22*a(lpos+int(j-1,8)*lda8+int(i,8))
2191 a(lpos+int(j-1,8)*lda8+int(i-1,8)) = mult1
2192 a(lpos+int(j-1,8)*lda8+int(i,8)) = mult2
2193 ENDDO
2194 ENDIF
2195 ENDDO
2196!$OMP END PARALLEL DO
2197 ENDDO
2198 END SUBROUTINE cmumps_fac_ldlt_copy2u_scalel
2199 SUBROUTINE cmumps_fac_ldlt_copyscale_u( IROWMAX, IROWMIN,
2200 & SIZECOPY, LDA, NCOLS, LIW, IW, OFFSET_IW,
2201 & LA, A, POSELT, A_LPOS, A_UPOS, A_DPOS )
2202!$ USE OMP_LIB
2203 INTEGER, INTENT(IN) :: IROWMAX, IROWMIN
2204 INTEGER, INTENT(IN) :: SIZECOPY
2205 INTEGER, INTENT(IN) :: LDA, NCOLS
2206 INTEGER, INTENT(IN) :: LIW
2207 INTEGER, INTENT(IN) :: IW(LIW)
2208 INTEGER, INTENT(IN) :: OFFSET_IW
2209 INTEGER(8), INTENT(IN) :: LA
2210 COMPLEX, INTENT(INOUT) :: A(LA)
2211 INTEGER(8), INTENT(IN) :: POSELT, A_LPOS, A_UPOS, A_DPOS
2212 INTEGER(8) :: LPOS, UPOS
2213 INTEGER(8) :: DPOS, POSPV1, POSPV2, OFFDAG
2214 INTEGER(8) :: LDA8
2215 INTEGER :: IROWEND, IROW, Block2
2216 INTEGER :: I, J
2217 COMPLEX :: MULT1, MULT2, A11, A22, A12
2218 INTEGER :: BLSIZECOPY
2219 COMPLEX :: ONE
2220 parameter(one=(1.0e0,0.0e0))
2221 INTEGER(8) :: LPOSI, UPOSI
2222 LOGICAL :: PIVOT_2X2
2223!$ LOGICAL :: OMP_FLAG
2224!$ INTEGER :: NOMP, CHUNK
2225 lda8 = int(lda,8)
2226 IF (sizecopy.NE.0) THEN
2227 blsizecopy = sizecopy
2228 ELSE
2229 blsizecopy = 250
2230 ENDIF
2231!$ NOMP = OMP_GET_MAX_THREADS()
2232!$ OMP_FLAG = .FALSE.
2233!$ CHUNK = (64/4)
2234!$ IF (NOMP .GT. 1 .AND. NCOLS .GE. 4*CHUNK) THEN
2235!$ OMP_FLAG = .TRUE.
2236!$ CHUNK = max(2*CHUNK, NCOLS/NOMP)
2237!$ ENDIF
2238 DO irowend = irowmax, irowmin, -blsizecopy
2239 block2 = min(blsizecopy, irowend)
2240 irow = irowend - block2 + 1
2241 lpos = a_lpos + int(irow-1,8)*lda8
2242 upos = a_upos + int(irow-1,8)
2243!$OMP PARALLEL DO PRIVATE(PIVOT_2X2, A11, DPOS,
2244!$OMP& POSPV1, POSPV2, OFFDAG, A22, A12, J, MULT1, MULT2
2245!$OMP& , LPOSI, UPOSI
2246!$OMP& ) FIRSTPRIVATE(Block2, LDA, LDA8, LPOS, UPOS, POSELT)
2247!$OMP& SCHEDULE(STATIC,CHUNK) IF(OMP_FLAG)
2248 DO i=1, ncols
2249 pivot_2x2 = .false.
2250 IF(iw(offset_iw+i-1) .LE. 0) THEN
2251 pivot_2x2 = .true.
2252 ELSE
2253 IF (i .GT. 1) THEN
2254 IF (iw(offset_iw+i-2) .LE. 0) THEN
2255 cycle
2256 ENDIF
2257 ENDIF
2258 ENDIF
2259 dpos = a_dpos + lda8*int(i-1,8) + int(i-1,8)
2260 IF(.not. pivot_2x2) THEN
2261 a11 = a(dpos)
2262 lposi = lpos+int(i-1,8)
2263 uposi = upos+int(i-1,8)*lda8
2264#if defined(__ve__)
2265!NEC$ IVDEP
2266#endif
2267 DO j = 1, block2
2268 a(uposi+int(j-1,8)) = a(lposi+int(j-1,8)*lda8)*a11
2269 END DO
2270 ELSE
2271 pospv1 = dpos
2272 pospv2 = dpos + int(lda+1,8)
2273 offdag = pospv1+1_8
2274 a11 = a(pospv1)
2275 a22 = a(pospv2)
2276 a12 = a(offdag)
2277#if defined(__ve__)
2278!NEC$ IVDEP
2279#endif
2280 DO j = 1,block2
2281 mult1 = a11*a(lpos+int(j-1,8)*lda8+int(i-1,8))
2282 & + a12*a(lpos+int(j-1,8)*lda8+int(i,8))
2283 mult2 = a12*a(lpos+int(j-1,8)*lda8+int(i-1,8))
2284 & + a22*a(lpos+int(j-1,8)*lda8+int(i,8))
2285 a(upos+int(i-1,8)*lda8+int(j-1,8)) = mult1
2286 a(upos+int(i,8)*lda8+int(j-1,8)) = mult2
2287 ENDDO
2288 ENDIF
2289 ENDDO
2290!$OMP END PARALLEL DO
2291 ENDDO
2292 RETURN
2293 END SUBROUTINE cmumps_fac_ldlt_copyscale_u
2294 SUBROUTINE cmumps_fac_t_ldlt(NFRONT,NASS,
2295 & IW,LIW,A,LA,
2296 & LDA,
2297 & IOLDPS,POSELT,KEEP,KEEP8,
2298 & POSTPONE_COL_UPDATE, ETATASS,
2299 & TYPEFile, LAFAC, MonBloc, NextPiv2beWritten,
2300 & LIWFAC, MYID, IFLAG, OFFSET_IW, INODE )
2301 USE cmumps_ooc
2302 IMPLICIT NONE
2303 INTEGER NFRONT, NASS,LIW
2304 INTEGER(8) :: LA
2305 COMPLEX A(LA)
2306 INTEGER IW(LIW)
2307 INTEGER KEEP(500)
2308 INTEGER(8) KEEP8(150)
2309 INTEGER(8) :: POSELT
2310 INTEGER LDA
2311 INTEGER IOLDPS, ETATASS
2312 LOGICAL POSTPONE_COL_UPDATE
2313 INTEGER(8) :: LAFAC
2314 INTEGER TYPEFile, NextPiv2beWritten
2315 INTEGER LIWFAC, MYID, IFLAG
2316 TYPE(io_block):: MonBloc
2317 INTEGER IDUMMY
2318 LOGICAL LAST_CALL
2319 INTEGER :: OFFSET_IW
2320 INTEGER, intent(in):: INODE
2321 INCLUDE 'mumps_headers.h'
2322 INTEGER(8) :: UPOS, APOS, LPOS
2323 INTEGER(8) :: LDA8
2324 INTEGER BLSIZE, BLSIZE2, Block, IROW, NPIV, IROWEND
2325 INTEGER I2, I2END, Block2
2326 COMPLEX ONE, ALPHA, BETA, ZERO
2327 parameter(one=(1.0e0,0.0e0), alpha=(-1.0e0,0.0e0))
2328 parameter(zero=(0.0e0,0.0e0))
2329 lda8 = int(lda,8)
2330 IF (etatass.EQ.1) THEN
2331 beta = zero
2332 ELSE
2333 beta = one
2334 ENDIF
2335 IF ( nfront - nass > keep(58) ) THEN
2336 IF ( nfront - nass > keep(57) ) THEN
2337 blsize = keep(58)
2338 ELSE
2339 blsize = (nfront - nass)/2
2340 END IF
2341 ELSE
2342 blsize = nfront - nass
2343 END IF
2344 blsize2 = keep(218)
2345 npiv = iw( ioldps + 1 + keep(ixsz))
2346 IF ( nfront - nass .GT. 0 ) THEN
2347 IF ( postpone_col_update ) THEN
2348 lpos = poselt + lda8 * int(nass,8)
2349 CALL ctrsm( 'L', 'U', 'T', 'U',
2350 & npiv, nfront-nass, one,
2351 & a( poselt ), lda,
2352 & a( lpos ), lda )
2353 ENDIF
2354#if defined(GEMMT_AVAILABLE)
2355 IF ( keep(421).EQ. -1) THEN
2356 lpos = poselt + int(nass,8)*lda8
2357 upos = poselt + int(nass,8)
2358 apos = poselt + int(nass,8)*lda8 + int(nass,8)
2359 IF (postpone_col_update) THEN
2360 CALL cmumps_fac_ldlt_copy2u_scalel( nfront - nass, 1,
2361 & keep(424), nfront, npiv,
2362 & liw, iw, offset_iw, la, a, poselt, lpos, upos,
2363 & poselt, .true. )
2364 ENDIF
2365 CALL cgemmt('U', 'N', 'N', nfront-nass, npiv,
2366 & alpha, a( upos ), lda,
2367 & a( lpos ), lda,
2368 & beta,
2369 & a( apos ), lda )
2370 ELSE
2371#endif
2372 DO irowend = nfront - nass, 1, -blsize
2373 block = min( blsize, irowend )
2374 irow = irowend - block + 1
2375 lpos = poselt + int(nass,8)*lda8 + int(irow-1,8) * lda8
2376 apos = poselt + int(nass,8)*lda8 + int(irow-1,8) * lda8 +
2377 & int(nass + irow - 1,8)
2378 upos = poselt + int(nass,8)
2379 IF (.NOT. postpone_col_update) THEN
2380 upos = poselt + int(nass + irow - 1,8)
2381 ENDIF
2382 IF (postpone_col_update) THEN
2383 CALL cmumps_fac_ldlt_copy2u_scalel( block, 1,
2384 & keep(424), nfront, npiv,
2385 & liw, iw, offset_iw, la, a, poselt, lpos, upos,
2386 & poselt, .true. )
2387 ENDIF
2388 DO i2end = block, 1, -blsize2
2389 block2 = min(blsize2, i2end)
2390 i2 = i2end - block2+1
2391 CALL cgemm('N', 'N', block2, block-i2+1, npiv, alpha,
2392 & a(upos+int(i2-1,8)), lda,
2393 & a(lpos+int(i2-1,8)*lda8), lda,
2394 & beta,
2395 & a(apos + int(i2-1,8) + int(i2-1,8)*lda8), lda)
2396 IF (keep(201).EQ.1) THEN
2397 IF (nextpiv2bewritten.LE.npiv) THEN
2398 last_call=.false.
2399 CALL cmumps_ooc_io_lu_panel(
2400 & strat_try_write, typefile,
2401 & a(poselt), lafac, monbloc,
2402 & nextpiv2bewritten, idummy,
2403 & iw(ioldps), liwfac, myid,
2404 & keep8(31),
2405 & iflag,last_call )
2406 IF (iflag .LT. 0 ) RETURN
2407 ENDIF
2408 ENDIF
2409 ENDDO
2410 IF ( nfront - nass - irow + 1 - block > 0 ) THEN
2411 CALL cgemm( 'N', 'N', block, nfront-nass-block-irow+1, npiv,
2412 & alpha, a( upos ), lda,
2413 & a( lpos + lda8 * int(block,8) ), lda,
2414 & beta,
2415 & a( apos + lda8 * int(block,8) ), lda )
2416 ENDIF
2417 END DO
2418#if defined(GEMMT_AVAILABLE)
2419 END IF
2420#endif
2421 IF ( (postpone_col_update).AND.(nass-npiv.GT.0) ) THEN
2422 lpos = poselt + int(npiv,8)*lda8
2423 upos = poselt + int(npiv,8)
2424 CALL cmumps_fac_ldlt_copyscale_u( nass-npiv, 1,
2425 & keep(424), nfront, npiv,
2426 & liw, iw, offset_iw, la, a, poselt, lpos, upos, poselt)
2427 lpos = poselt + lda8 * int(nass,8)
2428 CALL cgemm('N', 'N', nass-npiv, nfront-nass, npiv, alpha,
2429 & a( poselt + int(npiv,8)), lda,
2430 & a( lpos ), lda,
2431 & beta,
2432 & a( lpos + int(npiv,8) ), lda)
2433 ENDIF
2434 END IF
2435 RETURN
2436 END SUBROUTINE cmumps_fac_t_ldlt
2437 SUBROUTINE cmumps_store_perminfo( PIVRPTR, NBPANELS, PIVR, NASS,
2438 & K, P, LastPanelonDisk,
2439 & LastPIVRPTRIndexFilled )
2440 IMPLICIT NONE
2441 INTEGER, intent(in) :: NBPANELS, NASS, K, P
2442 INTEGER, intent(inout) :: PIVRPTR(NBPANELS), PIVR(NASS)
2443 INTEGER LastPanelonDisk, LastPIVRPTRIndexFilled
2444 INTEGER I
2445 IF ( LastPanelonDisk+1 > NBPANELS ) THEN
2446 WRITE(*,*) "INTERNAL ERROR IN CMUMPS_STORE_PERMINFO!"
2447 WRITE(*,*) "NASS=",nass,"PIVRPTR=",pivrptr(1:nbpanels)
2448 WRITE(*,*) "K=",k, "P=",p, "LastPanelonDisk=",lastpanelondisk
2449 WRITE(*,*) "LastPIVRPTRIndexFilled=", lastpivrptrindexfilled
2450 CALL mumps_abort()
2451 ENDIF
2452 pivrptr(lastpanelondisk+1) = k + 1
2453 IF (lastpanelondisk.NE.0) THEN
2454 pivr(k - pivrptr(1) + 1) = p
2455 DO i = lastpivrptrindexfilled + 1, lastpanelondisk
2456 pivrptr(i)=pivrptr(lastpivrptrindexfilled)
2457 ENDDO
2458 ENDIF
2459 lastpivrptrindexfilled = lastpanelondisk + 1
2460 RETURN
2461 END SUBROUTINE cmumps_store_perminfo
2462 SUBROUTINE cmumps_update_minmax_pivot
2463 & ( diag, dkeep, keep, nullpivot)
2464!$ USE OMP_LIB
2465 IMPLICIT NONE
2466 REAL, INTENT(IN) :: DIAG
2467 REAL, INTENT(INOUT) :: DKEEP(230)
2468 LOGICAL, INTENT(IN) :: NULLPIVOT
2469 INTEGER, INTENT(IN) :: KEEP(500)
2470 IF (KEEP(405).EQ.0) THEN
2471 DKEEP(21) = max(dkeep(21), diag)
2472 dkeep(19) = min(dkeep(19), diag)
2473 IF (.NOT.nullpivot) THEN
2474 dkeep(20) = min(dkeep(20), diag)
2475 ENDIF
2476 ELSE
2477!$OMP ATOMIC UPDATE
2478 dkeep(21) = max(dkeep(21), diag)
2479!$OMP END ATOMIC
2480!$OMP ATOMIC UPDATE
2481 dkeep(19) = min(dkeep(19), diag)
2482!$OMP END ATOMIC
2483 IF (.NOT.nullpivot) THEN
2484!$OMP ATOMIC UPDATE
2485 dkeep(20) = min(dkeep(20), diag)
2486!$OMP END ATOMIC
2487 ENDIF
2488 ENDIF
2489 RETURN
2490 END SUBROUTINE cmumps_update_minmax_pivot
2491 SUBROUTINE cmumps_get_size_schur_in_front (
2492 & N, NCB, SIZE_SCHUR, ROW_INDICES, PERM,
2493 & NVSCHUR
2494 & )
2495 IMPLICIT NONE
2496 INTEGER, intent(in) :: N, NCB, SIZE_SCHUR
2497 INTEGER, intent(in) :: ROW_INDICES(NCB), PERM(N)
2498 INTEGER, intent(out):: NVSCHUR
2499 INTEGER :: I, IPOS, IBEG_SCHUR
2500 IBEG_SCHUR = n - size_schur +1
2501 nvschur = 0
2502 ipos = ncb
2503 DO i= ncb,1,-1
2504 IF (abs(row_indices(i)).LE.n) THEN
2505 IF (perm(row_indices(i)).LT.ibeg_schur) EXIT
2506 ENDIF
2507 ipos = ipos -1
2508 ENDDO
2509 nvschur = ncb-ipos
2510 RETURN
2511 END SUBROUTINE cmumps_get_size_schur_in_front
2512 END MODULE cmumps_fac_front_aux_m
cmplx
float cmplx[2]
Definition pblas.h:136
mumps_abort
#define mumps_abort
Definition VE_Metis.h:25
cmumps_updatedeter
subroutine cmumps_updatedeter(piv, deter, nexp)
Definition cfac_determinant.F:15
cmumps_get_ooc_perm_ptr
subroutine cmumps_get_ooc_perm_ptr(typef, nbpanels, i_pivptr, i_piv, ipos, iw, liw)
Definition cooc_panel_piv.F:133
cswap
subroutine cswap(n, cx, incx, cy, incy)
CSWAP
Definition cswap.f:81
ccopy
subroutine ccopy(n, cx, incx, cy, incy)
CCOPY
Definition ccopy.f:81
cgeru
subroutine cgeru(m, n, alpha, x, incx, y, incy, a, lda)
CGERU
Definition cgeru.f:130
ctrsm
subroutine ctrsm(side, uplo, transa, diag, m, n, alpha, a, lda, b, ldb)
CTRSM
Definition ctrsm.f:180
cgemm
subroutine cgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
CGEMM
Definition cgemm.f:187
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21
cmumps_fac_front_aux_m
Definition cfac_front_aux.F:14
cmumps_fac_front_aux_m::cmumps_fac_sq_ldlt
subroutine cmumps_fac_sq_ldlt(ibeg_block, iend_block, npiv, nfront, nass, inode, a, la, lda, poselt, keep, keep8, first_row_trsm, last_row_trsm, last_col_gemm, last_row_gemm, call_trsm, call_gemm, lr_activated, iw, liw, offset_iw)
Definition cfac_front_aux.F:1918
cmumps_fac_front_aux_m::cmumps_fac_mq_ldlt
subroutine cmumps_fac_mq_ldlt(iend_block, nfront, nass, npiv, inode, a, la, lda, poselt, ifinb, pivsiz, maxfromm, is_maxfromm_avail, is_max_useful, parpiv_t1, last_row, iend_blr, nvschur_k253, lr_activated)
Definition cfac_front_aux.F:1599
cmumps_fac_front_aux_m::cmumps_get_size_schur_in_front
subroutine cmumps_get_size_schur_in_front(n, ncb, size_schur, row_indices, perm, nvschur)
Definition cfac_front_aux.F:2495
cmumps_fac_front_aux_m::cmumps_fac_n
subroutine cmumps_fac_n(nfront, nass, iw, liw, a, la, ioldps, poselt, ifinb, xsize, keep, maxfromn, is_maxfromn_avail, nvschur)
Definition cfac_front_aux.F:239
cmumps_fac_front_aux_m::cmumps_fac_t
subroutine cmumps_fac_t(a, la, npivb, nfront, npiv, nass, poselt)
Definition cfac_front_aux.F:413
cmumps_fac_front_aux_m::cmumps_fac_sq
subroutine cmumps_fac_sq(ibeg_block, iend_block, npiv, nfront, last_row, last_col, a, la, poselt, first_col, call_ltrsm, call_utrsm, call_gemm, with_comm_thread, lr_activated)
Definition cfac_front_aux.F:443
cmumps_fac_front_aux_m::cmumps_fac_h
subroutine cmumps_fac_h(nfront, nass, iw, liw, a, la, inopv, noffw, det_expw, det_mantw, det_signw, ioldps, poselt, uu, seuil, keep, keep8, dkeep, pp_first2swap_l, pp_lastpanelondisk_l, pp_lastpivrptrfilled_l, pp_first2swap_u, pp_lastpanelondisk_u, pp_lastpivrptrfilled_u, maxfromn, is_maxfromn_avail, inextpiv, ooc_effective_on_front, nvschur)
Definition cfac_front_aux.F:26
cmumps_fac_front_aux_m::cmumps_fac_pt_setlock427
subroutine cmumps_fac_pt_setlock427(k427_out, k427, k405, k222, nel1, nass)
Definition cfac_front_aux.F:339
cmumps_fac_front_aux_m::cmumps_fac_i
subroutine cmumps_fac_i(nfront, nass, last_row, ibeg_block, iend_block, n, inode, iw, liw, a, la, inopv, noffw, nbtinyw, det_expw, det_mantw, det_signw, iflag, ioldps, poselt, uu, seuil, keep, keep8, dkeep, pivnul_list, lpn_list pp_first2swap_l, pp_lastpanelondisk_l, pp_lastpivrptrfilled_l, pp_first2swap_u, pp_lastpanelondisk_u, pp_lastpivrptrfilled_u, pivot_option, lr_activated, iend_blr, inextpiv, ooc_effective_on_front, nvschur, parpiv_t1, tipiv)
Definition cfac_front_aux.F:720
cmumps_fac_front_aux_m::cmumps_store_perminfo
subroutine cmumps_store_perminfo(pivrptr, nbpanels, pivr, nass, k, p, lastpanelondisk, lastpivrptrindexfilled)
Definition cfac_front_aux.F:2440
cmumps_fac_front_aux_m::cmumps_fac_ldlt_copy2u_scalel
subroutine cmumps_fac_ldlt_copy2u_scalel(irowmax, irowmin, sizecopy, lda, ncols, liw, iw, offset_iw, la, a, poselt, a_lpos, a_upos, a_dpos, copy_needed)
Definition cfac_front_aux.F:2089
cmumps_fac_front_aux_m::cmumps_fac_ldlt_copyscale_u
subroutine cmumps_fac_ldlt_copyscale_u(irowmax, irowmin, sizecopy, lda, ncols, liw, iw, offset_iw, la, a, poselt, a_lpos, a_upos, a_dpos)
Definition cfac_front_aux.F:2202
cmumps_fac_front_aux_m::cmumps_fac_fr_update_cbrows
subroutine cmumps_fac_fr_update_cbrows(inode, nfront, nass, call_utrsm, a, la, lafac, poselt, iw, liw, ioldps, monbloc, myid, noffw, det_expw, det_mantw, det_signw, liwfac, pp_first2swap_l, pp_first2swap_u, lnextpiv2bewritten, unextpiv2bewritten, pp_lastpivrptrfilled_l, pp_lastpivrptrfilled_u xsize, seuil, uu, dkeep, keep8, keep, iflag, ooc_effective_on_front, nvschur)
Definition cfac_front_aux.F:628
cmumps_fac_front_aux_m::cmumps_fac_i_ldlt
subroutine cmumps_fac_i_ldlt(nfront, nass, inode, ibeg_block, iend_block, iw, liw, a, la, inopv, nnegw, nb22t1w, nbtinyw, det_expw, det_mantw, det_signw, iflag, ioldps, poselt, uu, seuil, keep, keep8, pivsiz, dkeep, pivnul_list, lpn_list, xsize, pp_first2swap_l, pp_lastpanelondisk, pp_lastpivrptrindexfilled, maxfromm, is_maxfromm_avail, pivot_option, iend_blr, inextpiv, ooc_effective_on_front, nvschur, parpiv_t1, lr_activated)
Definition cfac_front_aux.F:1102
cmumps_fac_front_aux_m::cmumps_fac_p
subroutine cmumps_fac_p(a, la, nfront, npiv, nass, poselt, call_utrsm, keep, inode, call_ooc, iwfac, liwfac, lafac, monbloc, myid, keep8, lnextpiv2bewritten, unextpiv2bewritten, iflag)
Definition cfac_front_aux.F:353
cmumps_fac_front_aux_m::cmumps_swap_ldlt
subroutine cmumps_swap_ldlt(a, la, iw, liw, ioldps, npivp1, ipiv, poselt, lastrow2swap, lda, nfront, level, parpiv, k50, xsize, ibeg_block_to_send)
Definition cfac_front_aux.F:2027
cmumps_fac_front_aux_m::cmumps_fac_t_ldlt
subroutine cmumps_fac_t_ldlt(nfront, nass, iw, liw, a, la, lda, ioldps, poselt, keep, keep8, postpone_col_update, etatass, typefile, lafac, monbloc, nextpiv2bewritten, liwfac, myid, iflag, offset_iw, inode)
Definition cfac_front_aux.F:2301
cmumps_fac_front_aux_m::cmumps_fac_mq
subroutine cmumps_fac_mq(ibeg_block, iend_block, nfront, nass, npiv, last_col, a, la, poselt, ifinb, lr_activated)
Definition cfac_front_aux.F:573
cmumps_fac_front_aux_m::cmumps_update_minmax_pivot
subroutine cmumps_update_minmax_pivot(diag, dkeep, keep, nullpivot)
Definition cfac_front_aux.F:2464
cmumps_ooc
Definition cmumps_ooc.F:14
cmumps_ooc::cmumps_ooc_io_lu_panel
subroutine, public cmumps_ooc_io_lu_panel(strat, typefile, afac, lafac, monbloc, lnextpiv2bewritten, unextpiv2bewritten, iw, liwfac, myid, filesize, ierr, last_call)
Definition cmumps_ooc.F:3032
cmumps_ooc::typef_both_lu
integer, parameter, public typef_both_lu
Definition cmumps_ooc.F:64
mumps_ooc_common
Definition mumps_ooc_common.F:14
mumps_ooc_common::typef_u
integer, public typef_u
Definition mumps_ooc_common.F:18
mumps_ooc_common::strat_try_write
integer, public strat_try_write
Definition mumps_ooc_common.F:45
mumps_ooc_common::typef_l
integer, public typef_l
Definition mumps_ooc_common.F:18