OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
stat_inimap1d_spmd.F File Reference
#include "implicit_f.inc"
#include "com01_c.inc"
#include "com04_c.inc"
#include "com08_c.inc"
#include "param_c.inc"
#include "scr03_c.inc"
#include "scr17_c.inc"
#include "task_c.inc"
#include "units_c.inc"
#include "chara_c.inc"

Go to the source code of this file.

Functions/Subroutines

subroutine stat_inimap1d_spmd (x, v, itab, ipart_state, nodtag, ipart, iparts, ipartq, iparttg, mat_param, igeo, iparg, ixs, ixq, ixtg, elbuf_tab, multi_fvm, bufmat, ipm)

Function/Subroutine Documentation

◆ stat_inimap1d_spmd()

subroutine stat_inimap1d_spmd ( dimension(3,numnod), intent(in) x,
dimension(3,numnod), intent(in) v,
integer, dimension(numnod), intent(in) itab,
integer, dimension(npart), intent(in) ipart_state,
integer, dimension(numnod), intent(inout) nodtag,
integer, dimension(lipart1,*), intent(in) ipart,
integer, dimension(*), target iparts,
integer, dimension(*), target ipartq,
integer, dimension(*), target iparttg,
type(matparam_struct_), dimension(nummat), intent(in) mat_param,
integer, dimension(npropgi,numgeo), intent(in) igeo,
integer, dimension(nparg,*), intent(in) iparg,
integer, dimension(nixs,numels), intent(in) ixs,
integer, dimension(nixq,numelq), intent(in) ixq,
integer, dimension(nixtg,numeltg), intent(in) ixtg,
type (elbuf_struct_), dimension(ngroup), intent(in), target elbuf_tab,
type(multi_fvm_struct), intent(in) multi_fvm,
dimension(*), target bufmat,
integer, dimension(npropmi,*), intent(inout) ipm )

Definition at line 40 of file stat_inimap1d_spmd.F.

