SPEED
WRITE_VTU_TIMEDATA.f90
Go to the documentation of this file.
1! Copyright (C) 2012 The SPEED FOUNDATION
2! Author: Ilario Mazzieri
3!
4! This file is part of SPEED.
5!
6! SPEED is free software; you can redistribute it and/or modify it
7! under the terms of the GNU Affero General Public License as
8! published by the Free Software Foundation, either version 3 of the
9! License, or (at your option) any later version.
10!
11! SPEED is distributed in the hope that it will be useful, but
12! WITHOUT ANY WARRANTY; without even the implied warranty of
13! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14! Affero General Public License for more details.
15!
16! You should have received a copy of the GNU Affero General Public License
17! along with SPEED. If not, see <http://www.gnu.org/licenses/>.
18
19
24
32 subroutine write_vtu_timedata(nn_loc, cs_nnz_loc, cs_loc, &
33 nmat, sdeg, prop_mat, tag_mat, &
34 xx_loc, yy_loc, zz_loc, mpi_id, nmpi, &
35 nn, vtk_numbering_map, &
36 its, strain, stress, disp )
37
38 implicit none
39
40 character*70 :: file_name_vtu, temp_char, file_name_pvtu, num_str1, num_str2
41 character*200 :: buffer
42
43 integer*4 :: nel_loc, nn_loc, cs_nnz_loc, mpi_id, nmat, nn, its, nmpi
44 integer*4 :: ie, inode, im, i
45 integer :: i1
46 integer*4 :: unit_mpi, pvtu_unit
47 integer*4, dimension(nmat) :: sdeg, tag_mat
48 integer*4, dimension(nn*nn*nn) :: vtk_numbering_map, node_numbering_vtkwrite, loc_nod_indx
49 integer*4, dimension(0:cs_nnz_loc) :: cs_loc
50
51 real*8, dimension(nmat,4) :: prop_mat
52 real*8, dimension(nn_loc) :: xx_loc, yy_loc, zz_loc
53 real*8, dimension(6*nn_loc) :: strain, stress
54 real*8, dimension(3*nn_loc) :: disp
55
56
57
58 nel_loc = cs_loc(0) - 1
59
60 !Coverting strain/stress to satisfy real*4 (Float32 in vtu) range
61 do inode=1,6*nn_loc
62 if (strain(inode).lt. 1.0e-29) strain(inode) = 1.0e-29;
63 if (strain(inode).gt. 1.0e+29) strain(inode) = 1.0e+29;
64 enddo
65
66
67
68 write(file_name_vtu,'(a,i8.8,a,i5.5,a)') './VTKOUT/SNAPSHOT_',its,'.',mpi_id,'.vtu'
69 !-------------- VTU FILE -------------------------------------------------------------------------------------
70 write(temp_char,*)'(i12,',nn*nn*nn,'i12)'
71
72 unit_mpi = 2500 + mpi_id
73 open(unit_mpi,file=file_name_vtu)
74
75 ! write header
76 buffer='<?xml version="1.0" encoding="utf-8"?>'
77 write(unit_mpi,'(a)')trim(buffer)
78 buffer='<VTKFile type="UnstructuredGrid" version="0.1" byte_order="LittleEndian">';
79 write(unit_mpi,'(a)')trim(buffer)
80 buffer='<UnstructuredGrid>'
81 write(unit_mpi,'(a)')trim(buffer)
82 write(num_str1,*)nn_loc; write(num_str2,*)nel_loc;
83 buffer='<Piece NumberOfPoints="'//trim(adjustl(num_str1))//'" NumberOfCells="'//trim(adjustl(num_str2))//'">'
84 write(unit_mpi,'(a)')trim(buffer)
85
86 ! Write POINTS Info (Node Coordinates)--------------------------------
87 buffer='<Points>'
88 write(unit_mpi,'(a)')trim(buffer)
89 buffer='<DataArray type="Float32" NumberOfComponents="3" format="ascii">'
90 write(unit_mpi,'(a)')trim(buffer)
91 do inode=1,nn_loc
92 write(unit_mpi,'(3e20.12)') xx_loc(inode),yy_loc(inode),zz_loc(inode)
93 enddo
94 buffer='</DataArray>'
95 write(unit_mpi,'(a)')trim(buffer)
96 buffer='</Points>'
97 write(unit_mpi,'(a)')trim(buffer)
98 !--------------------------------------------------------------------
99
100 ! Write CELL info-------------------------------------------------------
101 buffer='<Cells>'
102 write(unit_mpi,'(a)')trim(buffer)
103
104 buffer='<DataArray type="Int32" Name="connectivity" format="ascii">'
105 write(unit_mpi,'(a)')trim(buffer)
106 do ie=1,nel_loc
107 im = cs_loc(cs_loc(ie -1) + 0 )
108 ! nn = sdeg(im) +1
109 do i=1,(nn*nn*nn)
110 loc_nod_indx(i) = cs_loc(cs_loc(ie -1) + i) - 1
111 enddo
112 do i=1,(nn*nn*nn)
113 node_numbering_vtkwrite(i) = loc_nod_indx(vtk_numbering_map(i))
114 enddo
115 write(unit_mpi,temp_char) (node_numbering_vtkwrite(i), i=1,(nn*nn*nn))
116 enddo
117 buffer='</DataArray>'
118 write(unit_mpi,'(a)')trim(buffer)
119
120 buffer='<DataArray type="Int32" Name="offsets" format="ascii">'
121 write(unit_mpi,'(a)')trim(buffer)
122 write(unit_mpi,'(6i12)') (ie*(nn*nn*nn), ie=1,nel_loc)
123 buffer='</DataArray>'
124 write(unit_mpi,'(a)')trim(buffer)
125
126 buffer='<DataArray type="Int32" Name="types" format="ascii">'
127 write(unit_mpi,'(a)')trim(buffer)
128 ! vtkCellType: hexahedrons (ID = 12 for linear 8 noded hex)
129 ! ID = 29 for HEX27 (triquadratic Hex)
130 ! 72 for lagrange-higher order hex
131 ! Reference : https://vtk.org/doc/nightly/html/vtkCellType_8h_source.html
132 ! https://visit-sphinx-github-user-manual.readthedocs.io/en/v3.3.0/data_into_visit/VTKFormat.html
133 write(unit_mpi,'(6i12)') (72,ie=1,nel_loc)
134 buffer='</DataArray>'
135 write(unit_mpi,'(a)')trim(buffer)
136
137 buffer='</Cells>'
138 write(unit_mpi,'(a)')trim(buffer)
139 !--------------------------------------------------------------------
140
141 ! Write POINTS data --------------------------------------------
142 buffer='<PointData>'
143 write(unit_mpi,'(a)')trim(buffer)
144
145 !Strain
146 buffer='<DataArray type="Float32" NumberOfComponents="6" Name="Strain" format="ascii">'
147 write(unit_mpi,'(a)')trim(buffer)
148 do inode = 1,nn_loc
149 im = 6*(inode-1)
150 write(unit_mpi,'(6e20.12)') strain(im+1), strain(im+2), strain(im+3), strain(im+4), strain(im+5), strain(im+6)
151 enddo
152 buffer='</DataArray>'
153 write(unit_mpi,'(a)')trim(buffer)
154
155 !displacement
156 buffer='<DataArray type="Float32" NumberOfComponents="3" Name="Disp" format="ascii">'
157 write(unit_mpi,'(a)')trim(buffer)
158 do inode = 1,nn_loc
159 im = 3*(inode-1)
160 write(unit_mpi,'(3e20.12)') disp(im+1), disp(im+2), disp(im+3)
161 enddo
162 buffer='</DataArray>'
163 write(unit_mpi,'(a)')trim(buffer)
164
165
166 buffer='</PointData>'
167 write(unit_mpi,'(a)')trim(buffer)
168 !--------------------------------------------------------------------
169
170
171 ! Write CELLS data --------------------------------------------
172 ! buffer='<CellData>'
173 ! write(unit_mpi,'(a)')trim(buffer)
174 ! buffer='</CellData>'
175 ! write(unit_mpi,'(a)')trim(buffer)
176 !--------------------------------------------------------------------
177
178 buffer='</Piece>'
179 write(unit_mpi,'(a)')trim(buffer)
180 buffer='</UnstructuredGrid>'
181 write(unit_mpi,'(a)')trim(buffer)
182 buffer='</VTKFile>'
183 write(unit_mpi,'(a)')trim(buffer)
184 close(unit_mpi)
185
186
187 close(unit_mpi)
188 !-----------------------------------------------------------------------------------------------------------------
189
190
191
192
193
194 !-------------- PVTU FILE -------------------------------------------------------------------------------------
195 if (mpi_id.eq.0) then
196 write(file_name_pvtu,'(a,i8.8,a)') './VTKOUT/SNAPSHOT_',its,'.pvtu'
197 pvtu_unit = 2499
198 open(pvtu_unit,file=file_name_pvtu)
199
200 buffer='<?xml version="1.0" encoding="utf-8"?>'
201 write(pvtu_unit,'(a)')trim(buffer)
202 buffer='<VTKFile type="PUnstructuredGrid" version="0.1" byte_order="LittleEndian">';
203 write(pvtu_unit,'(a)')trim(buffer)
204 buffer='<PUnstructuredGrid GhostLevel="0">'
205 write(pvtu_unit,'(a)')trim(buffer)
206
207
208 buffer='<PPoints>'
209 write(pvtu_unit,'(a)')trim(buffer)
210 buffer='<PDataArray type="Float32" NumberOfComponents="3" format="ascii"/>'
211 write(pvtu_unit,'(a)')trim(buffer)
212 buffer='</PPoints>'
213 write(pvtu_unit,'(a)')trim(buffer)
214
215
216 buffer='<PCells>'
217 write(pvtu_unit,'(a)')trim(buffer)
218 buffer='<PDataArray type="Int32" Name="connectivity" NumberOfComponents="1" format="ascii"/>'
219 write(pvtu_unit,'(a)')trim(buffer)
220 buffer='<PDataArray type="Int32" Name="offsets" NumberOfComponents="1" format="ascii"/>'
221 write(pvtu_unit,'(a)')trim(buffer)
222 buffer='<PDataArray type="Int32" Name="types" NumberOfComponents="1" format="ascii"/>'
223 write(pvtu_unit,'(a)')trim(buffer)
224 buffer='</PCells>'
225 write(pvtu_unit,'(a)')trim(buffer)
226
227
228 buffer='<PPointData>'
229 write(pvtu_unit,'(a)')trim(buffer)
230 buffer='<PDataArray type="Float32" NumberOfComponents="6" Name="Strain" format="ascii"/>'
231 write(pvtu_unit,'(a)')trim(buffer)
232 buffer='<PDataArray type="Float32" NumberOfComponents="3" Name="Disp" format="ascii"/>'
233 write(pvtu_unit,'(a)')trim(buffer)
234 buffer='</PPointData>'
235 write(pvtu_unit,'(a)')trim(buffer)
236
237
238 ! buffer='<PCellData>'
239 ! write(pvtu_unit,'(a)')trim(buffer)
240 ! buffer='</PCellData>'
241 ! write(pvtu_unit,'(a)')trim(buffer)
242
243 do i=0,(nmpi-1)
244 write(file_name_vtu,'(a,i8.8,a,i5.5,a)') './SNAPSHOT_',its,'.',i,'.vtu'
245 buffer='<Piece Source="'//trim(file_name_vtu)//'"/>'
246 write(pvtu_unit,'(a)')trim(buffer)
247 enddo
248
249 buffer='</PUnstructuredGrid>'
250 write(pvtu_unit,'(a)')trim(buffer)
251 buffer='</VTKFile>'
252 write(pvtu_unit,'(a)')trim(buffer)
253 close(pvtu_unit);
254 endif
255
256 end subroutine write_vtu_timedata
257
258 subroutine write_pvd_timedata(nmpi, nts, ncount, dt)
259
260 integer*4, intent(in) :: nmpi, nts, ncount
261 integer*4 :: unit_mpi
262
263 real*8 :: dt
264
265 character*70 :: filename_pvtu, filename_pvd, str_time
266 character*200 :: buffer
267
268 write(filename_pvd,'(a,i8.8,a)') './VTKOUT/SNAPSHOT.pvd'
269 unit_mpi = 2499
270 open(unit_mpi,file=filename_pvd)
271
272 buffer='<?xml version="1.0" encoding="utf-8"?>'
273 write(unit_mpi,'(a)')trim(buffer)
274 buffer='<VTKFile type="Collection" version="0.1" byte_order="LittleEndian">'
275 write(unit_mpi,'(a)')trim(buffer)
276 buffer='<Collection>'
277 write(unit_mpi,'(a)')trim(buffer)
278
279 do i=0,nts,ncount
280 write(filename_pvtu,'(a,i8.8,a)') './SNAPSHOT_',i,'.pvtu'
281 write(str_time,'(f16.6)') i*dt
282 buffer='<DataSet timestep="'//trim(adjustl(str_time))//'" part="001" file="'//trim(filename_pvtu)//'"/>'
283 write(unit_mpi,'(a)')trim(buffer)
284 enddo
285
286 buffer='</Collection>'
287 write(unit_mpi,'(a)')trim(buffer)
288 buffer='</VTKFile>'
289 write(unit_mpi,'(a)')trim(buffer)
290 close(unit_mpi);
291
292 end subroutine WRITE_PVD_TIMEDATA
write_pvd_timedata
subroutine write_pvd_timedata(nmpi, nts, ncount, dt)
Definition WRITE_VTU_TIMEDATA.f90:259
write_vtu_timedata
subroutine write_vtu_timedata(nn_loc, cs_nnz_loc, cs_loc, nmat, sdeg, prop_mat, tag_mat, xx_loc, yy_loc, zz_loc, mpi_id, nmpi, nn, vtk_numbering_map, its, strain, stress, disp)
...Writing VTK file to visualise in Paraview
Definition WRITE_VTU_TIMEDATA.f90:37