FrontISTR  5.2.0
Large-scale structural analysis program with finit element method
hecmw_mpc_prepost.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 !-------------------------------------------------------------------------------
5 
7  use hecmw_util
13  implicit none
14 
15  private
16  public :: hecmw_mpc_mat_init
18  public :: hecmw_mpc_mat_finalize
20  public :: hecmw_mpc_mat_ass
21  public :: hecmw_mpc_trans_rhs
22  public :: hecmw_mpc_tback_sol
23  public :: hecmw_mpc_trans_mass
24  public :: hecmw_mpc_tback_eigvec
25  public :: hecmw_mpc_mark_slave
26 
27 contains
28 
29  !C
30  !C***
31  !C*** hecmw_mpc_mat_init
32  !C***
33  !C
34  subroutine hecmw_mpc_mat_init(hecMESH, hecMAT, hecMESHmpc, hecMATmpc)
35  implicit none
36  type (hecmwst_local_mesh), intent(inout), target :: hecmesh
37  type (hecmwst_matrix), intent(in), target :: hecmat
38  type (hecmwst_local_mesh), pointer :: hecmeshmpc
39  type (hecmwst_matrix), pointer :: hecmatmpc
40  integer(kind=kint) :: totalmpc, mpc_method, solver_type
41 
42  totalmpc = hecmesh%mpc%n_mpc
43  call hecmw_allreduce_i1 (hecmesh, totalmpc, hecmw_sum)
44 
45  if (totalmpc == 0) then
46  hecmeshmpc => hecmesh
47  hecmatmpc => hecmat
48  return
49  endif
50 
51  call hecmw_mpc_scale(hecmesh)
52 
53  mpc_method = hecmw_mat_get_mpc_method(hecmat)
54  if (mpc_method < 1 .or. 3 < mpc_method) then
55  solver_type = hecmw_mat_get_solver_type(hecmat)
56  if (solver_type > 1) then ! DIRECT SOLVER
57  mpc_method = 1 ! default: penalty
58  else ! ITERATIVE SOLVER
59  mpc_method = 3 ! default: elimination
60  endif
61  call hecmw_mat_set_mpc_method(hecmat, mpc_method)
62  endif
63 
64  if (mpc_method == 2) then
65  write(*,*) 'WARNING: MPCMETHOD=2 (MPCCG) is deprecated; may not work correctly'
66  ! MPC_METHOD = 3
67  ! call hecmw_mat_set_mpc_method(hecMAT, MPC_METHOD)
68  endif
69 
70  select case (mpc_method)
71  case (1) ! penalty
72  hecmeshmpc => hecmesh
73  hecmatmpc => hecmat
74  case (2) ! MPCCG
75  hecmeshmpc => hecmesh
76  hecmatmpc => hecmat
77  case (3) ! elimination
78  allocate(hecmeshmpc)
79  call hecmw_mpc_mesh_copy(hecmesh, hecmeshmpc)
80  allocate(hecmatmpc)
81  call hecmw_mat_init(hecmatmpc)
82  end select
83 
84  end subroutine hecmw_mpc_mat_init
85 
86  !C
87  !C***
88  !C*** hecmw_mpc_mat_init_explicit
89  !C***
90  !C
91  subroutine hecmw_mpc_mat_init_explicit(hecMESH, hecMAT, hecMATmpc)
92  implicit none
93  type (hecmwst_local_mesh), intent(inout), target :: hecmesh
94  type (hecmwst_matrix), intent(in), target :: hecmat
95  type (hecmwst_matrix), pointer :: hecmatmpc
96  integer(kind=kint) :: totalmpc, mpc_method
97 
98  totalmpc = hecmesh%mpc%n_mpc
99  call hecmw_allreduce_i1 (hecmesh, totalmpc, hecmw_sum)
100 
101  if (totalmpc == 0) then
102  hecmatmpc => hecmat
103  return
104  endif
105 
106  call hecmw_mpc_scale(hecmesh)
107 
108  ! Force MPC_METHOD=3
109  mpc_method = 3
110  call hecmw_mat_set_mpc_method(hecmat, mpc_method)
111 
112  allocate(hecmatmpc)
113  call hecmw_mat_init(hecmatmpc)
114 
115  hecmatmpc%N = hecmat%N
116  hecmatmpc%NP = hecmat%NP
117  hecmatmpc%NDOF = hecmat%NDOF
118  allocate(hecmatmpc%B(size(hecmat%B)))
119  allocate(hecmatmpc%X(size(hecmat%X)))
120  end subroutine hecmw_mpc_mat_init_explicit
121 
122  !C
123  !C***
124  !C*** hecmw_mpc_mat_finalize
125  !C***
126  !C
127  subroutine hecmw_mpc_mat_finalize(hecMESH, hecMAT, hecMESHmpc, hecMATmpc)
128  implicit none
129  type (hecmwst_local_mesh), intent(in) :: hecmesh
130  type (hecmwst_matrix), intent(in) :: hecmat
131  type (hecmwst_local_mesh), pointer :: hecmeshmpc
132  type (hecmwst_matrix), pointer :: hecmatmpc
133  integer(kind=kint) :: totalmpc, mpc_method
134 
135  totalmpc = hecmesh%mpc%n_mpc
136  call hecmw_allreduce_i1 (hecmesh, totalmpc, hecmw_sum)
137 
138  if (totalmpc == 0) then
139  nullify(hecmeshmpc)
140  nullify(hecmatmpc)
141  return
142  endif
143 
144  mpc_method = hecmw_mat_get_mpc_method(hecmat)
145 
146  select case (mpc_method)
147  case (1) ! penalty
148  nullify(hecmeshmpc)
149  nullify(hecmatmpc)
150  case (2) ! MPCCG
151  nullify(hecmeshmpc)
152  nullify(hecmatmpc)
153  case (3) ! elimination
154  call hecmw_mpc_mesh_free(hecmeshmpc)
155  deallocate(hecmeshmpc)
156  nullify(hecmeshmpc)
157  call hecmw_mat_finalize(hecmatmpc)
158  deallocate(hecmatmpc)
159  nullify(hecmatmpc)
160  end select
161 
162  end subroutine hecmw_mpc_mat_finalize
163 
164  !C
165  !C***
166  !C*** hecmw_mpc_mat_finalize_explicit
167  !C***
168  !C
169  subroutine hecmw_mpc_mat_finalize_explicit(hecMESH, hecMAT, hecMATmpc)
170  implicit none
171  type (hecmwst_local_mesh), intent(in) :: hecmesh
172  type (hecmwst_matrix), intent(in) :: hecmat
173  type (hecmwst_matrix), pointer :: hecmatmpc
174  integer(kind=kint) :: totalmpc, mpc_method
175 
176  totalmpc = hecmesh%mpc%n_mpc
177  call hecmw_allreduce_i1 (hecmesh, totalmpc, hecmw_sum)
178 
179  if (totalmpc == 0) then
180  nullify(hecmatmpc)
181  return
182  endif
183 
184  mpc_method = hecmw_mat_get_mpc_method(hecmat)
185 
186  select case (mpc_method)
187  case (1) ! penalty
188  nullify(hecmatmpc)
189  case (2) ! MPCCG
190  nullify(hecmatmpc)
191  case (3) ! elimination
192  call hecmw_mat_finalize(hecmatmpc)
193  deallocate(hecmatmpc)
194  nullify(hecmatmpc)
195  end select
196 
197  end subroutine hecmw_mpc_mat_finalize_explicit
198 
199  !C
200  !C***
201  !C*** hecmw_mpc_mat_ass
202  !C***
203  !C
204  subroutine hecmw_mpc_mat_ass(hecMESH, hecMAT, hecMESHmpc, hecMATmpc)
205  implicit none
206  type (hecmwst_local_mesh), intent(in) :: hecmesh
207  type (hecmwst_matrix), intent(inout) :: hecmat
208  type (hecmwst_local_mesh), pointer :: hecmeshmpc
209  type (hecmwst_matrix), pointer :: hecmatmpc
210  integer(kind=kint) :: totalmpc, mpc_method
211 
212  totalmpc = hecmesh%mpc%n_mpc
213  call hecmw_allreduce_i1 (hecmesh, totalmpc, hecmw_sum)
214 
215  if (totalmpc == 0) return
216 
217  mpc_method = hecmw_mat_get_mpc_method(hecmat)
218 
219  select case (mpc_method)
220  case (1) ! penalty
221  !if (hecMESH%my_rank.eq.0) write(0,*) "MPC Method: Penalty"
222  call hecmw_mat_ass_equation ( hecmesh, hecmat )
223  case (2) ! MPCCG
224  !if (hecMESH%my_rank.eq.0) write(0,*) "MPC Method: MPC-CG"
225  case (3) ! elimination
226  !if (hecMESH%my_rank.eq.0) write(0,*) "MPC Method: Elimination"
227  call hecmw_trimatmul_ttkt_mpc(hecmeshmpc, hecmat, hecmatmpc)
228  end select
229 
230  end subroutine hecmw_mpc_mat_ass
231 
232 
233  !C
234  !C***
235  !C*** hecmw_mpc_trans_rhs
236  !C***
237  !C
238  subroutine hecmw_mpc_trans_rhs(hecMESH, hecMAT, hecMATmpc)
239  implicit none
240  type (hecmwst_local_mesh), intent(inout) :: hecmesh
241  type (hecmwst_matrix), intent(inout) :: hecmat
242  type (hecmwst_matrix), pointer :: hecmatmpc
243  real(kind=kreal), allocatable :: btmp(:)
244  real(kind=kreal) :: time_dumm
245  integer(kind=kint) :: totalmpc, mpc_method, i
246 
247  totalmpc = hecmesh%mpc%n_mpc
248  call hecmw_allreduce_i1 (hecmesh, totalmpc, hecmw_sum)
249 
250  if (totalmpc == 0) return
251 
252  mpc_method = hecmw_mat_get_mpc_method(hecmat)
253 
254  select case (mpc_method)
255  case (1) ! penalty
256  call hecmw_mat_ass_equation_rhs ( hecmesh, hecmatmpc )
257  case (2) ! MPCCG
258  allocate(btmp(hecmat%NP*hecmat%NDOF))
259  do i = 1, hecmat%NP*hecmat%NDOF
260  btmp(i) = hecmat%B(i)
261  enddo
262  call hecmw_trans_b(hecmesh, hecmat, btmp, hecmatmpc%B, time_dumm)
263  deallocate(btmp)
264  case (3) ! elimination
265  call hecmw_trans_b(hecmesh, hecmat, hecmat%B, hecmatmpc%B, time_dumm)
266  hecmatmpc%Iarray=hecmat%Iarray
267  hecmatmpc%Rarray=hecmat%Rarray
268  end select
269 
270  end subroutine hecmw_mpc_trans_rhs
271 
272  !C
273  !C***
274  !C*** hecmw_mpc_tback_sol
275  !C***
276  !C
277  subroutine hecmw_mpc_tback_sol(hecMESH, hecMAT, hecMATmpc)
278  implicit none
279  type (hecmwst_local_mesh), intent(in) :: hecmesh
280  type (hecmwst_matrix), intent(inout) :: hecmat
281  type (hecmwst_matrix), pointer :: hecmatmpc
282  real(kind=kreal) :: time_dumm
283  integer(kind=kint) :: totalmpc, mpc_method, i
284 
285  totalmpc = hecmesh%mpc%n_mpc
286  call hecmw_allreduce_i1 (hecmesh, totalmpc, hecmw_sum)
287 
288  if (totalmpc == 0) return
289 
290  mpc_method = hecmw_mat_get_mpc_method(hecmat)
291 
292  select case (mpc_method)
293  case (1) ! penalty
294  ! do nothing
295  case (2) ! MPCCG
296  call hecmw_tback_x(hecmesh, hecmat%NDOF, hecmat%X, time_dumm)
297  case (3) ! elimination
298  do i = 1, size(hecmat%X)
299  hecmat%X(i) = hecmatmpc%X(i)
300  enddo
301  call hecmw_tback_x(hecmesh, hecmat%NDOF, hecmat%X, time_dumm)
302  hecmat%Iarray=hecmatmpc%Iarray
303  hecmat%Rarray=hecmatmpc%Rarray
304  end select
305  end subroutine hecmw_mpc_tback_sol
306 
307  !C
308  !C***
309  !C*** hecmw_mpc_trans_mass
310  !C***
311  !C
312  subroutine hecmw_mpc_trans_mass(hecMESH, hecMAT, mass)
313  implicit none
314  type (hecmwst_local_mesh), intent(inout) :: hecmesh
315  type (hecmwst_matrix), intent(inout) :: hecmat
316  real(kind=kreal), intent(inout) :: mass(:)
317 
318  real(kind=kreal), allocatable :: mtmp(:)
319  real(kind=kreal) :: time_dumm
320  integer(kind=kint) :: totalmpc, mpc_method, i
321 
322  totalmpc = hecmesh%mpc%n_mpc
323  call hecmw_allreduce_i1 (hecmesh, totalmpc, hecmw_sum)
324 
325  if (totalmpc == 0) return
326 
327  mpc_method = hecmw_mat_get_mpc_method(hecmat)
328 
329  select case (mpc_method)
330  case (1) ! penalty
331  ! do nothing
332  case (2,3) ! MPCCG or elimination
333  allocate(mtmp(hecmat%NP*hecmat%NDOF))
334  !C-- {Mt} = [T'] {w}
335  call hecmw_ttvec(hecmesh, hecmat%NDOF, mass, mtmp, time_dumm)
336  do i = 1, hecmat%NP*hecmat%NDOF
337  mass(i) = mtmp(i)
338  enddo
339  deallocate(mtmp)
340  end select
341 
342  end subroutine hecmw_mpc_trans_mass
343 
344  !C
345  !C***
346  !C*** hecmw_mpc_tback_eigvec
347  !C***
348  !C
349  subroutine hecmw_mpc_tback_eigvec(hecMESH, hecMAT, neig, eigvec)
350  implicit none
351  type (hecmwst_local_mesh), intent(in) :: hecmesh
352  type (hecmwst_matrix), intent(inout) :: hecmat
353  integer(kind=kint), intent(in) :: neig
354  real(kind=kreal), intent(inout) :: eigvec(:,:)
355 
356  real(kind=kreal) :: time_dumm
357  integer(kind=kint) :: totalmpc, mpc_method, i
358 
359  totalmpc = hecmesh%mpc%n_mpc
360  call hecmw_allreduce_i1 (hecmesh, totalmpc, hecmw_sum)
361 
362  if (totalmpc == 0) return
363 
364  mpc_method = hecmw_mat_get_mpc_method(hecmat)
365 
366  select case (mpc_method)
367  case (1) ! penalty
368  ! do nothing
369  case (2,3) ! MPCCG or elimination
370  do i = 1, neig
371  call hecmw_tback_x(hecmesh, hecmat%NDOF, eigvec(:,i), time_dumm)
372  !!! need normalization???
373  enddo
374  end select
375  end subroutine hecmw_mpc_tback_eigvec
376 
377  !C
378  !C***
379  !C*** hecmw_mpc_mark_slave
380  !C***
381  !C
382  subroutine hecmw_mpc_mark_slave(hecMESH, hecMAT, mark)
383  implicit none
384  type (hecmwst_local_mesh), intent(in) :: hecmesh
385  type (hecmwst_matrix), intent(inout) :: hecmat
386  integer(kind=kint), intent(out) :: mark(:)
387 
388  integer(kind=kint) :: ndof, i, j, k, kk
389 
390  ndof = hecmat%NDOF
391  mark(:) = 0
392  outer: do i = 1, hecmesh%mpc%n_mpc
393  do j = hecmesh%mpc%mpc_index(i-1)+1, hecmesh%mpc%mpc_index(i)
394  if (hecmesh%mpc%mpc_dof(j) > ndof) cycle outer
395  enddo
396  k = hecmesh%mpc%mpc_index(i-1)+1
397  kk = ndof * (hecmesh%mpc%mpc_item(k) - 1) + hecmesh%mpc%mpc_dof(k)
398  mark(kk) = 1
399  enddo outer
400  end subroutine hecmw_mpc_mark_slave
401 
402  !C
403  !C***
404  !C*** hecmw_mpc_scale
405  !C***
406  !C
407  subroutine hecmw_mpc_scale(hecMESH)
408  implicit none
409  type (hecmwst_local_mesh), intent(inout) :: hecmesh
410  integer(kind=kint) :: i, j, k
411  real(kind=kreal) :: wval
412 
413  !$omp parallel default(none),private(i,j,k,WVAL),shared(hecMESH)
414  !$omp do
415  do i = 1, hecmesh%mpc%n_mpc
416  k = hecmesh%mpc%mpc_index(i-1)+1
417  wval = 1.d0 / hecmesh%mpc%mpc_val(k)
418  hecmesh%mpc%mpc_val(k) = 1.d0
419  do j = hecmesh%mpc%mpc_index(i-1)+2, hecmesh%mpc%mpc_index(i)
420  hecmesh%mpc%mpc_val(j) = hecmesh%mpc%mpc_val(j) * wval
421  enddo
422  hecmesh%mpc%mpc_const(i) = hecmesh%mpc%mpc_const(i) * wval
423  enddo
424  !$omp end do
425  !$omp end parallel
426 
427  end subroutine hecmw_mpc_scale
428 
429 
430  !C
431  !C***
432  !C*** hecmw_trans_b
433  !C***
434  !C
435  subroutine hecmw_trans_b(hecMESH, hecMAT, B, BT, COMMtime)
436  implicit none
437  type (hecmwst_local_mesh), intent(in) :: hecmesh
438  type (hecmwst_matrix), intent(in) :: hecmat
439  real(kind=kreal), intent(in) :: b(:)
440  real(kind=kreal), intent(out), target :: bt(:)
441  real(kind=kreal), intent(inout) :: commtime
442 
443  real(kind=kreal), allocatable :: w(:)
444  real(kind=kreal), pointer :: xg(:)
445  integer(kind=kint) :: ndof, i, j, k, kk, flg_bak
446 
447  ndof = hecmat%NDOF
448 
449  allocate(w(hecmesh%n_node * ndof))
450 
451  !C===
452  !C +---------------------------+
453  !C | {bt}= [T']({b} - [A]{xg}) |
454  !C +---------------------------+
455  !C===
456  xg => bt
457  do i = 1, hecmat%N * ndof
458  xg(i) = 0.d0
459  enddo
460 
461  !C-- Generate {xg} from mpc_const
462  !$omp parallel default(none),private(i,k,kk),shared(hecMESH,XG),firstprivate(ndof)
463  !$omp do
464  outer: do i = 1, hecmesh%mpc%n_mpc
465  do j = hecmesh%mpc%mpc_index(i-1)+1, hecmesh%mpc%mpc_index(i)
466  if (hecmesh%mpc%mpc_dof(j) > ndof) cycle outer
467  enddo
468  k = hecmesh%mpc%mpc_index(i-1) + 1
469  kk = ndof * (hecmesh%mpc%mpc_item(k) - 1) + hecmesh%mpc%mpc_dof(k)
470  xg(kk) = hecmesh%mpc%mpc_const(i)
471  enddo outer
472  !$omp end do
473  !$omp end parallel
474 
475  !C-- {w} = {b} - [A]{xg}
476  flg_bak = hecmw_mat_get_flag_mpcmatvec(hecmat)
477  call hecmw_mat_set_flag_mpcmatvec(hecmat, 0)
478  call hecmw_matresid(hecmesh, hecmat, xg, b, w, commtime)
479  call hecmw_mat_set_flag_mpcmatvec(hecmat, flg_bak)
480 
481  !C-- {bt} = [T'] {w}
482  call hecmw_ttvec(hecmesh, ndof, w, bt, commtime)
483 
484  deallocate(w)
485  end subroutine hecmw_trans_b
486 
487 
488  !C
489  !C***
490  !C*** hecmw_tback_x
491  !C***
492  !C
493  subroutine hecmw_tback_x(hecMESH, ndof, X, COMMtime)
494  implicit none
495  type (hecmwst_local_mesh), intent(in) :: hecmesh
496  integer(kind=kint), intent(in) :: ndof
497  real(kind=kreal), intent(inout) :: x(:)
498  real(kind=kreal), intent(inout) :: commtime
499 
500  real(kind=kreal), allocatable :: w(:)
501  integer(kind=kint) :: i, j, k, kk
502 
503  allocate(w(hecmesh%n_node * ndof))
504 
505  !C-- {tx} = [T]{x}
506  call hecmw_tvec(hecmesh, ndof, x, w, commtime)
507 
508  !C-- {x} = {tx} + {xg}
509  !$omp parallel default(none),private(i,k,kk),shared(hecMESH,X,W),firstprivate(ndof)
510  !$omp do
511  do i= 1, hecmesh%nn_internal * ndof
512  x(i)= w(i)
513  enddo
514  !$omp end do
515 
516  !$omp do
517  outer: do i = 1, hecmesh%mpc%n_mpc
518  do j = hecmesh%mpc%mpc_index(i-1)+1, hecmesh%mpc%mpc_index(i)
519  if (hecmesh%mpc%mpc_dof(j) > ndof) cycle outer
520  enddo
521  k = hecmesh%mpc%mpc_index(i-1) + 1
522  kk = ndof * (hecmesh%mpc%mpc_item(k) - 1) + hecmesh%mpc%mpc_dof(k)
523  x(kk) = x(kk) + hecmesh%mpc%mpc_const(i)
524  enddo outer
525  !$omp end do
526  !$omp end parallel
527 
528  deallocate(w)
529 
530  if (ndof == 3) then
531  call hecmw_update_3_r(hecmesh, x, hecmesh%n_node)
532  else if (ndof == 2) then
533  call hecmw_update_2_r(hecmesh, x, hecmesh%n_node)
534  else
535  call hecmw_update_m_r(hecmesh, x, hecmesh%n_node, ndof)
536  endif
537  end subroutine hecmw_tback_x
538 
539  subroutine hecmw_mpc_mesh_copy(src, dst)
540  implicit none
541  type (hecmwst_local_mesh), intent(in) :: src
542  type (hecmwst_local_mesh), intent(out) :: dst
543  dst%zero = src%zero
544  dst%MPI_COMM = src%MPI_COMM
545  dst%PETOT = src%PETOT
546  dst%PEsmpTOT = src%PEsmpTOT
547  dst%my_rank = src%my_rank
548  dst%n_subdomain = src%n_subdomain
549  dst%n_node = src%n_node
550  dst%nn_internal = src%nn_internal
551  dst%n_elem = src%n_elem
552  dst%ne_internal = src%ne_internal
553  dst%n_elem_type = src%n_elem_type
554  dst%n_dof = src%n_dof
555  dst%n_neighbor_pe = src%n_neighbor_pe
556  if (src%n_neighbor_pe > 0) then
557  allocate(dst%neighbor_pe(dst%n_neighbor_pe))
558  dst%neighbor_pe(:) = src%neighbor_pe(:)
559  allocate(dst%import_index(0:dst%n_neighbor_pe))
560  dst%import_index(:)= src%import_index(:)
561  allocate(dst%export_index(0:dst%n_neighbor_pe))
562  dst%export_index(:)= src%export_index(:)
563  allocate(dst%import_item(dst%import_index(dst%n_neighbor_pe)))
564  dst%import_item(:) = src%import_item(:)
565  allocate(dst%export_item(dst%export_index(dst%n_neighbor_pe)))
566  dst%export_item(:) = src%export_item(:)
567  endif
568  allocate(dst%global_node_ID(dst%n_node))
569  dst%global_node_ID(1:dst%n_node) = src%global_node_ID(1:dst%n_node)
570  allocate(dst%node_ID(2*dst%n_node))
571  dst%node_ID(1:2*dst%n_node) = src%node_ID(1:2*dst%n_node)
572  allocate(dst%elem_type_item(dst%n_elem_type))
573  dst%elem_type_item(:) = src%elem_type_item(:)
574  !
575  dst%mpc%n_mpc = src%mpc%n_mpc
576  dst%mpc%mpc_index => src%mpc%mpc_index
577  dst%mpc%mpc_item => src%mpc%mpc_item
578  dst%mpc%mpc_dof => src%mpc%mpc_dof
579  dst%mpc%mpc_val => src%mpc%mpc_val
580  dst%mpc%mpc_const => src%mpc%mpc_const
581  !
582  dst%node_group%n_grp = src%node_group%n_grp
583  dst%node_group%n_bc = src%node_group%n_bc
584  dst%node_group%grp_name => src%node_group%grp_name
585  dst%node_group%grp_index => src%node_group%grp_index
586  dst%node_group%grp_item => src%node_group%grp_item
587  dst%node_group%bc_grp_ID => src%node_group%bc_grp_ID
588  dst%node_group%bc_grp_type => src%node_group%bc_grp_type
589  dst%node_group%bc_grp_index => src%node_group%bc_grp_index
590  dst%node_group%bc_grp_dof => src%node_group%bc_grp_dof
591  dst%node_group%bc_grp_val => src%node_group%bc_grp_val
592  !
593  dst%node => src%node
594  end subroutine hecmw_mpc_mesh_copy
595 
596  subroutine hecmw_mpc_mesh_free(hecMESH)
597  implicit none
598  type (hecmwst_local_mesh), intent(inout) :: hecmesh
599  if (hecmesh%n_neighbor_pe > 1) then
600  deallocate(hecmesh%neighbor_pe)
601  deallocate(hecmesh%import_index)
602  deallocate(hecmesh%export_index)
603  deallocate(hecmesh%import_item)
604  deallocate(hecmesh%export_item)
605  endif
606  deallocate(hecmesh%global_node_ID)
607  deallocate(hecmesh%node_ID)
608  deallocate(hecmesh%elem_type_item)
609  end subroutine hecmw_mpc_mesh_free
610 end module hecmw_mpc_prepost
subroutine, public hecmw_trimatmul_ttkt_mpc(hecMESH, hecMAT, hecTKT)
subroutine, public hecmw_mat_ass_equation_rhs(hecMESH, hecMAT)
subroutine, public hecmw_mat_ass_equation(hecMESH, hecMAT)
integer(kind=kint) function, public hecmw_mat_get_solver_type(hecMAT)
integer(kind=kint) function, public hecmw_mat_get_flag_mpcmatvec(hecMAT)
subroutine, public hecmw_mat_init(hecMAT)
subroutine, public hecmw_mat_finalize(hecMAT)
subroutine, public hecmw_mat_set_flag_mpcmatvec(hecMAT, flag_mpcmatvec)
integer(kind=kint) function, public hecmw_mat_get_mpc_method(hecMAT)
subroutine, public hecmw_mat_set_mpc_method(hecMAT, mpc_method)
subroutine, public hecmw_mpc_tback_sol(hecMESH, hecMAT, hecMATmpc)
subroutine, public hecmw_mpc_mat_ass(hecMESH, hecMAT, hecMESHmpc, hecMATmpc)
subroutine, public hecmw_mpc_mat_init(hecMESH, hecMAT, hecMESHmpc, hecMATmpc)
subroutine, public hecmw_mpc_mat_init_explicit(hecMESH, hecMAT, hecMATmpc)
subroutine, public hecmw_mpc_trans_mass(hecMESH, hecMAT, mass)
subroutine, public hecmw_mpc_mat_finalize(hecMESH, hecMAT, hecMESHmpc, hecMATmpc)
subroutine, public hecmw_mpc_mat_finalize_explicit(hecMESH, hecMAT, hecMATmpc)
subroutine, public hecmw_mpc_tback_eigvec(hecMESH, hecMAT, neig, eigvec)
subroutine, public hecmw_mpc_trans_rhs(hecMESH, hecMAT, hecMATmpc)
subroutine, public hecmw_mpc_mark_slave(hecMESH, hecMAT, mark)
subroutine, public hecmw_ttvec(hecMESH, ndof, X, Y, COMMtime)
subroutine, public hecmw_tvec(hecMESH, ndof, X, Y, COMMtime)
subroutine, public hecmw_matresid(hecMESH, hecMAT, X, B, R, COMMtime)
I/O and Utility.
Definition: hecmw_util_f.F90:7
integer(kind=kint), parameter hecmw_sum
integer(kind=4), parameter kreal
subroutine hecmw_allreduce_i1(hecMESH, s, ntag)
subroutine hecmw_update_2_r(hecMESH, val, n)
subroutine hecmw_update_3_r(hecMESH, val, n)
subroutine hecmw_update_m_r(hecMESH, val, n, m)