FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_solver_direct_MUMPS.f90
Go to the documentation of this file.
1!-------------------------------------------------------------------------------
2! Copyright (c) 2019 FrontISTR Commons
3! This software is released under the MIT License, see LICENSE.txt
4!-------------------------------------------------------------------------------
7 use hecmw_util
13
14 private
16
17 logical, save :: INITIALIZED = .false.
18 type (sparse_matrix), save :: spMAT
19
20contains
21
22 subroutine hecmw_solve_direct_mumps(hecMESH,hecMAT)
23 implicit none
24 type (hecmwst_local_mesh), intent(in) :: hecmesh
25 type (hecmwst_matrix ), intent(inout) :: hecmat
26 integer(kind=kint) :: spmat_type
27 integer(kind=kint) :: spmat_symtype
28 integer(kind=kint) :: mumps_job
29 integer(kind=kint) :: istat,myrank
30 real(kind=kreal) :: t1,t2,t3
31
32 call hecmw_mat_dump(hecmat, hecmesh)
33
34 t1=hecmw_wtime()
35 myrank=hecmw_comm_get_rank()
36
37 if (initialized .and. hecmat%Iarray(98) .eq. 1) then
38 mumps_job=-2
39 call hecmw_mumps_wrapper(spmat, mumps_job, istat)
40 if (istat < 0) then
41 write(*,*) 'ERROR: MUMPS returned with error', istat
42 stop
43 endif
44 call sparse_matrix_finalize(spmat)
45 initialized = .false.
46 endif
47
48 if (.not. initialized) then
49 spmat_type = sparse_matrix_type_coo
50 if (hecmat%symmetric) then
51 spmat_symtype = sparse_matrix_symtype_spd
52 else
53 spmat_symtype = sparse_matrix_symtype_asym
54 end if
55 call sparse_matrix_set_type(spmat, spmat_type, spmat_symtype)
56 mumps_job = -1
57 call hecmw_mumps_wrapper(spmat, mumps_job, istat)
58 if (istat < 0) then
59 write(*,*) 'ERROR: MUMPS returned with error', istat
60 stop
61 endif
62 initialized = .true.
63 hecmat%Iarray(98) = 1
64 endif
65
66 !* Flag to activate symbolic factorization: 1(yes) 0(no) hecMESH%Iarray(98)
67 !* Flag to activate numeric factorization: 1(yes) 0(no) hecMESH%Iarray(97)
68
69 if (hecmat%Iarray(98) .eq. 1) then
70 ! ANALYSIS and FACTORIZATION
71 call sparse_matrix_hec_init_prof(spmat, hecmat, hecmesh)
72 call sparse_matrix_hec_set_vals(spmat, hecmat)
73 !call sparse_matrix_dump(spMAT)
74 mumps_job=4
75 call hecmw_mumps_wrapper(spmat, mumps_job, istat)
76 if (istat < 0) then
77 write(*,*) 'ERROR: MUMPS returned with error', istat
78 stop
79 endif
80 if (myrank==0 .and. (spmat%iterlog > 0 .or. spmat%timelog > 0)) &
81 write(*,*) ' [MUMPS]: Analysis and Factorization completed.'
82 hecmat%Iarray(98) = 0
83 hecmat%Iarray(97) = 0
84 endif
85 if (hecmat%Iarray(97) .eq. 1) then
86 ! FACTORIZATION
87 call sparse_matrix_hec_set_vals(spmat, hecmat)
88 !call sparse_matrix_dump(spMAT)
89 mumps_job=2
90 call hecmw_mumps_wrapper(spmat, mumps_job, istat)
91 if (istat < 0) then
92 write(*,*) 'ERROR: MUMPS returned with error', istat
93 stop
94 endif
95 if (myrank==0 .and. (spmat%iterlog > 0 .or. spmat%timelog > 0)) &
96 write(*,*) ' [MUMPS]: Factorization completed.'
97 hecmat%Iarray(97) = 0
98 endif
99
100 t2=hecmw_wtime()
101
102 ! SOLUTION
103 call sparse_matrix_hec_set_rhs(spmat, hecmat)
104 mumps_job=3
105 call hecmw_mumps_wrapper(spmat, mumps_job, istat)
106 if (istat < 0) then
107 write(*,*) 'ERROR: MUMPS returned with error', istat
108 stop
109 endif
110 call sparse_matrix_hec_get_rhs(spmat, hecmat)
111 if (myrank==0 .and. (spmat%iterlog > 0 .or. spmat%timelog > 0)) &
112 write(*,*) ' [MUMPS]: Solution completed.'
113
114 t3=hecmw_wtime()
115 if (myrank==0 .and. spmat%timelog > 0) then
116 write(*,*) 'setup time : ',t2-t1
117 write(*,*) 'solve time : ',t3-t2
118 endif
119
120 call hecmw_mat_dump_solution(hecmat)
121
122 !call sparse_matrix_finalize(spMAT)
123 end subroutine hecmw_solve_direct_mumps
124
subroutine, public hecmw_mat_dump(hecmat, hecmesh)
subroutine, public hecmw_mat_dump_solution(hecmat)
This module provides linear equation solver interface for MUMPS.
subroutine, public hecmw_solve_direct_mumps(hecmesh, hecmat)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=4), parameter kreal
integer(kind=kint) function hecmw_comm_get_rank()
real(kind=kreal) function hecmw_wtime()
This module provides wrapper for parallel sparse direct solver MUMPS.
subroutine, public hecmw_mumps_wrapper(spmat, job, istat)
This module provides conversion routines between HEC data structure and DOF based sparse matrix struc...
subroutine, public sparse_matrix_hec_set_rhs(spmat, hecmat)
subroutine, public sparse_matrix_hec_set_vals(spmat, hecmat)
subroutine, public sparse_matrix_hec_init_prof(spmat, hecmat, hecmesh)
subroutine, public sparse_matrix_hec_get_rhs(spmat, hecmat)
This module provides DOF based sparse matrix data structure (CSR and COO)
integer(kind=kint), parameter, public sparse_matrix_type_coo
subroutine, public sparse_matrix_set_type(spmat, type, symtype)
subroutine, public sparse_matrix_finalize(spmat)
integer(kind=kint), parameter, public sparse_matrix_symtype_asym
integer(kind=kint), parameter, public sparse_matrix_symtype_spd