FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
fstr_rcap_io.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!-------------------------------------------------------------------------------
6 use m_fstr
7
8#ifdef WITH_REVOCAP
9 use rcapf
10#endif
11
12 public :: fstr_rcap_initialize ! call after fstr_setup
13 public :: fstr_rcap_finalize ! call before hecmw_finalize
14 public :: fstr_rcap_send
15 public :: fstr_rcap_get
16
17#ifdef WITH_REVOCAP
18 integer(kind=kint), private :: n_iter, n_rcap
19 real(kind=kreal), private :: t_start, t_end, t_rcap
20 real(kind=kreal), private :: t_start_all, t_end_all, t_rcap_all
21#endif
22
23contains
24
25 !------------------------------------------------------------------------------
26 subroutine fstr_rcap_initialize( hecMESH, fstrPARAM, fstrCPL )
27 implicit none
28 type( hecmwst_local_mesh ) :: hecmesh
29 type( fstr_param ) :: fstrparam
30 type( fstr_couple ) :: fstrcpl
31
32#ifdef WITH_REVOCAP
33 character( len=16) :: portfile
34 integer(kind=kint) :: myrank
35 integer(kind=kint) :: i,err,nid,nid_old
36 real(kind=kreal) :: t_s, t_e
37
38 write(idbg,*) "fstr_rcap_initialize: start"
39 t_start_all = hecmw_wtime()
40
41 if( fstrparam%fg_couple /= 1 ) return
42
43 write(idbg,*) "fstr_rcap_initialize: calling rcapf_init_solid_solver"
44
45 portfile='port'//char(0)
46
47 t_s = hecmw_wtime()
48 call rcapf_init_solid_solver( hecmesh%my_rank, portfile )
49 t_e = hecmw_wtime()
50 t_rcap_all = t_e - t_s
51
52 write(idbg,*) "fstr_rcap_initialize: calling rcapf_get_num_of_matching_node"
53
54 fstrcpl%dof = 3
55
56 t_s = hecmw_wtime()
57 call rcapf_get_num_of_matching_node( fstrcpl%coupled_node_n )
58 t_e = hecmw_wtime()
59 t_rcap_all = t_rcap_all + (t_e - t_s)
60
61 fstrcpl%ndof = fstrcpl%coupled_node_n * fstrcpl%dof
62 allocate( fstrcpl%coupled_node( fstrcpl%coupled_node_n ), stat=err)
63 if( err /= 0 ) goto 1000
64 allocate( fstrcpl%trac( fstrcpl%ndof ), stat=err)
65 if( err /= 0 ) goto 1000
66 allocate( fstrcpl%disp( fstrcpl%ndof ), stat=err)
67 if( err /= 0 ) goto 1000
68 allocate( fstrcpl%velo( fstrcpl%ndof ), stat=err)
69 if( err /= 0 ) goto 1000
70 allocate( fstrcpl%accel( fstrcpl%ndof ), stat=err)
71 if( err /= 0 ) goto 1000
72 allocate( fstrcpl%index( hecmesh%n_node_gross ), stat=err)
73 if( err /= 0 ) goto 1000
74 fstrcpl%trac = 0.0d0
75 fstrcpl%index = -1
76
77 write(idbg,*) "fstr_rcap_initialize: calling rcapf_get_mathcing_node_id"
78
79 t_s = hecmw_wtime()
80 call rcapf_get_matching_node_id( fstrcpl%coupled_node, fstrcpl%coupled_node_n )
81 t_e = hecmw_wtime()
82 t_rcap_all = t_rcap_all + (t_e - t_s)
83
84 write(idbg,*) "fstr_rcap_initialize: converting to local id: ", fstrcpl%coupled_node_n
85
86 do i=1, fstrcpl%coupled_node_n
87 nid = fstrcpl%coupled_node(i)
88 if( nid <= 0 .or. nid>hecmesh%n_node_gross ) then
89 write(*,*) '##Fatal error in fstr_rcap_initialize ', i, nid
90 call hecmw_abort( hecmw_comm_get_comm());
91 endif
92 if( hecmesh%n_refine > 0 ) then
93 nid_old = nid
94 if( associated( hecmesh%node_old2new ) ) then
95 nid = hecmesh%node_old2new( nid ) + 1
96 endif
97 write(idbg,*) nid_old, nid
98 endif
99 fstrcpl%index( nid ) = i
100 end do
101
102 n_iter = 0
103 n_rcap = 0
104
105 write(idbg,*) "fstr_rcap_initialize: end"
106
107 return
108 1000 write(*,*) "##Error : not enough memory"
109 call hecmw_abort( hecmw_comm_get_comm() )
110
111#else
112
113 if( fstrparam%fg_couple == 0 ) return
114
115 if( hecmw_comm_get_rank() == 0 ) then
116 write(*,*) "##Error : REVOCAP functions are not supported"
117 end if
118 call hecmw_abort( hecmw_comm_get_comm() )
119#endif
120 end subroutine fstr_rcap_initialize
121
122 !------------------------------------------------------------------------------
123 subroutine fstr_rcap_finalize( fstrPARAM, fstrCPL )
124 implicit none
125 type( fstr_param ) :: fstrparam
126 type( fstr_couple ) :: fstrcpl
127
128#ifdef WITH_REVOCAP
129 real(kind=kreal) :: t_tot, t_tot_avg, t_rcap_avg, tr
130 real(kind=kreal) :: t_tot_all, tr_all
131 real(kind=kreal) :: t_s, t_e
132
133 write(idbg,*) "fstr_rcap_finalize: start"
134
135 if( fstrparam%fg_couple /= 1 ) return
136 deallocate( fstrcpl%coupled_node )
137 deallocate( fstrcpl%trac )
138 deallocate( fstrcpl%disp )
139 deallocate( fstrcpl%velo )
140 deallocate( fstrcpl%accel )
141 deallocate( fstrcpl%index )
142
143 write(idbg,*) "fstr_rcap_finalize: calling rcapf_finalize"
144
145 t_s = hecmw_wtime()
146 call rcapf_finalize()
147 t_e = hecmw_wtime()
148 t_rcap_all = t_rcap_all + (t_e - t_s)
149
150 t_tot = t_end - t_start
151 t_tot_avg = t_tot / n_iter
152 t_rcap_avg = t_rcap / n_rcap
153 tr = t_rcap_avg / t_tot_avg * 100.d0
154
155 write(idbg,'(a,f11.3,a,i0,a,f11.3,a)') &
156 & 'fstr + rcap:', t_tot,'s for ',n_iter,' iters / avg. ', t_tot_avg,'s/iter'
157 write(idbg,'(a,f11.3,a,i0,a,f11.3,a)') &
158 & ' rcap:',t_rcap,'s for ',n_rcap,' calls / avg. ',t_rcap_avg,'s/call'
159 write(idbg,'(a,f11.3,a)') &
160 & ' rcap ratio:',tr,'%/iter'
161
162 t_end_all = hecmw_wtime()
163
164 t_tot_all = t_end_all - t_start_all
165 tr_all = t_rcap_all / t_tot_all * 100.d0
166
167 write(idbg,'(a,f11.3,a)') 'fstr total:',t_tot_all,'s'
168 write(idbg,'(a,f11.3,a)') 'rcap total:',t_rcap_all,'s'
169 write(idbg,'(a,f11.3,a)') 'rcap ratio:',tr_all,'%'
170
171 write(idbg,*) "fstr_rcap_finalize: end"
172
173#else
174
175 if( fstrparam%fg_couple == 0 ) return
176
177 if( hecmw_comm_get_rank() == 0 ) then
178 write(*,*) "##Error : REVOCAP functions are not supported"
179 end if
180 call hecmw_abort( hecmw_comm_get_comm() )
181#endif
182 end subroutine fstr_rcap_finalize
183 !------------------------------------------------------------------------------
184 subroutine fstr_rcap_send( fstrCPL )
185 implicit none
186 type( fstr_couple ) :: fstrcpl
187
188#ifdef WITH_REVOCAP
189 write(idbg,*) "fstr_rcap_send: start"
190
191 ! call rcapf_set_disp( fstrCPL%coupled_node, &
192 ! fstrCPL%coupled_node_n, &
193 ! fstrCPL%disp, fstrCPL%ndof )
194
195 call rcapf_set_velo( fstrcpl%coupled_node, &
196 fstrcpl%coupled_node_n, &
197 fstrcpl%velo, fstrcpl%ndof )
198
199 ! call rcapf_set_accel( fstrCPL%coupled_node, &
200 ! fstrCPL%coupled_node_n, &
201 ! fstrCPL%accel, fstrCPL%ndof )
202
203 write(idbg,*) "fstr_rcap_send: end"
204
205#else
206
207 if( hecmw_comm_get_rank() == 0 ) then
208 write(*,*) "##Error : REVOCAP functions are not supported"
209 end if
210 call hecmw_abort( hecmw_comm_get_comm() )
211#endif
212 end subroutine fstr_rcap_send
213 !------------------------------------------------------------------------------
214 subroutine fstr_rcap_get( fstrCPL )
215 implicit none
216 type( fstr_couple ) :: fstrcpl
217
218#ifdef WITH_REVOCAP
219 real(kind=kreal) :: t_s, t_e
220
221 write(idbg,*) "fstr_rcap_get: start"
222
223 if (n_rcap == 0) then
224 t_start = hecmw_wtime()
225 t_rcap = 0
226 else
227 t_end = hecmw_wtime()
228 n_iter = n_iter + 1
229 end if
230
231 write(idbg,*) "fstr_rcap_get: calling rcapf_get_trac"
232
233 t_s = hecmw_wtime()
234 call rcapf_get_trac( fstrcpl%coupled_node, &
235 fstrcpl%coupled_node_n, &
236 fstrcpl%trac, fstrcpl%ndof )
237 t_e = hecmw_wtime()
238
239 t_rcap = t_rcap + (t_e - t_s)
240 n_rcap = n_rcap + 1
241
242 t_rcap_all = t_rcap_all + (t_e - t_s)
243
244 write(idbg,*) "fstr_rcap_get: end"
245
246#else
247
248 if( hecmw_comm_get_rank() == 0 ) then
249 write(*,*) "##Error : REVOCAP functions are not supported"
250 end if
251 call hecmw_abort( hecmw_comm_get_comm() )
252#endif
253 end subroutine fstr_rcap_get
254 !------------------------------------------------------------------------------
255 subroutine fstr_get_convergence( revocap_flag )
256 implicit none
257 integer(kind=kint) :: revocap_flag
258
259#ifdef WITH_REVOCAP
260 write(idbg,*) "fstr_get_convergence: start"
261
262 call rcapf_get_convergence( revocap_flag )
263
264 write(idbg,*) "fstr_get_convergence: end"
265
266#else
267
268 if( hecmw_comm_get_rank() == 0 ) then
269 write(*,*) "##Error : REVOCAP functions are not supported"
270 end if
271 call hecmw_abort( hecmw_comm_get_comm() )
272#endif
273 end subroutine fstr_get_convergence
274 !------------------------------------------------------------------------------
275
276end module m_fstr_rcap_io
subroutine, public fstr_rcap_send(fstrcpl)
subroutine, public fstr_rcap_get(fstrcpl)
subroutine, public fstr_rcap_finalize(fstrparam, fstrcpl)
subroutine, public fstr_rcap_initialize(hecmesh, fstrparam, fstrcpl)
subroutine fstr_get_convergence(revocap_flag)
This module defined coomon data and basic structures for analysis.
Definition: m_fstr.f90:15
integer(kind=kint) myrank
PARALLEL EXECUTION.
Definition: m_fstr.f90:80
integer(kind=kint), parameter idbg
Definition: m_fstr.f90:95
Data for coupling analysis.
Definition: m_fstr.f90:580
FSTR INNER CONTROL PARAMETERS (fstrPARAM)
Definition: m_fstr.f90:138