SPEED
WRITE_VTK_MECH_PROP.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! propname = exaple 'MASS', 'FORCEZ', 'DAMPMATRIX'
33 subroutine write_vtk_mech_prop(nn_loc, cs_nnz_loc, cs_loc, &
34 nmat, sdeg, prop_mat, tag_mat, &
35 nmat_nlp, tag_mat_nlp, &
36 xx_loc, yy_loc, zz_loc, mpi_id, nn, vtk_numbering_map)
37
38 implicit none
39
40 character*70 :: file_name, prop_name, temp_char
41
42 integer*4 :: nel_loc, nn_loc, cs_nnz_loc, mpi_id, nmat, nmat_nlp, nn
43 integer*4 :: ie, inode, im, im_nlp, i, j, k
44 integer*4 :: unit_mpi
45 ! integer*4 :: ic1, ic2, ic3, ic4, ic5, ic6, ic7, ic8
46 integer*4, dimension(nmat) :: sdeg, tag_mat
47 integer*4, dimension(nn*nn*nn) :: vtk_numbering_map, node_numbering_vtkwrite, loc_nod_indx
48 integer*4, dimension(nmat_nlp) :: tag_mat_nlp
49 integer*4, dimension(0:cs_nnz_loc) :: cs_loc
50 integer*4, dimension(nn_loc) :: nlp_flag
51 ! integer*4, dimension(:), allocatable :: nlp_flag_el
52
53 real*8, dimension(nmat,4) :: prop_mat
54 real*8, dimension(nn_loc) :: rho, lambda, mu
55 real*8, dimension(nn_loc) :: xx_loc, yy_loc, zz_loc
56
57 nel_loc = cs_loc(0) - 1
58
59 ! allocate(nlp_flag_el(nel_loc))
60 rho = 0.d0; lambda = 0.d0; mu = 0.d0; nlp_flag = 0;
61
62 if (mpi_id.eq.0) write(*,'(A)')
63 if (mpi_id.eq.0) write(*,'(A)')'------Writing VTK file - SCALAR ----------'
64
65 prop_name = 'MECH_PROP_NLP'
66 write(file_name,'(a,i5.5,a)') trim(prop_name),mpi_id,'.vtk'
67 unit_mpi = 2500 + mpi_id
68
69 !----------------------------------------------------------------------
70 open(unit_mpi,file=file_name)
71 write(unit_mpi,'(a)') '# vtk DataFile Version 3.1'
72 write(unit_mpi,'(a)') 'material model VTK file'
73 write(unit_mpi,'(a)') 'ASCII'
74 write(unit_mpi,'(a)') 'DATASET UNSTRUCTURED_GRID'
75 write(unit_mpi, '(a,i12,a)') 'POINTS ', nn_loc, ' float'
76
77 ! Node Coordinates
78 do inode=1,nn_loc
79 write(unit_mpi,'(3e20.12)') xx_loc(inode),yy_loc(inode),zz_loc(inode)
80 enddo
81 write(unit_mpi,*) ''
82
83 ! Connectivity (note: node indices for vtk start at 0)
84 write(unit_mpi,'(a,i12,i12)') "CELLS ",nel_loc,nel_loc*(nn*nn*nn + 1)
85
86 write(temp_char,*)'(i12,',nn*nn*nn,'i12)'
87 do ie=1,nel_loc
88 im = cs_loc(cs_loc(ie -1) + 0 )
89
90 ! nn = sdeg(im) +1
91 do i=1,(nn*nn*nn)
92 loc_nod_indx(i) = cs_loc(cs_loc(ie -1) + i) - 1
93 enddo
94 do i=1,(nn*nn*nn)
95 node_numbering_vtkwrite(i) = loc_nod_indx(vtk_numbering_map(i))
96 enddo
97
98 ! Writing Element-node connectivity
99 ! write(unit_mpi,'(9i12)') 8, ic2, ic6, ic7, ic3, & ! Bottom Surface (front-left node, and then anticlosckwise)
100 ! ic1, ic5, ic8, ic4 ! Top Surface
101
102 write(unit_mpi,temp_char) nn*nn*nn, (node_numbering_vtkwrite(i), i=1,(nn*nn*nn))
103
104 ! dum_int = 0
105 ! do im_nlp=1,nmat_nlp
106 ! if (tag_mat_nlp(im_nlp) .eq. tag_mat(im)) dum_int = 1
107 ! enddo
108
109 do k=1,nn
110 do j=1,nn
111 do i=1,nn
112 inode = cs_loc( cs_loc(ie-1) + nn*nn*(k-1) + nn*(j-1) + i )
113 rho(inode) = prop_mat(im,1); lambda(inode) = prop_mat(im,2);
114 mu(inode) = prop_mat(im,3); !nlp_flag(inode) = nlp_flag_el(ie);
115 enddo
116 enddo
117 enddo
118
119 enddo
120 write(unit_mpi,*) ''
121
122 ! vtkCellType: hexahedrons (ID = 12 for linear 8 noded hex)
123 ! ID = 29 for HEX27 (triquadratic Hex)
124 ! 72 for lagrange-higher order hex
125 ! Reference : https://vtk.org/doc/nightly/html/vtkCellType_8h_source.html
126 ! https://visit-sphinx-github-user-manual.readthedocs.io/en/v3.3.0/data_into_visit/VTKFormat.html
127 write(unit_mpi,'(a,i12)') "CELL_TYPES ",nel_loc
128 write(unit_mpi,'(6i12)') (72,ie=1,nel_loc)
129 write(unit_mpi,*) ''
130
131 ! Writing node data------------------------------------------------------------------
132 write(unit_mpi,'(a,i12)') "POINT_DATA ",nn_loc
133
134 ! Density
135 write(unit_mpi,'(a)') 'SCALARS DENSITY float'
136 write(unit_mpi,'(a)') "LOOKUP_TABLE default"
137 do inode = 1,nn_loc
138 write(unit_mpi,*) rho(inode)
139 enddo
140 write(unit_mpi,*) ''
141
142 ! lambda
143 write(unit_mpi,'(a)') 'SCALARS Lambda float'
144 write(unit_mpi,'(a)') "LOOKUP_TABLE default"
145 do inode = 1,nn_loc
146 write(unit_mpi,*) lambda(inode)
147 enddo
148 write(unit_mpi,*) ''
149
150 ! mu
151 write(unit_mpi,'(a)') 'SCALARS Mu float'
152 write(unit_mpi,'(a)') "LOOKUP_TABLE default"
153 do inode = 1,nn_loc
154 write(unit_mpi,*) mu(inode)
155 enddo
156 write(unit_mpi,*) ''
157
158 ! ! nonlinear flag
159 ! write(unit_mpi,'(a)') 'SCALARS NLP_FLAG integer'
160 ! write(unit_mpi,'(a)') "LOOKUP_TABLE default"
161 ! do inode = 1,nn_loc
162 ! write(unit_mpi,*) nlp_flag(inode)
163 ! enddo
164 ! write(unit_mpi,*) ''
165 !-----------------------------------------------------------------------------------
166
167 ! Writing Element data------------------------------------------------------------------
168 ! write(unit_mpi,'(a,i12)') "CELL_DATA ",nel_loc
169
170 ! ! nonlinear flag
171 ! write(unit_mpi,'(a)') 'SCALARS NLP_FLAG_EL integer'
172 ! write(unit_mpi,'(a)') "LOOKUP_TABLE default"
173 ! do ie = 1,nel_loc
174 ! write(unit_mpi,*) nlp_flag_el(ie)
175 ! enddo
176 ! write(unit_mpi,*) ''
177 !-----------------------------------------------------------------------------------
178
179
180 close(unit_mpi)
181 !------------------------------------------------------------------
182
183 if (mpi_id.eq.0) write(*,'(A)')'Completed.'
184
185 end subroutine write_vtk_mech_prop
186
187
188
189
190
191 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
192 ! Maps the conversion -> speed node numbering to VTK nodenumbering for VTK_LAGRANGE_HEXAHEDRAN of order n
193 ! Reference: https://www.kitware.com/main/wp-content/uploads/2018/09/Source_Issue_43.pdf
194 ! https://gitlab.kitware.com/vtk/vtk/-/issues/17746
195 subroutine vtk_node_numbering_map(nn, vtk_numbering)
196
197 integer*4, intent(in) :: nn
198 integer*4, dimension(nn*nn*nn), intent(out) :: vtk_numbering
199
200 integer*4 :: i, j, k, i1, i2, i3, ncount
201
202 vtk_numbering = 0;
203
204 if (nn.le.1) call exit()
205
206 ! 8 Vertices (first 4 belong to bottom face), next 4 belong to top face
207 k= 1; j=nn; i=nn; vtk_numbering(1) = nn*nn*(k-1) + nn*(j-1) + i;
208 k=nn; j=nn; i=nn; vtk_numbering(2) = nn*nn*(k-1) + nn*(j-1) + i;
209 k=nn; j=nn; i= 1; vtk_numbering(3) = nn*nn*(k-1) + nn*(j-1) + i;
210 k= 1; j=nn; i= 1; vtk_numbering(4) = nn*nn*(k-1) + nn*(j-1) + i;
211
212 k= 1; j= 1; i=nn; vtk_numbering(5) = nn*nn*(k-1) + nn*(j-1) + i;
213 k=nn; j= 1; i=nn; vtk_numbering(6) = nn*nn*(k-1) + nn*(j-1) + i;
214 k=nn; j= 1; i= 1; vtk_numbering(7) = nn*nn*(k-1) + nn*(j-1) + i;
215 k= 1; j= 1; i= 1; vtk_numbering(8) = nn*nn*(k-1) + nn*(j-1) + i;
216
217 if (nn.le.2) return
218
219 ! 12 Edges
220 do i1 = 2,(nn-1)
221 !Bottom face
222 ! Edge 1 - Bottom - X1
223 k=i1; j=nn; i=nn; vtk_numbering(8 + i1 -1) = nn*nn*(k-1) + nn*(j-1) + i;
224
225 ! Edge 2 - Bottom - Y2
226 k=nn; j=nn; i=nn-i1+1; vtk_numbering(8 + (nn-2) + i1 -1) = nn*nn*(k-1) + nn*(j-1) + i;
227
228 ! Edge 3 - Bottom - X2
229 k=i1; j=nn; i= 1; vtk_numbering(8 + 2*(nn-2) + i1 -1) = nn*nn*(k-1) + nn*(j-1) + i;
230
231 ! Edge 4 - Bottom - Y2
232 k= 1; j=nn; i=nn-i1+1; vtk_numbering(8 + 3*(nn-2) + i1-1) = nn*nn*(k-1) + nn*(j-1) + i;
233
234
235 !Top Face
236 ! Edge 5 - Top - X1
237 k=i1; j= 1; i=nn; vtk_numbering(8 + 4*(nn-2) + i1-1) = nn*nn*(k-1) + nn*(j-1) + i;
238
239 ! Edge 6 - Top - Y2
240 k=nn; j= 1; i=nn-i1+1; vtk_numbering(8 + 5*(nn-2) + i1-1) = nn*nn*(k-1) + nn*(j-1) + i;
241
242 ! Edge 7 - Top - X2
243 k=i1; j= 1; i= 1; vtk_numbering(8 + 6*(nn-2) + i1-1) = nn*nn*(k-1) + nn*(j-1) + i;
244
245 ! Edge 8 - Top - Y2
246 k= 1; j= 1; i=nn-i1+1; vtk_numbering(8 + 7*(nn-2) + i1-1) = nn*nn*(k-1) + nn*(j-1) + i;
247
248
249 !Vertical Edges
250 ! Edge 9 - front-left
251 k= 1; j=nn-i1+1; i=nn; vtk_numbering(8 + 8*(nn-2) + i1-1) = nn*nn*(k-1) + nn*(j-1) + i;
252
253 ! Edge 10 - front-right
254 k=nn; j=nn-i1+1; i=nn; vtk_numbering(8 + 9*(nn-2) + i1-1) = nn*nn*(k-1) + nn*(j-1) + i;
255
256 ! Edge 11 - back-left
257 k= 1; j=nn-i1+1; i= 1; vtk_numbering(8 + 10*(nn-2) + i1-1) = nn*nn*(k-1) + nn*(j-1) + i;
258
259 ! Edge 12 - back-right
260 k=nn; j=nn-i1+1; i= 1; vtk_numbering(8 + 11*(nn-2) + i1-1) = nn*nn*(k-1) + nn*(j-1) + i;
261 enddo
262
263
264 !Other Exterior nodes present on 6 faces
265 ncount = 0;
266 do i1 = 2,(nn-1)
267 do i2 = 2,(nn-1)
268 ncount = ncount+1;
269
270 !1st Face - side left (Xmin)
271 k= 1; j=nn-i1+1; i=nn-i2+1; vtk_numbering(8 + 12*(nn-2) + ncount) = nn*nn*(k-1) + nn*(j-1) + i;
272
273 !2nd Face - side left (Xmax)
274 k=nn; j=nn-i1+1; i=nn-i2+1; vtk_numbering(8 + 12*(nn-2) + 1*(nn-2)*(nn-2) + ncount) = nn*nn*(k-1) + nn*(j-1) + i;
275
276 !3rd Face - front (Ymin)
277 k=i2; j=nn-i1+1; i=nn; vtk_numbering(8 + 12*(nn-2) + 2*(nn-2)*(nn-2) + ncount) = nn*nn*(k-1) + nn*(j-1) + i;
278
279 !4th Face - Back (Ymax)
280 k=i2; j=nn-i1+1; i= 1; vtk_numbering(8 + 12*(nn-2) + 3*(nn-2)*(nn-2) + ncount) = nn*nn*(k-1) + nn*(j-1) + i;
281
282 !5th Face - Bottom (Zmin)
283 k=i2; j=nn; i=nn-i1+1; vtk_numbering(8 + 12*(nn-2) + 4*(nn-2)*(nn-2) + ncount) = nn*nn*(k-1) + nn*(j-1) + i;
284
285 !6th Face - Top (Zmax)
286 k=i2; j=1; i=nn-i1+1; vtk_numbering(8 + 12*(nn-2) + 5*(nn-2)*(nn-2) + ncount) = nn*nn*(k-1) + nn*(j-1) + i;
287 enddo
288 enddo
289
290 ! Volume nodes
291 ncount = 0;
292 do i1=2,(nn-1)
293 do i2=2,(nn-1)
294 do i3=2,(nn-1)
295 ncount=ncount+1;
296 k=i3; j=nn-i1+1; i=nn-i2+1; vtk_numbering(8 + 12*(nn-2) + 6*(nn-2)*(nn-2) + ncount) = nn*nn*(k-1) + nn*(j-1) + i;
297 enddo
298 enddo
299 enddo
300
301 end subroutine vtk_node_numbering_map
write_vtk_mech_prop
subroutine write_vtk_mech_prop(nn_loc, cs_nnz_loc, cs_loc, nmat, sdeg, prop_mat, tag_mat, nmat_nlp, tag_mat_nlp, xx_loc, yy_loc, zz_loc, mpi_id, nn, vtk_numbering_map)
...Writing VTK file to visualise in Paraview
Definition WRITE_VTK_MECH_PROP.f90:37
vtk_node_numbering_map
subroutine vtk_node_numbering_map(nn, vtk_numbering)
Definition WRITE_VTK_MECH_PROP.f90:196