45C-----------------------------------------------
46C Description
47C-----------------------------------------------
48C This subroutine is generating mapping data to be used with /INIMAP1D option.
49C include file is incremented starting from ROOT_1D_0001.inimap
50C It contains 1D fonctions for submaterial data (volume fraction, mass density, energy density)
51C and also function for global velocity
52C User can use the generated file in a second run using #include command in the new Starter input file (target mesh)
53C-----------------------------------------------
54C M o d u l e s
55C-----------------------------------------------
56 USE message_mod
57 USE elbufdef_mod
59 USE multi_fvm_mod
60 USE multimat_param_mod , ONLY : m51_n0phas, m51_nvphas
61 USE matparam_def_mod, ONLY : matparam_struct_
62C-----------------------------------------------
63C I m p l i c i t T y p e s
64C-----------------------------------------------
65#include "implicit_f.inc"
66C-----------------------------------------------
67C C o m m o n B l o c k s
68C-----------------------------------------------
69#include "com01_c.inc"
70#include "com04_c.inc"
71#include "com08_c.inc"
72#include "param_c.inc"
73#include "scr03_c.inc"
74#include "scr17_c.inc"
75#include "task_c.inc"
76#include "units_c.inc"
77#include "chara_c.inc"
78C-----------------------------------------------
79C D u m m y A r g u m e n t s
80C-----------------------------------------------
81 INTEGER,INTENT(IN) :: ITAB(NUMNOD), IPART(LIPART1,*),IGEO(NPROPGI,NUMGEO),IXS(NIXS,NUMELS), IPART_STATE(NPART)
82 INTEGER,INTENT(IN) :: IPARG(NPARG,*),IXQ(NIXQ,NUMELQ),IXTG(NIXTG,NUMELTG)
83 INTEGER, INTENT(INOUT) :: NODTAG(NUMNOD),IPM(NPROPMI,*)
84 INTEGER, TARGET :: IPARTS(*),IPARTQ(*),IPARTTG(*)
85 my_real,INTENT(IN) :: x(3,numnod),v(3,numnod)
86 my_real, INTENT(IN), TARGET :: bufmat(*)
87 TYPE (ELBUF_STRUCT_), DIMENSION(NGROUP), TARGET, INTENT(IN) :: ELBUF_TAB
88 TYPE(MULTI_FVM_STRUCT), INTENT(IN) :: MULTI_FVM
89 TYPE(MATPARAM_STRUCT_) ,DIMENSION(NUMMAT) ,INTENT(IN) :: MAT_PARAM
90C-----------------------------------------------
91C L o c a l V a r i a b l e s
92C-----------------------------------------------
93 INTEGER I, N, JJ,J, IPRT0, IPRT, K, STAT_NUMELS_1, KK, INOD
94 INTEGER NG, NEL, NFT, LFT, LLT, ITY, ISOLNOD, IOFF,NPT
95 INTEGER NUM_CENTROIDS, IPOS,MLW,IFORM,NBMAT,NB2,ISUBMAT,NNOD,NNOD2
96 INTEGER NUVAR,NUM_CELL
97 TYPE(G_BUFEL_) ,POINTER :: GBUF
98 my_real p0(3),p0_inf(3),p0_sup(3),shift_c,shift_n,length
99 my_real max_xc,max_yc,max_zc,min_xc,min_yc,min_zc
100 my_real lx,ly,lz
101 my_real dx,dy,dz
102 my_real dotprod,tol,xyz(3),vel(3),dist,vect(3)
103 CHARACTER FILNAM*100, CHSTAT*4
104 INTEGER XMIN_CELL_ID,YMIN_CELL_ID,ZMIN_CELL_ID,XMAX_CELL_ID,YMAX_CELL_ID,ZMAX_CELL_ID
105 INTEGER IS_ITY_1, IS_ITY_2, IS_ITY_7, LAST_CELL,FIRST_CELL,IMAT,NPAR,IADBUF
106 INTEGER, ALLOCATABLE, DIMENSION(:) :: IDX
107 INTEGER, POINTER,DIMENSION(:) :: IPART_PTR
108 my_real, POINTER,DIMENSION(:) :: uparam
109 TYPE(BUF_MAT_) ,POINTER :: MBUF
110 my_real, ALLOCATABLE, DIMENSION(:,:) :: map_nodes
111 INTEGER, ALLOCATABLE, DIMENSION(:,:) :: GET_CELL_FOM_CENTROID !(NG,I+NFT)
112 INTEGER NPTS(NSPMD),NCELL(NSPMD),NPTS_TOT,NCELL_TOT
113 my_real, DIMENSION(:,:), ALLOCATABLE :: work
114 INTEGER, DIMENSION(:), ALLOCATABLE :: WORK_INDX
115 my_real :: len_(nspmd),len_tot,shift_c_min,shift_n_min
116 INTEGER :: IDX1(21),IDX2(21),IDX3(21)
117C-----------------------------------------------
118C S o u r c e L i n e s
119C-----------------------------------------------
120
121 !---INITIALIZATION
122 !
123 state_inimap_call_number = state_inimap_call_number +1 !number of written files ROOT_INIMAP_00**.rad
124 num_centroids = 0
125 mlw=0
126 !box containing user domain :
127 min_xc = ep20
128 min_yc = ep20
129 min_zc = ep20
130 max_xc = -ep20
131 max_yc = -ep20
132 max_zc = -ep20
133 !detecting elem types to manager error messages :
134 is_ity_1 = 0
135 is_ity_2 = 0
136 is_ity_7 = 0
137
138 !---ALLOCATIONS
139 !
140 IF(.NOT.(ALLOCATED(state_inimap_buf))) THEN
141 IF(ispmd/=0)THEN
142 ALLOCATE(state_inimap_buf(1))
143 ELSE
144 ALLOCATE(state_inimap_buf(nspmd)) !process 0 will gather all data
145 ENDIF
146 ENDIF
147
148 !---ENUMARATION : ELEM TYPES AND BOX DIMENSION
149 !
150 DO ng=1,ngroup
151 ity =iparg(5,ng)
152 isolnod = iparg(28,ng)
153 nel =iparg(2,ng)
154 nft =iparg(3,ng)
155 gbuf => elbuf_tab(ng)%GBUF
156 mlw = iparg(1,ng)
157 lft=1
158 llt=nel
159 npt=0
160 IF(ity == 1) THEN
161 !---bricks
162 is_ity_1=1
163 npt=isolnod
164 ipart_ptr => iparts(1:numels)
165 ELSEIF(ity == 2)THEN
166 !---quads
167 is_ity_2=1
168 npt=4
169 ipart_ptr => ipartq(1:numelq)
170 ELSEIF(ity == 7 .AND. n2d /= 0)THEN
171 !---Triangles
172 is_ity_7=1
173 npt=3
174 ipart_ptr => iparttg(1:numeltg)
175 ENDIF
176 IF(npt /= 0)THEN
177 DO i=lft,llt
178 n = i + nft
179 iprt=ipart_ptr(n)
180 imat =ipart(1,iprt)
181 IF(ipart_state(iprt)==0)cycle
182 num_centroids = num_centroids +1
183 IF(is_ity_1==1)THEN
184 p0(1) = sum( x(1,ixs(2:9,n)) ) / npt
185 p0(2) = sum( x(2,ixs(2:9,n)) ) / npt
186 p0(3) = sum( x(3,ixs(2:9,n)) ) / npt
187 ELSEIF(is_ity_2==1)THEN
188 p0(1) = sum( x(1,ixq(2:5,n)) ) / npt
189 p0(2) = sum( x(2,ixq(2:5,n)) ) / npt
190 p0(3) = sum( x(3,ixq(2:5,n)) ) / npt
191 ENDIF
192 IF(min_xc>p0(1))THEN
193 min_xc=p0(1)
194 xmin_cell_id = n
195 ENDIF
196 IF(min_yc>p0(2))THEN
197 min_yc=p0(2)
198 ymin_cell_id = n
199 ENDIF
200 IF(min_zc>p0(3))THEN
201 min_zc=p0(3)
202 zmin_cell_id = n
203 ENDIF
204 IF(max_xc<p0(1))THEN
205 max_xc=p0(1)
206 xmax_cell_id = n
207 ENDIF
208 IF(max_yc<p0(2))THEN
209 max_yc=p0(2)
210 ymax_cell_id = n
211 ENDIF
212 IF(max_zc<p0(3))THEN
213 max_zc=p0(3)
214 zmax_cell_id = n
215 ENDIF
216 END DO
217 ELSE
218 !no related cells: bricks, quads, and triangles only
219 END IF
220 END do! next NG
221
222 !---NUMBERING :
223 ! +--CENTROIDS
224 state_inimap_buf(1)%NUM_CENTROIDS = num_centroids
225 state_inimap_buf(1)%NUM_POINTS = 0
226 ! +--FIRST ESTIMATION OF NODE NUMBER (DUPLICATED PROJECTIONS AT THIS STEP)
227 nnod=0
228 DO i=1,numnod
229 IF(nodtag(i) == 1)THEN
230 nnod=nnod+1
231 ENDIF
232 ENDDO
233
234 !---ANOTHER ALLOCATIONS
235 !
236 IF(.NOT.ALLOCATED(map_nodes))ALLOCATE(map_nodes(3,nnod))
237 IF(.NOT.ALLOCATED(state_inimap_buf(1)%POS_CENTROIDS))ALLOCATE(state_inimap_buf(1)%POS_CENTROIDS(num_centroids))
238 IF(.NOT.ALLOCATED(get_cell_fom_centroid))THEN
239 ALLOCATE(get_cell_fom_centroid(2,num_centroids))
240 ENDIF
241
242 lx=zero
243 ly=zero
244 lz=zero
245 IF(num_centroids > 0)THEN
246 !---BOX CONTAINING USER-1D-DOMAIN & BUILDING 1D AXIS
247 ! use centroids otherwise with spmd decomposition 1d portion may be smaller than orthogonal directions. (S >> dl)
248 lx=max_xc-min_xc
249 ly=max_yc-min_yc
250 lz=max_zc-min_zc
251 vect(1:3)=(/lx,ly,lz/)
252 !---first and last cell along this user 1D domain
253 ipos = maxloc(vect(1:3),1)
254 SELECT CASE(ipos)
255 CASE(1)
256 first_cell = xmin_cell_id
257 last_cell = xmax_cell_id
258 CASE(2)
259 first_cell = ymin_cell_id
260 last_cell = ymax_cell_id
261 CASE(3)
262 first_cell = zmin_cell_id
263 last_cell = zmax_cell_id
264 END SELECT
265 !---first and last centroid position
266 IF(is_ity_1==1)THEN
267 p0_inf(1) = sum( x(1,ixs(2:9,first_cell)) ) / npt
268 p0_inf(2) = sum( x(2,ixs(2:9,first_cell)) ) / npt
269 p0_inf(3) = sum( x(3,ixs(2:9,first_cell)) ) / npt
270 p0_sup(1) = sum( x(1,ixs(2:9,last_cell)) ) / npt
271 p0_sup(2) = sum( x(2,ixs(2:9,last_cell)) ) / npt
272 p0_sup(3) = sum( x(3,ixs(2:9,last_cell)) ) / npt
273 ELSEIF(is_ity_2==1)THEN
274 p0_inf(1) = sum( x(1,ixq(2:5,first_cell)) ) / npt
275 p0_inf(2) = sum( x(2,ixq(2:5,first_cell)) ) / npt
276 p0_inf(3) = sum( x(3,ixq(2:5,first_cell)) ) / npt
277 p0_sup(1) = sum( x(1,ixq(2:5,last_cell)) ) / npt
278 p0_sup(2) = sum( x(2,ixq(2:5,last_cell)) ) / npt
279 p0_sup(3) = sum( x(3,ixq(2:5,last_cell)) ) / npt
280 ENDIF
281 !---first and last centroid position are determining 1D axis :
282 vect(1:3)=(/p0_sup(1)-p0_inf(1),p0_sup(2)-p0_inf(2),p0_sup(3)-p0_inf(3)/)
283 lx=vect(1)
284 ly=vect(2)
285 lz=vect(3)
286 length = sqrt(vect(1)*vect(1) + vect(2)*vect(2) + vect(3)*vect(3)) !1d axis length (needed to normalize following dot products)
287 shift_c=zero
288 IF(length > zero)shift_c = (p0_inf(1)*lx + p0_inf(2)*ly + p0_inf(3)*lz) / length !distance from origin at first centroid position
289 state_inimap_buf(1)%SHIFT_Cy = shift_c
290 state_inimap_buf(1)%SHIFT_Cz = zero
291 state_inimap_buf(1)%LENGTH = length
292 ELSE
293 state_inimap_buf(1)%SHIFT_Cy = zero
294 state_inimap_buf(1)%SHIFT_Cz = zero
295 state_inimap_buf(1)%LENGTH = zero
296 ENDIF
297 !---ERROR MESSAGES
298 !
299 IF(is_ity_7 > 0)THEN
300 CALL ancmsg(msgid=284,anmode=aninfo,c1=" -- 1D DOMAIN MUST BE MESHED WITH BRICKS OR QUADS ONLY")
301 return
302 ENDIF
303
304 !---STARTING POINT OF FIRST NODE IN GENERAL FRAME (SHIFT_N)
305 ! SHIFT_N : nodal posisiton (first node along axis)
306 ! SHIFT_C : cell centroid position (first centroid along axis)
307 ! by nature : SHIFT_N < SHIFT_C
308 !
309 IF(is_ity_1==1)THEN
310 dx = x(1,ixs(2,first_cell))
311 dy = x(2,ixs(2,first_cell))
312 dz = x(3,ixs(2,first_cell))
313 DO jj=3,npt
314 IF(x(1,ixs(jj,first_cell)) < dx)dx=x(1,ixs(jj,first_cell))
315 IF(x(2,ixs(jj,first_cell)) < dy)dy=x(2,ixs(jj,first_cell))
316 IF(x(3,ixs(jj,first_cell)) < dz)dz=x(3,ixs(jj,first_cell))
317 ENDDO
318 ELSEIF(is_ity_2==1)THEN
319 dx = x(1,ixq(2,first_cell))
320 dy = x(2,ixq(2,first_cell))
321 dz = x(3,ixq(2,first_cell))
322 DO jj=3,npt
323 IF(x(1,ixq(jj,first_cell)) < dx)dx=x(1,ixq(jj,first_cell))
324 IF(x(2,ixq(jj,first_cell)) < dx)dy=x(2,ixq(jj,first_cell))
325 IF(x(3,ixq(jj,first_cell)) < dx)dz=x(3,ixq(jj,first_cell))
326 ENDDO
327 ENDIF
328 !first point projection on 1d axis
329 shift_n = zero
330 IF(length > zero)shift_n=(dx*lx + dy*ly + dz*lz) / length
331 state_inimap_buf(1)%SHIFT_Ny = shift_n
332 state_inimap_buf(1)%SHIFT_Nz = zero
333
334 !---ABSCISSA : CENTROIDS POSITION LIST ( STATE_INIMAP_BUF(1)%POS_CENTROIDS(1:NUM_CENTROIDS) )
335 !
336 k=1
337 DO ng=1,ngroup
338 ity =iparg(5,ng)
339 isolnod = iparg(28,ng)
340 nel =iparg(2,ng)
341 nft =iparg(3,ng)
342 gbuf => elbuf_tab(ng)%GBUF
343 mlw = iparg(1,ng)
344 lft=1
345 llt=nel
346 IF(npt /= 0)THEN
347 DO i=lft,llt
348 n = i + nft
349 iprt=ipart_ptr(n)
350 IF(ipart_state(iprt)==0)cycle
351 IF(is_ity_1==1)THEN
352 p0(1) = sum( x(1,ixs(2:9,n)) ) / npt
353 p0(2) = sum( x(2,ixs(2:9,n)) ) / npt
354 p0(3) = sum( x(3,ixs(2:9,n)) ) / npt
355 ELSEIF(is_ity_2==1)THEN
356 p0(1) = sum( x(1,ixq(2:5,n)) ) / npt
357 p0(2) = sum( x(2,ixq(2:5,n)) ) / npt
358 p0(3) = sum( x(3,ixq(2:5,n)) ) / npt
359 ENDIF
360 dx = p0(1)-p0_inf(1)
361 dy = p0(2)-p0_inf(2)
362 dz = p0(3)-p0_inf(3)
363 !DIST = SQRT(DX*DX+DY*DY+DZ*DZ)
364 dotprod = zero
365 IF(length > zero)dotprod = (lx*dx + ly*dy + lz*dz) / length
366 state_inimap_buf(1)%POS_CENTROIDS(k) = dotprod + state_inimap_buf(1)%SHIFT_Cy
367 get_cell_fom_centroid(1,k) = ng
368 get_cell_fom_centroid(2,k) = i
369 k=k+1
370 END DO
371 END IF
372 END do! next NG
373 !---ABSCISSA : CENTROIDS POSITION - SORTING
374 IF(.NOT.ALLOCATED(idx))ALLOCATE(idx(num_centroids))
375 DO k=1,num_centroids ; idx(k)=k; ENDDO
376 IF(num_centroids>0)CALL quicksort(state_inimap_buf(1)%POS_CENTROIDS(:), idx, 1, num_centroids)
377
378 !---CELL DATA TREATMENT
379 ! storing submaterial data : vfrac,rho,E
380 !
381 IF(num_centroids > 0)THEN
382 IF(mlw==151)THEN
383 nbmat = multi_fvm%NBMAT
384 ELSEIF(mlw==51)THEN
385 nbmat = 4
386 ELSE
387 nbmat = 1
388 ENDIF
389 state_inimap_buf(1)%MLW = mlw
390 state_inimap_buf(1)%NSUBMAT = nbmat
391 ALLOCATE(state_inimap_buf(1)%SUBMAT(nbmat))
392 DO i=1,nbmat
393 ALLOCATE(state_inimap_buf(1)%SUBMAT(i)%VFRAC(num_centroids))
394 ALLOCATE(state_inimap_buf(1)%SUBMAT(i)%RHO(num_centroids))
395 ALLOCATE(state_inimap_buf(1)%SUBMAT(i)%E(num_centroids))
396 ALLOCATE(state_inimap_buf(1)%SUBMAT(i)%PRES(num_centroids))
397 ENDDO
398 IF(mlw==151)THEN
399 !velocities
400 state_inimap_buf(1)%NUM_POINTS = num_centroids
401 IF(.NOT.ALLOCATED(state_inimap_buf(1)%POS_NODES))ALLOCATE(state_inimap_buf(1)%POS_NODES(num_centroids))
402 IF(.NOT.ALLOCATED(state_inimap_buf(1)%VEL_NODES))ALLOCATE(state_inimap_buf(1)%VEL_NODES(num_centroids))
403 DO k=1, num_centroids
404 ng = get_cell_fom_centroid(1,idx(k))
405 i = get_cell_fom_centroid(2,idx(k))
406 nft = iparg(3,ng)
407 state_inimap_buf(1)%POS_NODES(k) = state_inimap_buf(1)%POS_CENTROIDS(k)
408 xyz(1:3) = multi_fvm%VEL(1:3,i+nft)
409 dotprod=zero
410 IF(length > zero)dotprod = (lx*xyz(1) + ly*xyz(2) + lz*xyz(3)) / length
411 state_inimap_buf(1)%VEL_NODES(k) = dotprod
412 ENDDO
413 !submat
414 DO isubmat=1,nbmat
415 DO k=1, num_centroids
416 ng = get_cell_fom_centroid(1,idx(k))
417 i = get_cell_fom_centroid(2,idx(k))
418 nft = iparg(3,ng)
419 state_inimap_buf(1)%SUBMAT(isubmat)%VFRAC(k) = multi_fvm%PHASE_ALPHA(isubmat,i+nft)
420 state_inimap_buf(1)%SUBMAT(isubmat)%RHO(k) = multi_fvm%PHASE_RHO(isubmat,i+nft)
421 state_inimap_buf(1)%SUBMAT(isubmat)%E(k) = multi_fvm%PHASE_EINT(isubmat,i+nft)
422 state_inimap_buf(1)%SUBMAT(isubmat)%PRES(k) = multi_fvm%PRES(i+nft)
423 ENDDO
424 ENDDO
425 ELSEIF(mlw==51)THEN
426 nb2=0
427 DO isubmat=1,nbmat
428 DO k=1, num_centroids
429 ng = get_cell_fom_centroid(1,idx(k))
430 i = get_cell_fom_centroid(2,idx(k))
431 nft = iparg(3,ng)
432 nel = iparg(2,ng)
433 n = i + nft
434 iprt=ipart_ptr(n)
435 imat =ipart(1,iprt)
436 IF(mat_param(imat)%MULTIMAT%MID(isubmat) == 0)EXIT !submaterial not defined.
437 nb2=max(nb2,ipm(5,imat))
438 iadbuf = ipm(7,imat)
439 npar = ipm(9,imat)
440 nuvar = ipm(8,imat)
441 uparam => bufmat(iadbuf:iadbuf+npar)
442 kk = m51_n0phas + (uparam(276+isubmat)-1)*m51_nvphas
443 mbuf => elbuf_tab(ng)%BUFLY(1)%MAT(1,1,1)
444 state_inimap_buf(1)%SUBMAT(isubmat)%VFRAC(k) = mbuf%VAR(nel*(01+kk-1)+i)
445 state_inimap_buf(1)%SUBMAT(isubmat)%RHO(k) = mbuf%VAR(nel*(12+kk-1)+i)
446 state_inimap_buf(1)%SUBMAT(isubmat)%E(k) = mbuf%VAR(nel*(08+kk-1)+i)
447 state_inimap_buf(1)%SUBMAT(isubmat)%PRES(k) = mbuf%VAR(nel*(18+kk-1)+i)
448 ENDDO
449 ENDDO
450 state_inimap_buf(1)%NSUBMAT = nb2
451 ELSE
452 DO k=1, num_centroids
453 ng = get_cell_fom_centroid(1,idx(k))
454 nel = iparg(2,ng)
455 i = get_cell_fom_centroid(2,idx(k))
456 gbuf => elbuf_tab(ng)%GBUF
457 state_inimap_buf(1)%SUBMAT(1)%VFRAC(k) = 1.d00
458 state_inimap_buf(1)%SUBMAT(1)%RHO(k) = gbuf%RHO(i)
459 state_inimap_buf(1)%SUBMAT(1)%E(k) = gbuf%EINT(i)
460 state_inimap_buf(1)%SUBMAT(1)%PRES(k) = -third * (gbuf%SIG(i) + gbuf%SIG(i + nel) + gbuf%SIG(i + 2 * nel))
461 ENDDO
462 ENDIF
463 ENDIF
464
465 !---VELOCITY TREATMENT FOR STAGGERED SCHEME
466 ! law 51 only
467 !
468 ! MAP_NODES(1,...) : node_id
469 ! MAP_NODES(2,...) : 1d-position
470 ! MAP_NODES(3,...) : vel_1d-value
471 !
472 IF(num_centroids > 0)THEN
473 IF(mlw /= 151)THEN
474 k=1
475 DO i=1,numnod
476 IF(nodtag(i) == 1)THEN
477 map_nodes(1,k)=i
478 xyz(1:3)=x(1:3,i)
479 dotprod=zero
480 IF(length > zero)dotprod = (lx*xyz(1) + ly*xyz(2) + lz*xyz(3)) / length
481 map_nodes(2,k)=dotprod
482 k=k+1
483 ENDIF
484 ENDDO
485 IF(ALLOCATED(idx))DEALLOCATE(idx)
486 ALLOCATE(idx(nnod))
487 DO k=1,nnod ; idx(k)=k; ENDDO
488 CALL quicksort(map_nodes(2,:), idx, 1, nnod)
489 tol=em10*length
490 !---check duplicates and remove
491 nnod2=nnod
492 DO i=2,nnod
493 dist = abs(map_nodes(2,i)-map_nodes(2,i-1))
494 IF(dist <= tol) map_nodes(1,idx(i)) = zero
495 ENDDO
496 k=0
497 DO i=1,nnod
498 IF(map_nodes(1,idx(i)) /= zero)THEN
499 k=k+1
500 ENDIF
501 ENDDO
502 !---VELOCITY : WRITE IN DOMAIN BUFFER (SORTED)
503 IF(.NOT.ALLOCATED(state_inimap_buf(1)%POS_NODES))ALLOCATE(state_inimap_buf(1)%POS_NODES(k))
504 IF(.NOT.ALLOCATED(state_inimap_buf(1)%VEL_NODES))ALLOCATE(state_inimap_buf(1)%VEL_NODES(k))
505 k=0
506 DO i=1,nnod
507 IF(map_nodes(1,idx(i)) /= zero)THEN
508 k=k+1
509 state_inimap_buf(1)%POS_NODES(k) = map_nodes(2,i) !already sorted k<=i nothing is erased
510 vel(1:3)=v(1:3,int(map_nodes(1,idx(i))))
511 dotprod=zero
512 IF(length > zero)dotprod = (lx*vel(1) + ly*vel(2) + lz*vel(3)) / length
513 state_inimap_buf(1)%VEL_NODES(k)=dotprod
514 ENDIF
515 ENDDO
516 state_inimap_buf(1)%NUM_POINTS=k
517 ELSE
518
519 ENDIF
520 ENDIF
521
522C-----------------------------------------------
523C S P M D E x c h a n g e
524C-----------------------------------------------
525 IF(nspmd > 1)THEN
528 !
529 IF(ispmd == 0)THEN
530 shift_c_min = state_inimap_buf(1)%SHIFT_Cy
531 shift_n_min = state_inimap_buf(1)%SHIFT_Ny
532 DO i=2,nspmd
533 IF(state_inimap_buf(i)%NUM_CENTROIDS==0)cycle
534 shift_c_min =min(shift_c_min, state_inimap_buf(i)%SHIFT_Cy)
535 shift_n_min =min(shift_n_min, state_inimap_buf(i)%SHIFT_Ny)
536 ENDDO
537 state_inimap_buf(1)%SHIFT_Cy = shift_c_min
538 state_inimap_buf(1)%SHIFT_Ny = shift_n_min
539 ENDIF
540 ENDIF
541
542C-----------------------------------------------
543C S P M D - G a t h e r i n g & S o r t i n g
544C-----------------------------------------------
545 ncell_tot = state_inimap_buf(1)%NUM_CENTROIDS
546 len_tot = state_inimap_buf(1)%LENGTH
547 IF(ispmd == 0 .AND. nspmd > 1)THEN
548 !--cumulated dimensions
549 !
550 npts_tot = 0
551 ncell_tot = 0
552 len_tot = zero
553 DO i=1,nspmd
554 IF(state_inimap_buf(i)%NUM_CENTROIDS==0)THEN
555 npts(i)=0
556 len_(i)=zero
557 ncell(i)=0
558 cycle
559 ENDIF
560 npts(i)= state_inimap_buf(i)%NUM_POINTS ;
561 npts_tot=npts_tot+npts(i)
562 len_(i)=state_inimap_buf(i)%LENGTH ;
563 len_tot=len_tot+len_(i)
564 ncell(i)= state_inimap_buf(i)%NUM_CENTROIDS ;
565 ncell_tot = ncell_tot + ncell(i)
566 ENDDO
567 ALLOCATE(work(npts_tot,3),work_indx(npts_tot))
568 !stat_inimap1d_mp.F
569 !--gathering velocity into working_array
570 !
571 j=0
572 DO i=1,nspmd
573 DO k=1,npts(i)
574 j=j+1
575 work(j,1) = state_inimap_buf(i)%POS_NODES(k)
576 work(j,3) = state_inimap_buf(i)%VEL_NODES(k) !tmp
577 ENDDO
578 ENDDO
579 !
580 !sorting velocity
581 !
582 work_indx(1:npts_tot) = (/(j,j=1,npts_tot)/)
583 CALL quicksort(work(:,1), work_indx, 1, npts_tot)
584 !sort velocity consequently
585 DO i=1,npts_tot
586 work(i,2)=work(work_indx(i),3) !sorted vel
587 ENDDO
588 tol=em10*len_tot
589 work_indx(1:npts_tot) = 0
590 !
591 !---remove duplicates (possible common nodes on adjacent domains)
592 ! STAGGERED SCHEME ONLY
593 !
594 IF(mlw /= 151)THEN
595 DO i=2,npts_tot
596 dist = abs(work(i,1)-work(i-1,1))
597 IF(dist <= tol) THEN
598 work_indx(i) = 1
599 ENDIF
600 ENDDO
601 k=0
602 DO i=1,npts_tot
603 IF(work_indx(i) ==0 )THEN
604 k=k+1
605 work(k,1)=work(i,1) !abscicca
606 work(k,2)=work(i,2) !ordinates
607 ENDIF
608 ENDDO
609 DO i=k+1,npts_tot ; work(i,1:2)=zero ; ENDDO
610 npts_tot=k
611 ENDIF
612 !
613 !---store in relevant buffer (reallocate)
614 !
615 IF(ALLOCATED(state_inimap_buf(1)%VEL_NODES))DEALLOCATE(state_inimap_buf(1)%VEL_NODES)
616 IF(ALLOCATED(state_inimap_buf(1)%POS_NODES))DEALLOCATE(state_inimap_buf(1)%POS_NODES)
617 ALLOCATE(state_inimap_buf(1)%VEL_NODES(npts_tot))
618 ALLOCATE(state_inimap_buf(1)%POS_NODES(npts_tot))
619 state_inimap_buf(1)%NUM_POINTS=npts_tot
620 state_inimap_buf(1)%POS_NODES(1:npts_tot)=work(1:npts_tot,1)
621 state_inimap_buf(1)%VEL_NODES(1:npts_tot)=work(1:npts_tot,2)
622 IF(ALLOCATED(work))DEALLOCATE(work)
623 IF(ALLOCATED(work_indx))DEALLOCATE(work_indx)
624
625
626 nbmat=1
627 DO i=1,nspmd
628 nbmat=max(nbmat,state_inimap_buf(i)%NSUBMAT)
629 ENDDO
630 ALLOCATE(work(ncell_tot,1+4*nbmat))
631 ALLOCATE(work_indx(ncell_tot))
632 !
633 !--gathering submaterial data into working_array (duplicates are not possible with centroids)
634 !
635 j=0
636 DO i=1,nspmd
637 DO k=1,ncell(i)
638 j=j+1
639 work(j,1) = state_inimap_buf(i)%POS_CENTROIDS(k)
640 nbmat = state_inimap_buf(i)%NSUBMAT
641 DO jj=1,nbmat
642 work(j,1+ 4*(jj-1)+1) = state_inimap_buf(i)%SUBMAT(jj)%VFRAC(k)
643 work(j,1+ 4*(jj-1)+2) = state_inimap_buf(i)%SUBMAT(jj)%RHO(k)
644 work(j,1+ 4*(jj-1)+3) = state_inimap_buf(i)%SUBMAT(jj)%E(k)
645 work(j,1+ 4*(jj-1)+4) = state_inimap_buf(i)%SUBMAT(jj)%PRES(k)
646 ENDDO
647 ENDDO
648 ENDDO
649 !
650 !sorting
651 !
652 work_indx(1:ncell_tot) = (/(j,j=1,ncell_tot)/)
653 CALL quicksort(work(:,1), work_indx, 1, ncell_tot)
654 !
655 !---store in relevant buffer (reallocate)
656 !
657 IF(ALLOCATED(state_inimap_buf(1)%POS_CENTROIDS))DEALLOCATE(state_inimap_buf(1)%POS_CENTROIDS)
658 nbmat = state_inimap_buf(1)%NSUBMAT
659 DO jj=1,nbmat
660 IF(ALLOCATED(state_inimap_buf(1)%SUBMAT(jj)%VFRAC))DEALLOCATE(state_inimap_buf(1)%SUBMAT(jj)%VFRAC)
661 IF(ALLOCATED(state_inimap_buf(1)%SUBMAT(jj)%RHO))DEALLOCATE(state_inimap_buf(1)%SUBMAT(jj)%RHO)
662 IF(ALLOCATED(state_inimap_buf(1)%SUBMAT(jj)%E))DEALLOCATE(state_inimap_buf(1)%SUBMAT(jj)%E)
663 IF(ALLOCATED(state_inimap_buf(1)%SUBMAT(jj)%PRES))DEALLOCATE(state_inimap_buf(1)%SUBMAT(jj)%PRES)
664 ENDDO
665 ALLOCATE(state_inimap_buf(1)%POS_CENTROIDS(ncell_tot))
666 DO jj=1,nbmat
667 ALLOCATE(state_inimap_buf(1)%SUBMAT(jj)%VFRAC(ncell_tot))
668 ALLOCATE(state_inimap_buf(1)%SUBMAT(jj)%RHO(ncell_tot))
669 ALLOCATE(state_inimap_buf(1)%SUBMAT(jj)%E(ncell_tot))
670 ALLOCATE(state_inimap_buf(1)%SUBMAT(jj)%PRES(ncell_tot))
671 ENDDO
672 DO j=1,ncell_tot
673 state_inimap_buf(1)%POS_CENTROIDS(j)=work(j,1) !-SHIFT_C_MIN
674 DO jj=1,nbmat
675 state_inimap_buf(1)%SUBMAT(jj)%VFRAC(j)=work(work_indx(j),1+ 4*(jj-1)+1)
676 state_inimap_buf(1)%SUBMAT(jj)%RHO(j) =work(work_indx(j),1+ 4*(jj-1)+2)
677 state_inimap_buf(1)%SUBMAT(jj)%E(j) =work(work_indx(j),1+ 4*(jj-1)+3)
678 state_inimap_buf(1)%SUBMAT(jj)%PRES(j)=work(work_indx(j),1+ 4*(jj-1)+4)
679 ENDDO
680 ENDDO
681 state_inimap_buf(1)%NUM_CENTROIDS = ncell_tot
682 state_inimap_buf(1)%LENGTH = len_tot
683
684
685 endif! IF(ISPMD == 0 .AND. NSPMD > 1)THEN
686
687 IF(ispmd == 0)THEN
688 IF(ncell_tot == 0 .OR. len_tot == zero)THEN
689 print *, "** ERROR WITH /STATE/INIMAP"
690 print *, " -- SITUATION NOT EXPECTED"
691 print *, " -- 1D DOMAIN IS NOT DETECTED."
692 return
693 ENDIF
694 ENDIF
695
696C-----------------------------------------------
697C O u t p u t F i l e
698C-----------------------------------------------
699
700 IF(ispmd == 0)THEN
701 WRITE(chstat,'(I4.4)')state_inimap_call_number
702 filnam=rootnam(1:rootlen)//'_1D_'//chstat//'.inimap'
703 OPEN(unit=220582,file=filnam(1:len(trim(filnam))),access='SEQUENTIAL',form='FORMATTED',status='UNKNOWN')
704 WRITE(unit=220582,fmt='(A,I10)') '#state file for mappgin with /INIMAP1D, iteration = ',state_inimap_call_number
705 WRITE(unit=220582,fmt='(A,A)') '# ROOTNAME = ',rootnam(1:rootlen)
706 WRITE(unit=220582,fmt='(A,I0)') '# VERSION = ',st_invers
707 WRITE(unit=220582,fmt='(A,F20.13)')'# TIME = ',tt
708 WRITE(unit=220582,fmt='(A,I10)') '# NCYCLE = ',ncycle
709 WRITE(unit=220582,fmt='(A,I10)') '# NCELL = ',ncell_tot
710 !WRITE(UNIT=220582,FMT='(A)')'#//SUBMODEL/1'
711 !WRITE(UNIT=220582,FMT='(A)')'#MAPPING DATA (FUNCTIONS)'
712 !WRITE(UNIT=220582,FMT='(A)')"## off_def off_nod off_ele off_part off_mat off_type off_sub"
713 !WRITE(UNIT=220582,FMT='(A1,7I10)')'#',1000, 0, 0, 0, 0, 0, 0
714 ENDIF
715
716 IF(ispmd == 0)THEN
717 !--- OUTPUT FUNCTION FROM CELL DATA BUFFER ---!
718 nbmat = state_inimap_buf(1)%NSUBMAT
719 shift_c = state_inimap_buf(1)%SHIFT_Cy
720 shift_n = state_inimap_buf(1)%SHIFT_Ny
721 num_centroids = state_inimap_buf(1)%NUM_CENTROIDS
722 !---volume fractions
723 ipos=0 !1:1+21
724 DO isubmat = 1,nbmat
725 WRITE(unit=220582,fmt=2001)ipos+isubmat,isubmat
726 DO k=1, num_centroids
727 WRITE(unit=220582,fmt='(2E20.12)') state_inimap_buf(1)%POS_CENTROIDS(k)
728 . ,state_inimap_buf(1)%SUBMAT(isubmat)%VFRAC(k)
729 ENDDO
730 WRITE(unit=220582,fmt=1500)ipos+isubmat,ipos+isubmat,1.00,1.00,-shift_n,0.00
731 ENDDO
732 !---mass densities
733 ipos=100 !31:31+21
734 DO isubmat = 1,nbmat
735 WRITE(unit=220582,fmt=2002)ipos+isubmat,isubmat
736 DO k=1, num_centroids
737 WRITE(unit=220582,fmt='(2E20.12)') state_inimap_buf(1)%POS_CENTROIDS(k)
738 . ,state_inimap_buf(1)%SUBMAT(isubmat)%RHO(k)
739 ENDDO
740 WRITE(unit=220582,fmt=1500)ipos+isubmat,ipos+isubmat,1.00,1.00,-shift_n,0.00
741 ENDDO
742 !---energy density
743 ipos=200 !31:31+21
744 DO isubmat = 1,nbmat
745 WRITE(unit=220582,fmt=2003)ipos+isubmat,isubmat
746 DO k=1, num_centroids
747 WRITE(unit=220582,fmt='(2E20.12)') state_inimap_buf(1)%POS_CENTROIDS(k)
748 . ,state_inimap_buf(1)%SUBMAT(isubmat)%E(k)
749 ENDDO
750 WRITE(unit=220582,fmt=1500)ipos+isubmat,ipos+isubmat,1.00,1.00,-shift_n,0.00
751 ENDDO
752 !---pressure fraction
753 ipos=300 !91:91+21
754 DO isubmat = 1,nbmat
755 WRITE(unit=220582,fmt=2004)ipos+isubmat,isubmat
756 DO k=1, num_centroids
757 WRITE(unit=220582,fmt='(2E20.12)') state_inimap_buf(1)%POS_CENTROIDS(k)
758 . ,state_inimap_buf(1)%SUBMAT(isubmat)%PRES(k)
759 ENDDO
760 WRITE(unit=220582,fmt=1500)ipos+isubmat,ipos+isubmat,1.00,1.00,-shift_n,0.00
761 ENDDO
762 !--- OUTPUT VELOCITY FUNCTION ---!
763 !
764 WRITE(unit=220582,fmt=3000)
765 DO jj=1,state_inimap_buf(1)%NUM_POINTS
766 WRITE(unit=220582,fmt='(2E20.12)')state_inimap_buf(1)%POS_NODES(jj) , state_inimap_buf(1)%VEL_NODES(jj)
767 ENDDO
768 WRITE(unit=220582,fmt=1500)400,400,1.00,1.00,-shift_n,0.00
769 ENDIF
770
771C-----------------------------------------------
772C D e a l l o c a t e & C l o s e
773C-----------------------------------------------
774 IF(ispmd == 0)THEN
775 !---DEALLOCATE
776 IF(ALLOCATED(map_nodes))DEALLOCATE(map_nodes)
777 IF(ALLOCATED(get_cell_fom_centroid))DEALLOCATE(get_cell_fom_centroid)
778 IF(ALLOCATED(idx))DEALLOCATE(idx)
779 DO jj=1,nspmd
780 nbmat = state_inimap_buf(jj)%NSUBMAT
781 IF(ALLOCATED(state_inimap_buf(jj)%SUBMAT))THEN
782 DO i=1,nbmat
783 IF(ALLOCATED(state_inimap_buf(jj)%SUBMAT(i)%VFRAC))DEALLOCATE(state_inimap_buf(jj)%SUBMAT(i)%VFRAC)
784 IF(ALLOCATED(state_inimap_buf(jj)%SUBMAT(i)%RHO))DEALLOCATE(state_inimap_buf(jj)%SUBMAT(i)%RHO)
785 IF(ALLOCATED(state_inimap_buf(jj)%SUBMAT(i)%E))DEALLOCATE(state_inimap_buf(jj)%SUBMAT(i)%E)
786 IF(ALLOCATED(state_inimap_buf(jj)%SUBMAT(i)%PRES))DEALLOCATE(state_inimap_buf(jj)%SUBMAT(i)%PRES)
787 ENDDO
788 ENDIF
789 IF(ALLOCATED(state_inimap_buf(jj)%SUBMAT ))DEALLOCATE(state_inimap_buf(jj)%SUBMAT)
790 IF(ALLOCATED(state_inimap_buf(jj)%POS_NODES ))DEALLOCATE(state_inimap_buf(jj)%POS_NODES)
791 IF(ALLOCATED(state_inimap_buf(jj)%VEL_NODES ))DEALLOCATE(state_inimap_buf(jj)%VEL_NODES)
792 IF(ALLOCATED(state_inimap_buf(jj)%POS_CENTROIDS))DEALLOCATE(state_inimap_buf(jj)%POS_CENTROIDS)
793 IF(ALLOCATED(state_inimap_buf(jj)%SUBMAT))DEALLOCATE(state_inimap_buf(jj)%SUBMAT)
794 ENDDO
795
796
797 !---OUTPUT FILE : FOOTER & CLOSE
798 WRITE(unit=220582,fmt=1000)
799 !WRITE(UNIT=220582,FMT='(A)')'#//ENDSUB'
800
801 idx1=(/1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21/) !vfrac
802 idx2=100+idx1 !rho
803 idx3=300+idx1 !pressure
804 IF(is_stat_inimap_vp)THEN
805 WRITE(unit=220582,fmt='(A)') '#/INIMAP1D/VP/1'
806 ELSE
807 WRITE(unit=220582,fmt='(A)') '#/INIMAP1D/VE/1'
808 ENDIF
809 WRITE(unit=220582,fmt='(A)') '#default input to update from /STATE/INIMAP1D'
810 WRITE(unit=220582,fmt='(A)') '## Type'
811 WRITE(unit=220582,fmt='(A)') '# 1'
812 WRITE(unit=220582,fmt='(A)') '## Node1 Node2'
813 WRITE(unit=220582,fmt='(A)') '# 0 0'
814 WRITE(unit=220582,fmt='(A)') '## Grbric Grquad Grtria'
815 WRITE(unit=220582,fmt='(A)') '# 0 0 0'
816 WRITE(unit=220582,fmt='(A)') '## Fct_v Fscale_v'
817 WRITE(unit=220582,fmt='(a)') '# 400 1.0'
818 DO imat=1,min(21,nbmat)
819 WRITE(unit=220582,fmt='(A)') '## Fct_vf Fct_rho Fscale_rho Fct_p Fscale_p'
820 WRITE(unit=220582,fmt='(A1,I10,2(I10,F20.0))')'#', idx1(imat),idx2(imat),1.0,idx3(imat),1.0
821 ENDDO
822 WRITE(unit=220582,fmt=1000)
823
824 WRITE (iout,500) filnam(1:len(trim(filnam)))
825 WRITE (istdo,500) filnam(1:len(trim(filnam)))
826
827 CLOSE(unit=220582)
828
829 ENDIF
830
831 IF(ALLOCATED(state_inimap_buf))DEALLOCATE(state_inimap_buf)
832
833
834
835
836
837C-----------------------------------------------
838C O u t p u t F o r m a t
839C-----------------------------------------------
840
841 500 FORMAT (4x,' STATE FILE:',1x,a,' WRITTEN')
842
843 1000 FORMAT('#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|')
844
845 1500 FORMAT('#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|',/,
846 . '/MOVE_FUNCT/',i0,/,
847 . 'move_function__',i0,/,
848 . '# ASCALEx FSCALEy ASHIFTx FSHIFTy',/,
849 . 4(6x,e14.7) )
850
851 2001 FORMAT('#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|',/,
852 . '/FUNCT/',i0,/,
853 . 'volume fraction submaterial_',i0,/,
854 . '#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|')
855 2002 FORMAT('#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|',/,
856 . '/FUNCT/',i0,/,
857 . 'mass density submaterial_',i0,/,
858 . '#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|')
859 2003 FORMAT('#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|',/,
860 . '/FUNCT/',i0,/,
861 . 'energy density submaterial_',i0,/,
862 . '#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|')
863 2004 FORMAT('#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|',/,
864 . '/FUNCT/',i0,/,
865 . 'pressure submaterial_',i0,/,
866 . '#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|')
867
868 3000 FORMAT('#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|',/,
869 . '/FUNCT/400',/,
870 . 'velocity_function'/,
871 . '#---1----|----2----|----3----|----4----|----5----|----6----|----7----|----8----|----9----|---10----|')
872C-----------------------------------------------
873 RETURN
#define my_real
Definition cppsort.cpp:32
#define min(a, b)
Definition macros.h:20
#define max(a, b)
Definition macros.h:21
integer state_inimap_call_number
type(map_struct), dimension(:), allocatable state_inimap_buf
recursive subroutine quicksort(a, idx, first, last)
Definition quicksort.F:34
subroutine spmd_state_inimap1d_exch_data()
subroutine spmd_state_inimap_exch_siz()
subroutine ancmsg(msgid, msgtype, anmode, i1, i2, i3, i4, i5, i6, i7, i8, i9, i10, i11, i12, i13, i14, i15, i16, i17, i18, i19, i20, r1, r2, r3, r4, r5, r6, r7, r8, r9, c1, c2, c3, c4, c5, c6, c7, c8, c9, prmode)
Definition message.F:889