OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
imp_mumps.F
Go to the documentation of this file.
1Copyright> OpenRadioss
2Copyright> Copyright (C) 1986-2025 Altair Engineering Inc.
3Copyright>
4Copyright> This program is free software: you can redistribute it and/or modify
5Copyright> it under the terms of the GNU Affero General Public License as published by
6Copyright> the Free Software Foundation, either version 3 of the License, or
7Copyright> (at your option) any later version.
8Copyright>
9Copyright> This program is distributed in the hope that it will be useful,
10Copyright> but WITHOUT ANY WARRANTY; without even the implied warranty of
11Copyright> MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12Copyright> GNU Affero General Public License for more details.
13Copyright>
14Copyright> You should have received a copy of the GNU Affero General Public License
15Copyright> along with this program. If not, see <https://www.gnu.org/licenses/>.
16Copyright>
17Copyright>
18Copyright> Commercial Alternative: Altair Radioss Software
19Copyright>
20Copyright> As an alternative to this open-source version, Altair also offers Altair Radioss
21Copyright> software under a commercial license. Contact Altair to discuss further if the
22Copyright> commercial version may interest you: https://www.altair.com/radioss/.
23!||====================================================================
24!|| imp_mumps1 ../engine/source/implicit/imp_mumps.F
25!||--- called by ------------------------------------------------------
26!|| imp_solv ../engine/source/implicit/imp_solv.F
27!||--- calls -----------------------------------------------------
28!|| mumps_set2 ../engine/source/implicit/imp_mumps.F
29!|| print_stiff_mat ../engine/source/implicit/imp_mumps.F
30!|| spmd_cddl ../engine/source/mpi/implicit/imp_spmd.F
31!|| spmd_inf_g ../engine/source/mpi/implicit/imp_spmd.F
32!|| spmd_mumps_count ../engine/source/mpi/implicit/imp_spmd.F
33!|| spmd_mumps_deal ../engine/source/mpi/implicit/imp_spmd.F
34!|| spmd_mumps_gath ../engine/source/mpi/implicit/imp_spmd.F
35!|| spmd_mumps_ini ../engine/source/mpi/implicit/imp_spmd.F
36!|| tmpenvf ../engine/source/system/tmpenv_c.c
37!||--- uses -----------------------------------------------------
38!|| imp_intm ../engine/share/modules/imp_intm.f
39!|| imp_kbcs ../engine/share/modules/impbufdef_mod.F
40!|| message_mod ../engine/share/message_module/message_mod.F
41!|| spmd_comm_world_mod ../engine/source/mpi/spmd_comm_world.f90
42!||====================================================================
43 SUBROUTINE imp_mumps1(NDDL0 , NNZK0 , NDDL , NNZK , NNMAX ,
44 . NODGLOB, IDDL , NDOF , INLOC, IKC ,
45 . IADK , JDIK , DIAG_K, LT_K , IAD_ELEM,
46 . FR_ELEM, MUMPS_PAR, CDDLP , IADI , JDII ,
47 . ITOK , DIAG_I , LT_I , NDDLI, NNZI ,
48 . IMPRINT, IT)
49C-----------------------------------------------
50C M o d u l e s
51C-----------------------------------------------
52 USE imp_intm
53 USE imp_kbcs
54 USE message_mod
55C-----------------------------------------------
56C I m p l i c i t T y p e s
57C-----------------------------------------------
58 USE spmd_comm_world_mod, ONLY : spmd_comm_world
59#include "implicit_f.inc"
60C-----------------------------------------------
61C C o m m o n B l o c k s
62C-----------------------------------------------
63#if defined(MUMPS5)
64#include "dmumps_struc.h"
65#endif
66#include "impl1_c.inc"
67#include "task_c.inc"
68#include "units_c.inc"
69#include "com01_c.inc"
70C-----------------------------------------------
71C D u m m y A r g u m e n t s
72C-----------------------------------------------
73 INTEGER NDDL0, NNZK0, NDDL, NNZK, NNMAX, NODGLOB(*), IDDL(*),
74 . NDOF(*), INLOC(*), IKC(*), IADK(*), JDIK(*),
75 . IAD_ELEM(2,*), FR_ELEM(*), CDDLP(*), IADI(*), JDII(*),
76 . ITOK(*), NDDLI, NNZI,IMPRINT,TLEN, IT
77 my_real
78 . diag_k(*), lt_k(*), diag_i(*), lt_i(*)
79#ifdef MUMPS5
80 TYPE(dmumps_struc) MUMPS_PAR
81#else
82 ! Fake declaration as DMUMPS_STRUC is shipped with MUMPS
83 INTEGER MUMPS_PAR
84#endif
85#ifdef MUMPS5
86C----------------------------------------------
87C L o c a l V a r i a b l e s
88C-----------------------------------------------
89 INTEGER I, J , N
90 INTEGER NDDLG0, NNZKG0, NDDLG, NNZKG, NNMAXG,
91 . NDDL0P(NSPMD), NNZK0P(NSPMD), NDDLP(NSPMD),
92 . NNZKP(NSPMD), NNMAXP(NSPMD), NKLOC,
93 . NKFRONT, NKFLOC, NZLOC, NNZ, NZP(NSPMD-1), IACTI(NDDL),
94 . nnzt
95 INTEGER, DIMENSION(:,:), ALLOCATABLE :: ITK
96 my_real
97 . , DIMENSION(:), ALLOCATABLE :: RTK
98C
99 IF (NDDLI==0) nnzi=0
100C Desallocations si necessaire
101
102 CALL spmd_mumps_deal(mumps_par)
103C
104 CALL spmd_mumps_ini(mumps_par, 1)
105C
106
107 IF (ncycle==1.AND.inconv==1) THEN
108 mumps_par%ICNTL(3) = iout
109 ELSE
110 mumps_par%ICNTL(3) = -1
111 ENDIF
112C--Level of info to be printed with user input /MUMPS/MSGLV/n
113 IF(m_msg > 0) THEN
114 mumps_par%ICNTL(3) = iout ! standard output
115 mumps_par%ICNTL(4) = m_msg ! max lev of info
116 ENDIF
117C
118 IF (m_order==0) THEN
119 mumps_par%ICNTL(7) = 7 ! Automatic choice of ordering
120 ELSE
121 mumps_par%ICNTL(7) = m_order
122 END IF
123
124 mumps_par%ICNTL(13) = 1 ! Disable scalapack for the root matrix
125C uncomment to set out of core
126C MUMPS_PAR%ICNTL(22)=1
127
128 IF (m_ocore > 0) THEN
129 CALL tmpenvf(mumps_par%OOC_TMPDIR,tlen)
130 mumps_par%ICNTL(22)=1
131 ENDIF
132
133 IF (nspmd>1) THEN
134C LMEMV is the memory on the host (i.e. node)
135C MUMPS_PAR%ICNTL(23)=LMEMV/NSPMD_PER_NODE
136 IF (imumpsd==1) THEN
137 mumps_par%ICNTL(18)=3
138 ELSEIF (imumpsd==2) THEN
139 mumps_par%ICNTL(18)=0
140 ENDIF
141 IF (idtc==3) mumps_par%ICNTL(13)=1
142C
143 nddlg0 = nddl0
144 nnzkg0 = nnzk0
145 nddlg = nddl
146 nnzkg = nnzk
147 nnmaxg = nnmax
148 CALL spmd_inf_g(
149 1 nddlg0 ,nnzkg0 ,nddlg ,nnzkg ,nnmaxg ,
150 2 nddl0p ,nnzk0p ,nddlp ,nnzkp ,nnmaxp )
151C
152 CALL spmd_cddl(nddl, nodglob, iddl, ndof, cddlp,
153 . inloc, ikc, nddlg, nddlp)
154C
155 nnzt = nddl+nnzk+nnzi+nz_sl+nz_si
156 ALLOCATE(itk(2,nnzt),rtk(nnzt))
157C
158 DO i=1,nddl
159 iacti(i)=i
160 ENDDO
161C
162 CALL mumps_set2(
163 . iadk, jdik, diag_k, lt_k, cddlp,
164 . nkloc, nkfront, itk, rtk, iddl,
165 . inloc, iad_elem, fr_elem, ndof, ikc,
166 . nddl, nnzk, iacti, nddli, nnzi,
167 . iadi, jdii, itok, diag_i, lt_i )
168C
169c CALL SPMD_MUMPS_FRONT(
170c . ITK, RTK, NKFRONT, NKFLOC, NKLOC,
171c . NDDLG, IMPRINT )
172C
173 nkfloc = 0
174 nzloc=nkloc+nkfloc
175 IF (imumpsd==1) THEN
176 ALLOCATE(mumps_par%A_LOC(nzloc),
177 . mumps_par%IRN_LOC(nzloc),
178 . mumps_par%JCN_LOC(nzloc))
179 IF (ispmd==0) THEN
180 ALLOCATE(mumps_par%RHS(nddlg))
181 ELSE
182 ALLOCATE(mumps_par%RHS(0))
183 ENDIF
184 mumps_par%N=nddlg
185 mumps_par%NZ_LOC=nzloc
186C
187 DO i=1,nzloc
188 mumps_par%IRN_LOC(i)=itk(1,i)
189 mumps_par%JCN_LOC(i)=itk(2,i)
190 mumps_par%A_LOC(i)=rtk(i)
191 ENDDO
192 ELSEIF (imumpsd==2) THEN
193 CALL spmd_mumps_count(nzloc, nzp, nnz)
194C
195 IF (ispmd==0) THEN
196 ALLOCATE(mumps_par%A(nnz),
197 . mumps_par%IRN(nnz),
198 . mumps_par%JCN(nnz),
199 . mumps_par%RHS(nddlg))
200 mumps_par%N=nddlg
201 mumps_par%NZ=nnz
202 ELSE
203 ALLOCATE(mumps_par%A(0),
204 . mumps_par%IRN(0),
205 . mumps_par%JCN(0),
206 . mumps_par%RHS(0))
207 ENDIF
208C
209 CALL spmd_mumps_gath(
210 . itk, rtk, nzloc, mumps_par%A, mumps_par%IRN,
211 . mumps_par%JCN, nzp)
212C
213 ENDIF
214 DEALLOCATE(itk, rtk)
215 ELSE
216 mumps_par%ICNTL(18)=0
217C
218 DO i=1,nddl
219 cddlp(i)=i
220 ENDDO
221 nnzt = nnzk
222 nnzk = nnzk + nddli + nnzi
223C
224 ALLOCATE(mumps_par%A(nddl+nnzk),
225 . mumps_par%IRN(nddl+nnzk),
226 . mumps_par%JCN(nddl+nnzk),
227 . mumps_par%RHS(nddl))
228C
229 nnz=0
230 DO i=1,nddli
231 j=itok(i)
232 nnz=nnz+1
233 mumps_par%IRN(nnz)=j
234 mumps_par%JCN(nnz)=j
235 mumps_par%A(nnz)=diag_i(i)
236 DO n=iadi(i),iadi(i+1)-1
237 nnz=nnz+1
238 mumps_par%IRN(nnz)=j
239 mumps_par%JCN(nnz)=itok(jdii(n))
240 mumps_par%A(nnz)=lt_i(n)
241 ENDDO
242 ENDDO
243 DO i=1,nddl
244 nnz=nnz+1
245 mumps_par%IRN(nnz)=i
246 mumps_par%JCN(nnz)=i
247 mumps_par%A(nnz)=diag_k(i)
248 DO j=iadk(i),iadk(i+1)-1
249 nnz=nnz+1
250 mumps_par%IRN(nnz)=i
251 mumps_par%JCN(nnz)=jdik(j)
252 mumps_par%A(nnz)=lt_k(j)
253 ENDDO
254 ENDDO
255C
256 IF (imprint/=0) THEN
257 WRITE(istdo,*)
258 WRITE(istdo,'(A21,I10,A8,I10)')
259 .' MUMPS DIM : NNZ =',nnzk+nddl,' NNZFR =',0
260 ENDIF
261C
262 mumps_par%N=nddl
263 mumps_par%NZ=nnzk+nddl
264 nnzk = nnzt
265 ENDIF
266C
267c WRITE(IOUT,*) "NCYCLE,IT=",NCYCLE,IT
268 IF (prstifmat == 1 .AND. (iline==1 .OR. (prstifmat_nc == ncycle
269 . .AND. prstifmat_it == it))) THEN
270 IF (ispmd == 0) THEN
271 WRITE(iout,1000)
272 WRITE(istdo,1000)
273 WRITE(iout,*)
274 WRITE(istdo,*)
275 ENDIF
276 CALL print_stiff_mat(mumps_par, nddl, nodglob, iddl, ndof,
277 . cddlp, inloc, ikc, nddlg, nddlp)
278
279
280 ENDIF
2811000 FORMAT(5x,'--STIFFNESS MATRIX IS PRINTED--')
282 RETURN
283#endif
284 END
285!||====================================================================
286!|| print_stiff_mat ../engine/source/implicit/imp_mumps.F
287!||--- called by ------------------------------------------------------
288!|| imp_mumps1 ../engine/source/implicit/imp_mumps.F
289!||--- calls -----------------------------------------------------
290!|| spmd_int_allreduce_max ../engine/source/mpi/implicit/imp_spmd.F
291!||--- uses -----------------------------------------------------
292!|| spmd_comm_world_mod ../engine/source/mpi/spmd_comm_world.F90
293!||====================================================================
294 SUBROUTINE print_stiff_mat(MUMPS_PAR, NDDL, NODGLOB, IDDL, NDOF,
295 . CDDLP, INLOC, IKC, NDDLG, NDDLP)
296C-----------------------------------------------
297C I m p l i c i t T y p e s
298C-----------------------------------------------
299 USE spmd_comm_world_mod, ONLY : spmd_comm_world
300#include "implicit_f.inc"
301C-----------------------------------------------
302C C o m m o n B l o c k s
303C-----------------------------------------------
304#if defined(MUMPS5)
305#include "dmumps_struc.h"
306#endif
307#include "task_c.inc"
308#include "com01_c.inc"
309#include "com04_c.inc"
310#include "impl2_c.inc"
311C-----------------------------------------------
312C M e s s a g e P a s s i n g
313C-----------------------------------------------
314#include "spmd.inc"
315C-----------------------------------------------
316C D u m m y A r g u m e n t s
317C-----------------------------------------------
318#ifdef MUMPS5
319 TYPE(dmumps_struc) mumps_par
320#else
321 ! Fake declaration as DMUMPS_STRUC is shipped with MUMPS
322 INTEGER mumps_par
323#endif
324 INTEGER nddl, nodglob(*), IDDL(*), NDOF(*), CDDLP(*), INLOC(*), IKC(*), NDDLG, NDDLP(*)
325#ifdef MUMPS5
326C-----------------------------------------------
327C L o c a l V a r i a b l e s
328C-----------------------------------------------
329 INTEGER I,J,II,JJ,NROWS,NCOLS,NENTRIES,IItmp,JJtmp,
330 . tmpirows,tmpicols,length,l,ind,count_dup,
331 . nkc,tddl(2,mumps_par%N),n,id,nd,ierror,offst
332 my_real
333 . tmpk,sumk,tmpproc
334 INTEGER, DIMENSION(:), ALLOCATABLE :: IROWS,ICOLS,NENTRIEStmp
335 my_real, DIMENSION(:), ALLOCATABLE :: k
336 LOGICAL SWITCH
337 CHARACTER FILNAME*100,FILNAME2*100,CSPMD
338C-----------------------------------------------
339C Define offset for UNIT file
340 offst = 100
341C Automatic write of stiffness coefficients in MatrixMarket format by MUMPS: MUMPS_PAR%WRITE_PROBLEM = 'string'
342C MUMPS_PAR%WRITE_PROBLEM ="./stiffness_matrix_MUMPS"
343 WRITE(cspmd,'(I1)') ispmd
344 filname = 'local_stiffness_matrix_domain'//cspmd
345 OPEN(unit=offst+ispmd,file=filname(1:30),access="SEQUENTIAL",
346 . action="WRITE",status="UNKNOWN")
347C Manual write of stiffness coefficients in MatrixMarket format
348C TDDL: local DOF (at MPI domain level) to global node and direction (DX,DY,DZ,RX,RY,RZ)
349 nkc = 0
350 tddl(1,:) = 0
351 tddl(2,:) = 0
352 DO n=1,numnod
353 i=inloc(n)
354 DO j=1,ndof(i)
355 nd=iddl(i)+j
356 id=nd-nkc
357 IF (ikc(nd)<1) THEN
358 tddl(1,cddlp(id))=nodglob(i)
359 tddl(2,cddlp(id))=j
360 ELSE
361 nkc=nkc+1
362 ENDIF
363 ENDDO
364 ENDDO
365C Communications between processes for TDDL(1:2,1:MUMPS_PAR%N)
366 IF (nspmd > 1) THEN
367 CALL spmd_int_allreduce_max(tddl(1,:),tddl(1,:),
368 . mumps_par%N)
369 CALL spmd_int_allreduce_max(tddl(2,:),tddl(2,:),
370 . mumps_par%N)
371 END IF
372C All processes: write local stiffness coefficients (at MPI domain level) in MatrixMarket format
373 IF (nspmd == 1) THEN
374 WRITE(offst+ispmd,1002) mumps_par%N,mumps_par%N,mumps_par%NZ
375 DO i=1,mumps_par%NZ
376 iitmp = mumps_par%IRN(i)
377 jjtmp = mumps_par%JCN(i)
378 ii = 6*(tddl(1,iitmp)-1)+tddl(2,iitmp)
379 jj = 6*(tddl(1,jjtmp)-1)+tddl(2,jjtmp)
380 IF (jj > ii) THEN
381 iitmp = ii
382 ii = jj
383 jj = iitmp
384 ENDIF
385 WRITE(offst+ispmd,1003) ii,jj,mumps_par%A(i)
386 ENDDO
387 ELSE
388 WRITE(offst+ispmd,1002) mumps_par%N,mumps_par%N,mumps_par%NZ_LOC
389 DO i=1,mumps_par%NZ_LOC
390 iitmp = mumps_par%IRN_LOC(i)
391 jjtmp = mumps_par%JCN_LOC(i)
392 ii = 6*(tddl(1,iitmp)-1)+tddl(2,iitmp)
393 jj = 6*(tddl(1,jjtmp)-1)+tddl(2,jjtmp)
394 IF (jj > ii) THEN
395 iitmp = ii
396 ii = jj
397 jj = iitmp
398 ENDIF
399 WRITE(offst+ispmd,1003) ii,jj,mumps_par%A_LOC(i)
400 ENDDO
401 ENDIF
402 CLOSE(unit=offst+ispmd)
403#ifdef MPI
404 CALL mpi_barrier(spmd_comm_world,ierror)
405#endif
406C Process 0: read stiffness coefficients from all processes, sort (bubble) and write stiffness coefficients in MatrixMarket format
407 IF (ispmd == 0) THEN
408 ALLOCATE(nentriestmp(nspmd))
409 WRITE(cspmd,'(I1)') nspmd
410 nentries = 0
411 filname2 = 'stiffness_matrix_'//cspmd//'_SPMD'
412 OPEN(unit=offst+nspmd,file=filname2(1:23),access="SEQUENTIAL",
413 . action="WRITE",status="UNKNOWN")
414 DO l = 0,nspmd-1
415 WRITE(cspmd,'(I1)') l
416 filname = 'local_stiffness_matrix_domain'//cspmd
417 OPEN(unit=offst+l,file=filname(1:30),access="SEQUENTIAL",
418 . action="READ",status="UNKNOWN")
419 READ(unit=offst+l,fmt=*) nrows,ncols,nentriestmp(l+1)
420 nentries = nentries + nentriestmp(l+1)
421 ENDDO
422 ALLOCATE(irows(nentries))
423 ALLOCATE(icols(nentries))
424 ALLOCATE(k(nentries))
425 ind = 0
426 sumk = zero
427 DO l = 0,nspmd-1
428 DO i = 1,nentriestmp(l+1)
429 ind = ind + 1
430 READ(unit=offst+l,fmt=*) irows(ind),icols(ind),k(ind)
431 ENDDO
432 ENDDO
433C Bubble sort in ascending order of ICOLS and then IROWS
434 i = nentries
435 switch = .true.
436 DO WHILE ((i>0) .AND. (switch))
437 switch = .false.
438 DO j = 1,i-1
439 IF (icols(j) > icols(j+1) .OR. (icols(j) == icols(j+1)
440 . .AND. irows(j) > irows(j+1))) THEN
441 tmpirows = irows(j)
442 irows(j) = irows(j+1)
443 irows(j+1) = tmpirows
444 tmpicols = icols(j)
445 icols(j) = icols(j+1)
446 icols(j+1) = tmpicols
447 tmpk = k(j)
448 k(j) = k(j+1)
449 k(j+1) = tmpk
450 switch = .true.
451 ENDIF
452 ENDDO
453 i = i - 1
454 ENDDO
455C Write stiffness coefficients in only one file (suppress duplications)
456 ind = 1
457 DO WHILE (ind <= nentries)
458 tmpk = k(ind)
459 DO WHILE (ind <= nentries .AND. irows(ind)==irows(ind+1)
460 . .AND. icols(ind)==icols(ind+1))
461 ind = ind + 1
462 tmpk = tmpk + k(ind)
463 ENDDO
464 IF (abs(tmpk)>=prstifmat_tol) THEN
465 WRITE(offst+nspmd,1003) irows(ind),icols(ind),tmpk
466 sumk = sumk + abs(tmpk)
467 ENDIF
468 ind = ind + 1
469 ENDDO
470 WRITE(offst+nspmd,1001) sumk
471 DO l = 0,nspmd
472 CLOSE(unit=offst+l)
473 ENDDO
474 DEALLOCATE(nentriestmp)
475 DEALLOCATE(irows)
476 DEALLOCATE(icols)
477 DEALLOCATE(k)
478 ENDIF
4791000 FORMAT(i10,i10,i10,i10,e10.2)
4801001 FORMAT('Sum ABS(K_ij) = ',e10.2)
4811002 FORMAT(i10,i10,i10)
4821003 FORMAT(i10,i10,e10.2)
483C
484#endif
485 RETURN
486 END
487!||====================================================================
488!|| imp_mumps2 ../engine/source/implicit/imp_mumps.F
489!||--- called by ------------------------------------------------------
490!|| lin_solvp2 ../engine/source/implicit/lin_solv.F
491!||--- calls -----------------------------------------------------
492!|| spmd_mumps_exec ../engine/source/mpi/implicit/imp_spmd.F
493!|| spmd_mumps_rhs ../engine/source/mpi/implicit/imp_spmd.F
494!||--- uses -----------------------------------------------------
495!|| spmd_comm_world_mod ../engine/source/mpi/spmd_comm_world.F90
496!||====================================================================
497 SUBROUTINE imp_mumps2(MUMPS_PAR, CDDLP, F, X, NDDL)
498C-----------------------------------------------
499C I m p l i c i t T y p e s
500C-----------------------------------------------
501 USE spmd_comm_world_mod, ONLY : spmd_comm_world
502#include "implicit_f.inc"
503C-----------------------------------------------
504C C o m m o n B l o c k s
505C-----------------------------------------------
506#if defined(MUMPS5)
507#include "dmumps_struc.h"
508#endif
509#include "impl1_c.inc"
510#include "filescount_c.inc"
511C-----------------------------------------------
512C D u m m y A r g u m e n t s
513C-----------------------------------------------
514 INTEGER CDDLP(*), NDDL
515 my_real F(*), X(*)
516#ifdef MUMPS5
517 TYPE(DMUMPS_STRUC) MUMPS_PAR
518#else
519 ! Fake declaration as DMUMPS_STRUC is shipped with MUMPS
520 INTEGER MUMPS_PAR
521#endif
522
523#ifdef MUMPS5
524C----------------------------------------------
525C L o c a l V a r i a b l e s
526C-----------------------------------------------
527 INTEGER I, NDDLG,LENT
528C
529 IF (MUMPS_PAR%N<=0) return
530 IF (idsc==1) CALL spmd_mumps_exec(mumps_par, 1)
531C---------For licence
532 IF (idsc==1) isolv_d = 1
533C
534 nddlg=mumps_par%N
535 CALL spmd_mumps_rhs(f, cddlp, mumps_par%RHS, nddl, 1,
536 . nddlg)
537C
538 CALL spmd_mumps_exec(mumps_par, 2)
539C
540 CALL spmd_mumps_rhs(x, cddlp, mumps_par%RHS, nddl, 2,
541 . nddlg)
542C----FLAG for MUMPS: IF (IMPL_S>0.AND.ISOLV==3)
543 lent = 1024*mumps_par%INFO(16)
544 mumpsfilesize = max(mumpsfilesize,lent)
545C
546 RETURN
547#endif
548 END
549
550!||====================================================================
551!|| mumps_set ../engine/source/implicit/imp_mumps.F
552!||--- calls -----------------------------------------------------
553!|| spmd_ifri ../engine/source/mpi/implicit/imp_spmd.F
554!||--- uses -----------------------------------------------------
555!|| imp_intm ../engine/share/modules/imp_intm.F
556!|| spmd_comm_world_mod ../engine/source/mpi/spmd_comm_world.F90
557!||====================================================================
558 SUBROUTINE mumps_set(IADK , JDIK , DIAG_K , LT_K , CDDLP,
559 . NKLOC, NKFRONT , ITK , RTK , IDDL ,
560 . INLOC, IAD_ELEM, FR_ELEM, NDOF , IKC ,
561 . NDDL , NNZK , IACTI , NDDLI , NNZI ,
562 . IADI , JDII , ITOK , DIAG_I, LT_I )
563C-----------------------------------------------
564C M o d u l e s
565C-----------------------------------------------
566 USE imp_intm
567C-----------------------------------------------
568C I m p l i c i t T y p e s
569C-----------------------------------------------
570 USE spmd_comm_world_mod, ONLY : spmd_comm_world
571#include "implicit_f.inc"
572C-----------------------------------------------
573C C o m m o n B l o c k s
574C-----------------------------------------------
575#include "com01_c.inc"
576#include "com04_c.inc"
577#include "com08_c.inc"
578#include "task_c.inc"
579C-----------------------------------------------
580C D u m m y A r g u m e n t s
581C-----------------------------------------------
582 INTEGER IADK(*), JDIK(*), CDDLP(*), NKLOC, NKFRONT, ITK(2,*),
583 . IDDL(*), INLOC(*), IAD_ELEM(2,*), FR_ELEM(*), NDOF(*),
584 . IKC(*), NDDL, NNZK, IACTI(*), NDDLI, NNZI, IADI(*),
585 . JDII(*), ITOK(*)
586 my_real
587 . DIAG_K(*), LT_K(*), RTK(*), DIAG_I(*), LT_I(*)
588#ifdef MUMPS5
589C----------------------------------------------
590C L o c a l V a r i a b l e s
591C-----------------------------------------------
592 INTEGER I, IDDL_FRONT(NSPMD+1,NDDL), NKC, N, TNKC(NUMNOD),
593 . J, ND, NOD, INOD, KK, K, ID, NN, ILOC, JJ,
594 . ITAG(2,NSPMD), INDEX, II, IDIAG(NDDL), IADL(NDDL),
595 . IADLFRONT(NDDL), IFOUND, CDDLP_REM(NDDL_SI),
596 . iddl_rem(nddl_si)
597 INTEGER, DIMENSION(:,:), ALLOCATABLE :: ITKFRONT
598 my_real, DIMENSION(:), ALLOCATABLE :: RTKFRONT
599C
600 ALLOCATE(ITKFRONT(2,NDDL+NNZK+NNZI+NZ_SI),
601 . RTKFRONT(NDDL+NNZK+NNZI+NZ_SI))
602C
603 DO i=1,nddl
604 iddl_front(1,i)=1
605 iddl_front(2,i)=ispmd+1
606 ENDDO
607C
608 nkc=0
609 DO n=1,numnod
610 i=inloc(n)
611 tnkc(i)=nkc
612 DO j=1,ndof(i)
613 nd=iddl(i)+j
614 id=nd-nkc
615 IF (ikc(nd)>=1) nkc=nkc+1
616 ENDDO
617 ENDDO
618C
619 DO i=1,nspmd
620C IF (I==ISPMD+1) CYCLE
621C
622 DO j=iad_elem(1,i),iad_elem(1,i+1)-1
623 inod=fr_elem(j)
624 nkc=tnkc(inod)
625 DO k=1,ndof(inod)
626 nd=iddl(inod)+k
627 id=nd-nkc
628 IF (ikc(nd)<1) THEN
629 iddl_front(1,id)=iddl_front(1,id)+1
630 nn=iddl_front(1,id)
631 iddl_front(nn+1,id)=i
632 ELSE
633 nkc=nkc+1
634 ENDIF
635 ENDDO
636 ENDDO
637 ENDDO
638C
639 nkloc=0
640 nkfront=0
641 DO i=1,nddl
642 IF (iacti(i)==0) cycle
643 ii=iacti(i)
644 iadl(ii)=nkloc+1
645 iadlfront(ii)=nkfront+1
646 IF (iddl_front(1,ii)==1) THEN
647 nkloc=nkloc+1
648 idiag(ii)=nkloc
649 itk(1,nkloc)=cddlp(ii)
650 itk(2,nkloc)=cddlp(ii)
651 rtk(nkloc)=diag_k(i)
652 ELSE
653 nkfront=nkfront+1
654 idiag(ii)=nkfront
655 itkfront(1,nkfront)=cddlp(ii)
656 itkfront(2,nkfront)=cddlp(ii)
657 rtkfront(nkfront)=diag_k(i)
658 ENDIF
659C
660 DO j=iadk(i),iadk(i+1)-1
661 iloc=1
662 jj=iacti(jdik(j))
663 IF (jj==0) cycle
664 IF (iddl_front(1,ii)==1.OR.iddl_front(1,jj)==1) THEN
665 iloc=0
666 ELSE
667 DO k=1,nspmd
668 itag(1,k)=0
669 itag(2,k)=0
670 ENDDO
671 DO k=1,iddl_front(1,ii)
672 kk=iddl_front(1+k,ii)
673 itag(1,kk)=1
674 ENDDO
675 DO k=1,iddl_front(1,jj)
676 kk=iddl_front(1+k,jj)
677 itag(2,kk)=1
678 ENDDO
679 index=0
680 DO k=1,nspmd
681 index=index+itag(1,k)*itag(2,k)
682 ENDDO
683 IF (index==1) iloc=0
684 ENDIF
685C
686 IF (iloc==0) THEN
687 nkloc=nkloc+1
688 itk(1,nkloc)=cddlp(ii)
689 itk(2,nkloc)=cddlp(jj)
690 rtk(nkloc)=lt_k(j)
691 ELSEIF (iloc==1) THEN
692 nkfront=nkfront+1
693 itkfront(1,nkfront)=cddlp(ii)
694 itkfront(2,nkfront)=cddlp(jj)
695 rtkfront(nkfront)=lt_k(j)
696 ENDIF
697 ENDDO
698 ENDDO
699 IF (nddli>0) THEN
700C Matrice de rigidite d'interface
701 DO n=1,nddli
702 i=itok(n)
703 IF (iacti(i)==0) cycle
704 ii=iacti(i)
705 IF (iddl_front(1,ii)==1) THEN
706 j=idiag(ii)
707 rtk(j)=rtk(j)+diag_i(n)
708 ELSE
709 j=idiag(ii)
710 rtkfront(j)=rtkfront(j)+diag_i(n)
711 ENDIF
712C
713 DO j=iadi(n),iadi(n+1)-1
714 iloc=1
715 jj=itok(jdii(j))
716 jj=iacti(jj)
717 IF (jj==0) cycle
718 IF (iddl_front(1,ii)==1.OR.iddl_front(1,jj)==1) THEN
719 iloc=0
720 ELSE
721 DO k=1,nspmd
722 itag(1,k)=0
723 itag(2,k)=0
724 ENDDO
725 DO k=1,iddl_front(1,ii)
726 kk=iddl_front(1+k,ii)
727 itag(1,kk)=1
728 ENDDO
729 DO k=1,iddl_front(1,jj)
730 kk=iddl_front(1+k,jj)
731 itag(2,kk)=1
732 ENDDO
733 index=0
734 DO k=1,nspmd
735 index=index+itag(1,k)*itag(2,k)
736 ENDDO
737 IF (index==1) iloc=0
738 ENDIF
739C
740 IF (iloc==0) THEN
741 ifound=0
742 k=iadl(ii)
743 DO WHILE (ifound==0.AND.k<=iadl(ii+1)-1)
744 IF (cddlp(ii)==itk(1,k)
745 . .AND.cddlp(jj)==itk(2,k)) ifound=k
746 k=k+1
747 ENDDO
748 IF (ifound/=0) THEN
749 rtk(ifound)=rtk(ifound)+lt_i(j)
750 ELSE
751 nkloc=nkloc+1
752 itk(1,nkloc)=cddlp(ii)
753 itk(2,nkloc)=cddlp(jj)
754 rtk(nkloc)=lt_i(j)
755 ENDIF
756 ELSEIF (iloc==1) THEN
757 ifound=0
758 k=iadlfront(ii)
759 DO WHILE (ifound==0.AND.k<=iadlfront(ii+1)-1)
760 IF (cddlp(ii)==itkfront(1,k)
761 . .AND.cddlp(jj)==itkfront(2,k)) ifound=k
762 k=k+1
763 ENDDO
764 IF (ifound/=0) THEN
765 rtkfront(ifound)=rtkfront(ifound)+lt_i(j)
766 ELSE
767 nkfront=nkfront+1
768 itkfront(1,nkfront)=cddlp(ii)
769 itkfront(2,nkfront)=cddlp(jj)
770 rtkfront(nkfront)=lt_i(j)
771 ENDIF
772 ENDIF
773 ENDDO
774 ENDDO
775 ENDIF
776C Complement to the interface stiffness matrix for remote seconds
777 DO i=1,nddl_sl
778 ii=iddl_sl(i)
779 DO j=iad_ss(i),iad_ss(i+1)-1
780 iloc=1
781 jj=jdi_sl(j)
782 IF (iddl_front(1,ii)==1.OR.iddl_front(1,jj)==1) THEN
783 iloc=0
784 ELSE
785 DO k=1,nspmd
786 itag(1,k)=0
787 itag(2,k)=0
788 ENDDO
789 DO k=1,iddl_front(1,ii)
790 kk=iddl_front(1+k,ii)
791 itag(1,kk)=1
792 ENDDO
793 DO k=1,iddl_front(1,jj)
794 kk=iddl_front(1+k,jj)
795 itag(2,kk)=1
796 ENDDO
797 index=0
798 DO k=1,nspmd
799 index=index+itag(1,k)*itag(2,k)
800 ENDDO
801 IF (index==1) iloc=0
802 ENDIF
803 IF (iloc==0) THEN
804 nkloc=nkloc+1
805 itk(1,nkloc)=cddlp(ii)
806 itk(2,nkloc)=cddlp(jj)
807 rtk(nkloc)=lt_sl(j)
808 ELSEIF (iloc==1) THEN
809 nkfront=nkfront+1
810 itkfront(1,nkfront)=cddlp(ii)
811 itkfront(2,nkfront)=cddlp(jj)
812 rtkfront(nkfront)=lt_sl(j)
813 ENDIF
814 ENDDO
815 ENDDO
816C----- It lacks diag_sl ---------
817 DO n=1,nddl_sl
818 i=iddl_sl(n)
819 IF (iacti(i)==0) cycle
820 ii=iacti(i)
821 j=idiag(ii)
822 IF (iddl_front(1,ii)==1) THEN
823 rtk(j)=rtk(j)+diag_sl(n)
824 ELSE
825 rtkfront(j)=rtkfront(j)+diag_sl(n)
826 ENDIF
827 ENDDO
828 IF ((nddl_si+nddl_sl)>0) THEN
829 CALL spmd_ifri(cddlp, cddlp_rem)
830 DO i=1,nddl
831 iadl(i) = iddl_front(1,i)
832 ENDDO
833 CALL spmd_ifri(iadl, iddl_rem)
834 ENDIF
835 DO i=1,nddl_si
836 DO j=iad_si(i),iad_si(i+1)-1
837 jj=jdi_si(j)
838 IF (iddl_front(1,jj)==1.OR.iddl_rem(i)==1) THEN
839 nkloc=nkloc+1
840 itk(1,nkloc)=cddlp_rem(i)
841 itk(2,nkloc)=cddlp(jj)
842 rtk(nkloc)=lt_si(j)
843 ELSE
844 nkfront=nkfront+1
845 itkfront(1,nkfront)=cddlp_rem(i)
846 itkfront(2,nkfront)=cddlp(jj)
847 rtkfront(nkfront)=lt_si(j)
848 ENDIF
849 ENDDO
850 ENDDO
851C
852 DO i=1,nkfront
853 itk(1,nkloc+i)=itkfront(1,i)
854 itk(2,nkloc+i)=itkfront(2,i)
855 rtk(nkloc+i)=rtkfront(i)
856 ENDDO
857C
858 DEALLOCATE(itkfront, rtkfront)
859C
860 RETURN
861#endif
862 END
863
864
865!||====================================================================
866!|| mumps_set2 ../engine/source/implicit/imp_mumps.F
867!||--- called by ------------------------------------------------------
868!|| imp_buck ../engine/source/implicit/imp_buck.F
869!|| imp_mumps1 ../engine/source/implicit/imp_mumps.F
870!||--- calls -----------------------------------------------------
871!|| spmd_ifri ../engine/source/mpi/implicit/imp_spmd.F
872!||--- uses -----------------------------------------------------
873!|| imp_intm ../engine/share/modules/imp_intm.F
874!|| spmd_comm_world_mod ../engine/source/mpi/spmd_comm_world.F90
875!||====================================================================
876 SUBROUTINE mumps_set2(IADK , JDIK , DIAG_K , LT_K , CDDLP,
877 . NKLOC, NKFRONT , ITK , RTK , IDDL ,
878 . INLOC, IAD_ELEM, FR_ELEM, NDOF , IKC ,
879 . NDDL , NNZK , IACTI , NDDLI , NNZI ,
880 . IADI , JDII , ITOK , DIAG_I, LT_I )
881C-----------------------------------------------
882C M o d u l e s
883C-----------------------------------------------
884 USE imp_intm
885C-----------------------------------------------
886C I m p l i c i t T y p e s
887C-----------------------------------------------
888 USE spmd_comm_world_mod, ONLY : spmd_comm_world
889#include "implicit_f.inc"
890C-----------------------------------------------
891C C o m m o n B l o c k s
892C-----------------------------------------------
893#include "com01_c.inc"
894#include "com04_c.inc"
895C-----------------------------------------------
896C D u m m y A r g u m e n t s
897C-----------------------------------------------
898 INTEGER IADK(*), JDIK(*), CDDLP(*), NKLOC, NKFRONT, ITK(2,*),
899 . IDDL(*), INLOC(*), IAD_ELEM(2,*), FR_ELEM(*), NDOF(*),
900 . IKC(*), NDDL, NNZK, IACTI(*), NDDLI, NNZI, IADI(*),
901 . JDII(*), ITOK(*)
902 my_real
903 . diag_k(*), lt_k(*), rtk(*), diag_i(*), lt_i(*)
904#ifdef MUMPS5
905C----------------------------------------------
906C L o c a l V a r i a b l e s
907C-----------------------------------------------
908 INTEGER I, NKC, N,
909 . J, ND, NOD, INOD, KK, K, ID, NN, ILOC, JJ,
910 . ITAG(2,NSPMD), INDEX, II, IDIAG(NDDL), IADL(NDDL+1),
911 . IFOUND, CDDLP_REM(NDDL_SI),
912 . iddl_rem(nddl_si)
913 nkc=0
914 DO n=1,numnod
915 i=inloc(n)
916 DO j=1,ndof(i)
917 nd=iddl(i)+j
918 id=nd-nkc
919 IF (ikc(nd)>=1) nkc=nkc+1
920 ENDDO
921 ENDDO
922
923 nkloc=0
924 DO i=1,nddl
925 IF (iacti(i)==0) cycle
926 ii=iacti(i)
927 iadl(ii)=nkloc+1
928 nkloc=nkloc+1
929 idiag(ii)=nkloc
930 itk(1,nkloc)=cddlp(ii)
931 itk(2,nkloc)=cddlp(ii)
932 rtk(nkloc)=diag_k(i)
933 DO j=iadk(i),iadk(i+1)-1
934 iloc=1
935 jj=iacti(jdik(j))
936 IF (jj==0) cycle
937 nkloc=nkloc+1
938 itk(1,nkloc)=cddlp(ii)
939 itk(2,nkloc)=cddlp(jj)
940 rtk(nkloc)=lt_k(j)
941 ENDDO
942 ENDDO
943 iadl(nddl+1) = nkloc+1
944 IF (nddli>0) THEN
945C Matrice de rigidite d'interface
946 DO n=1,nddli
947 i=itok(n)
948 IF (iacti(i)==0) cycle
949 ii=iacti(i)
950 j=idiag(ii)
951 rtk(j)=rtk(j)+diag_i(n)
952 DO j=iadi(n),iadi(n+1)-1
953 iloc=1
954 jj=itok(jdii(j))
955 jj=iacti(jj)
956 IF (jj==0) cycle
957 ifound=0
958 k=iadl(ii)
959 DO WHILE (ifound==0.AND.k<=iadl(ii+1)-1)
960 IF (cddlp(ii)==itk(1,k)
961 . .AND.cddlp(jj)==itk(2,k)) ifound=k
962 k=k+1
963 ENDDO
964 IF (ifound/=0) THEN
965 rtk(ifound)=rtk(ifound)+lt_i(j)
966 ELSE
967 nkloc=nkloc+1
968 itk(1,nkloc)=cddlp(ii)
969 itk(2,nkloc)=cddlp(jj)
970 rtk(nkloc)=lt_i(j)
971 ENDIF
972 ENDDO
973 ENDDO
974 ENDIF
975C Complement to the interface stiffness matrix for remote seconds
976 DO i=1,nddl_sl
977 ii=iddl_sl(i)
978 DO j=iad_ss(i),iad_ss(i+1)-1
979 iloc=1
980 jj=jdi_sl(j)
981 nkloc=nkloc+1
982 itk(1,nkloc)=cddlp(ii)
983 itk(2,nkloc)=cddlp(jj)
984 rtk(nkloc)=lt_sl(j)
985 ENDDO
986 ENDDO
987C----- il manque DIAG_SL--------
988 DO n=1,nddl_sl
989 i=iddl_sl(n)
990 IF (iacti(i)==0) cycle
991 ii=iacti(i)
992 j=idiag(ii)
993 rtk(j)=rtk(j)+diag_sl(n)
994 ENDDO
995 IF ((nddl_si+nddl_sl)>0) THEN
996 CALL spmd_ifri(cddlp, cddlp_rem)
997 ENDIF
998 DO i=1,nddl_si
999 DO j=iad_si(i),iad_si(i+1)-1
1000 jj=jdi_si(j)
1001 nkloc=nkloc+1
1002 itk(1,nkloc)=cddlp_rem(i)
1003 itk(2,nkloc)=cddlp(jj)
1004 rtk(nkloc)=lt_si(j)
1005 ENDDO
1006 ENDDO
1007C
1008 RETURN
1009#endif
1010 END
my_real
#define my_real
Definition cppsort.cpp:32
mumps_set
subroutine mumps_set(iadk, jdik, diag_k, lt_k, cddlp, nkloc, nkfront, itk, rtk, iddl, inloc, iad_elem, fr_elem, ndof, ikc, nddl, nnzk, iacti, nddli, nnzi, iadi, jdii, itok, diag_i, lt_i)
Definition imp_mumps.F:563
print_stiff_mat
subroutine print_stiff_mat(mumps_par, nddl, nodglob, iddl, ndof, cddlp, inloc, ikc, nddlg, nddlp)
Definition imp_mumps.F:296
mumps_set2
subroutine mumps_set2(iadk, jdik, diag_k, lt_k, cddlp, nkloc, nkfront, itk, rtk, iddl, inloc, iad_elem, fr_elem, ndof, ikc, nddl, nnzk, iacti, nddli, nnzi, iadi, jdii, itok, diag_i, lt_i)
Definition imp_mumps.F:881
imp_mumps1
subroutine imp_mumps1(nddl0, nnzk0, nddl, nnzk, nnmax, nodglob, iddl, ndof, inloc, ikc, iadk, jdik, diag_k, lt_k, iad_elem, fr_elem, mumps_par, cddlp, iadi, jdii, itok, diag_i, lt_i, nddli, nnzi, imprint, it)
Definition imp_mumps.F:49
imp_mumps2
subroutine imp_mumps2(mumps_par, cddlp, f, x, nddl)
Definition imp_mumps.F:498
spmd_ifri
subroutine spmd_ifri(ig, il)
Definition imp_spmd.F:4262
spmd_mumps_rhs
subroutine spmd_mumps_rhs(v, cddlp, rhs, nddl, isens, nddlg)
Definition imp_spmd.F:611
spmd_mumps_ini
subroutine spmd_mumps_ini(mumps_par, sym)
Definition imp_spmd.F:498
spmd_mumps_gath
subroutine spmd_mumps_gath(itk, rtk, nzloc, a, irn, jcn, nzp)
Definition imp_spmd.F:408
spmd_mumps_deal
subroutine spmd_mumps_deal(mumps_par)
Definition imp_spmd.F:558
spmd_mumps_count
subroutine spmd_mumps_count(nzloc, nzp, nnz)
Definition imp_spmd.F:350
spmd_int_allreduce_max
subroutine spmd_int_allreduce_max(sendbuf, recvbuf, count)
Definition imp_spmd.F:5131
spmd_inf_g
subroutine spmd_inf_g(nddl0, nzzk0, nddl, nzzk, nnmax, nddl0p, nzzk0p, nddlp, nzzkp, nnmaxp)
Definition imp_spmd.F:1514
spmd_mumps_exec
subroutine spmd_mumps_exec(mumps_par, itask)
Definition imp_spmd.F:724
spmd_cddl
subroutine spmd_cddl(nddl, nodglob, iddl, ndof, cddlp, inloc, ikc, nddlg, nddlp)
Definition imp_spmd.F:3146
max
#define max(a, b)
Definition macros.h:21
mpi_barrier
subroutine mpi_barrier(comm, ierr)
Definition mpi.f:188
id
initmumps id
Definition multiplerhs_example.m:4
imp_intm
Definition imp_intm.F:118
imp_intm::jdi_si
integer, dimension(:), allocatable jdi_si
Definition imp_intm.F:174
imp_intm::iddl_sl
integer, dimension(:), allocatable iddl_sl
Definition imp_intm.F:178
imp_intm::iad_ss
integer, dimension(:), allocatable iad_ss
Definition imp_intm.F:175
imp_intm::nz_sl
integer nz_sl
Definition imp_intm.F:173
imp_intm::iad_si
integer, dimension(:), allocatable iad_si
Definition imp_intm.F:174
imp_intm::jdi_sl
integer, dimension(:), allocatable jdi_sl
Definition imp_intm.F:175
imp_intm::nddl_sl
integer nddl_sl
Definition imp_intm.F:173
imp_intm::nz_si
integer nz_si
Definition imp_intm.F:173
imp_kbcs
Definition impbufdef_mod.F:471
message_mod
Definition message_mod.F:1257
tmpenvf
void tmpenvf(char *tmpdir, int *tmplen)
Definition tmpenv_c.c:149