Computes max values for Peak Ground Map.
51
52
53 implicit none
54
55 integer*4 :: nnod_loc, ne_loc, cs_nnz_loc, nm
56 integer*4 :: nn, im, imon, iaz
57 integer*4 :: i, j, k, is, in, id
58 integer*4 :: ie, ielem
59 integer*4, dimension(0:cs_nnz_loc) :: cs_loc
60 integer*4, dimension(nm) :: sdeg_mat
61 integer*4, dimension(ne_loc) :: local_el_num
62
63 real*8 :: uxm, uym, uzm
64 real*8 :: vxm, vym, vzm
65 real*8 :: axm, aym, azm
66
67
68 real*8, dimension(:), allocatable :: ct,ww
69 real*8, dimension(:,:), allocatable :: dd
70 real*8, dimension(:,:,:), allocatable :: ux_el, uy_el, uz_el
71
72 real*8, dimension(3*nnod_loc) :: u_1,u0,u1
73
74 integer*4 :: imonpgm, nmonitpgm,ndt_monitpgm
75 integer*4, dimension(nmonitpgm) :: node_mpgm
76 integer*4, dimension(nmonitpgm) :: elem_mpgm
77 real*8 :: tmp, dt
78 real*8 :: rotang_monitpgm
79 real*8 :: upm,unm,vpm,vnm,apm,anm
80 real*8, dimension(nmonitpgm) :: xr_mpgm,yr_mpgm,zr_mpgm
81 real*8, dimension(nmonitpgm,9) :: max_u,max_v,max_a
82 real*8, dimension(nmonitpgm,3) :: max_o
83
84
85 real*8 :: variable1m,variable2m,variable3m
86 real*8 :: variable4m,variable5m,variable6m
87
88 real*8, dimension(3*nnod_loc) :: omega
89
90 real*8, dimension(:,:,:), allocatable :: variable1_el,variable2_el
91 real*8, dimension(:,:,:), allocatable :: variable4_el,variable5_el
92
93
94
95 upm = 0.0d0; unm = 0.0d0; vpm = 0.0d0; vnm = 0.0d0; apm = 0.0d0; anm
96
97
98 do imonpgm = 1,nmonitpgm
99
100 ielem = elem_mpgm(imonpgm)
102
103 if (ie .ne. 0) then
104
105 im = cs_loc(cs_loc(ie -1) +0)
106 nn = sdeg_mat(im) +1
107 allocate(ct(nn),ww(nn),dd(nn,nn))
108 allocate(ux_el(nn,nn,nn),uy_el(nn,nn,nn),uz_el(nn,nn,nn))
109 allocate(variable1_el(nn,nn,nn),variable2_el(nn,nn,nn),variable3_el
110 allocate(variable4_el(nn,nn,nn),variable5_el(nn,nn,nn),variable6_el
112
113 do k = 1,nn
114 do j = 1,nn
115 do i = 1,nn
116 is = nn*nn*(k -1) +nn*(j -1) +i
117 in = cs_loc(cs_loc(ie -1) + is)
118
119 iaz = 3*(in -1) +1
120 ux_el(i,j,k) = u1(iaz)
121 iaz = 3*(in -1) +2
122 uy_el(i,j,k) = u1(iaz)
123 iaz = 3*(in -1) +3
124 uz_el(i,j,k) = u1(iaz)
125 enddo
126 enddo
127 enddo
128
130 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm
132 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm
134 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm
135
136 do k = 1,nn
137 do j = 1,nn
138 do i = 1,nn
139 is = nn*nn*(k -1) +nn*(j -1) +i
140 in = cs_loc(cs_loc(ie -1) + is)
141
142 iaz = 3*(in -1) +1
143 ux_el(i,j,k) = (3.0d0*u1(iaz) - 4.0d0*u0(iaz) +
144 iaz = 3*(in -1) +2
145 uy_el(i,j,k) = (3.0d0*u1(iaz) - 4.0d0*u0(iaz) +
146 iaz = 3*(in -1) +3
147 uz_el(i,j,k) = (3.0d0*u1(iaz) - 4.0d0*u0(iaz) +
148 enddo
149 enddo
150 enddo
151
153 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm),vxm
155 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm),vym
157 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm),vzm
158
159 do k = 1,nn
160 do j = 1,nn
161 do i = 1,nn
162 is = nn*nn*(k -1) +nn*(j -1) +i
163 in = cs_loc(cs_loc(ie -1) + is)
164
165 iaz = 3*(in -1) +1
166 ux_el(i,j,k) = (u1(iaz) -2.0*u0(iaz) +u_1(iaz))
167 iaz = 3*(in -1) +2
168 uy_el(i,j,k) = (u1(iaz) -2.0*u0(iaz) +u_1(iaz))
169 iaz = 3*(in -1) +3
170 uz_el(i,j,k) = (u1(iaz) -2.0*u0(iaz) +u_1(iaz))
171 enddo
172 enddo
173 enddo
174
176 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm
178 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm
180 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm
181
182 if (dabs(uxm).lt.1.0e-30) uxm = 0.0e+00
183 if (dabs(uym).lt.1.0e-30) uym = 0.0e+00
184 if (dabs(uzm).lt.1.0e-30) uzm = 0.0e+00
185 if (dabs(vxm).lt.1.0e-30) vxm = 0.0e+00
186 if (dabs(vym).lt.1.0e-30) vym = 0.0e+00
187 if (dabs(vzm).lt.1.0e-30) vzm = 0.0e+00
188 if (dabs(axm).lt.1.0e-30) axm = 0.0e+00
189 if (dabs(aym).lt.1.0e-30) aym = 0.0e+00
190 if (dabs(azm).lt.1.0e-30) azm = 0.0e+00
191
192 if (rotang_monitpgm.ne.0.0d0) then
193
194
195
196
197 upm = uxm * cos(rotang_monitpgm) + uym * sin(rotang_monitpgm
198 unm = -uxm * sin(rotang_monitpgm) + uym * cos(rotang_monitpgm
199 vpm = vxm * cos(rotang_monitpgm) + vym * sin(rotang_monitpgm
200 vnm = -vxm * sin(rotang_monitpgm) + vym * cos(rotang_monitpgm
201 apm = axm * cos(rotang_monitpgm) + aym * sin(rotang_monitpgm
202 anm = -axm * sin(rotang_monitpgm) + aym * cos(rotang_monitpgm
203
204 uxm = upm
205 uym = unm
206 vxm = vpm
207 vym = vnm
208 axm = apm
209 aym = anm
210 endif
211
212 if (dabs(uxm).gt.max_u(imonpgm,1)) max_u(imonpgm,1) = dabs
213 if (dabs(uym).gt.max_u(imonpgm,2)) max_u(imonpgm,2) = dabs
214 if (dabs(uzm).gt.max_u(imonpgm,3)) max_u(imonpgm,3) = dabs
215
216 tmp = dsqrt(dabs(uxm*uym))
217 if (tmp.gt.max_u(imonpgm,4)) max_u(imonpgm,4) = tmp
218
219 tmp = ((dabs(uxm)+dabs(uym))/2)
220 if (tmp.gt.max_u(imonpgm,5)) max_u(imonpgm,5) = tmp
221
222 tmp = (dabs(uxm*uym*uzm))**(0.3333333)
223 if (tmp.gt.max_u(imonpgm,6)) max_u(imonpgm,6) = tmp
224
225 tmp = ((dabs(uxm)+dabs(uym)+dabs(uzm))/3)
226 if (tmp.gt.max_u(imonpgm,7)) max_u(imonpgm,7) = tmp
227
228 tmp = dsqrt((uxm)**2 + (uym)**2)
229 if (tmp.gt.max_u(imonpgm,8)) max_u(imonpgm,8) = tmp
230
231 tmp = dsqrt((uxm)**2 + (uym)**2 + (uzm)**2)
232 if (tmp.gt.max_u(imonpgm,9)) max_u(imonpgm,9) = tmp
233
234
235 if (dabs(vxm).gt.max_v(imonpgm,1)) max_v(imonpgm,1) = dabs
236 if (dabs(vym).gt.max_v(imonpgm,2)) max_v(imonpgm,2) = dabs
237 if (dabs(vzm).gt.max_v(imonpgm,3)) max_v(imonpgm,3) = dabs
238
239 tmp = dsqrt(dabs(vxm*vym))
240 if (tmp.gt.max_v(imonpgm,4)) max_v(imonpgm,4) = tmp
241
242 tmp = ((dabs(vxm)+dabs(vym))/2)
243 if (tmp.gt.max_v(imonpgm,5)) max_v(imonpgm,5) = tmp
244
245 tmp = (dabs(vxm*vym*vzm))**(0.3333333)
246 if (tmp.gt.max_v(imonpgm,6)) max_v(imonpgm,6) = tmp
247
248 tmp = ((dabs(vxm)+dabs(vym)+dabs(vzm))/3)
249 if (tmp.gt.max_v(imonpgm,7)) max_v(imonpgm,7) = tmp
250
251 tmp = dsqrt((vxm)**2 + (vym)**2)
252 if (tmp.gt.max_v(imonpgm,8)) max_v(imonpgm,8) = tmp
253
254 tmp = dsqrt((vxm)**2 + (vym)**2 + (vzm)**2)
255 if (tmp.gt.max_v(imonpgm,9)) max_v(imonpgm,9) = tmp
256
257 if (dabs(axm).gt.max_a(imonpgm,1)) max_a(imonpgm,1) = dabs
258 if (dabs(aym).gt.max_a(imonpgm,2)) max_a(imonpgm,2) = dabs
259 if (dabs(azm).gt.max_a(imonpgm,3)) max_a(imonpgm,3) = dabs
260
261 tmp = dsqrt(dabs(axm*aym))
262 if (tmp.gt.max_a(imonpgm,4)) max_a(imonpgm,4) = tmp
263
264 tmp = ((dabs(axm)+dabs(aym))/2)
265 if (tmp.gt.max_a(imonpgm,5)) max_a(imonpgm,5) = tmp
266
267 tmp = (dabs(axm*aym*azm))**(0.3333333)
268 if (tmp.gt.max_a(imonpgm,6)) max_a(imonpgm,6) = tmp
269
270 tmp = ((dabs(axm)+dabs(aym)+dabs(azm))/3)
271 if (tmp.gt.max_a(imonpgm,7)) max_a(imonpgm,7) = tmp
272
273 tmp = dsqrt((axm)**2 + (aym)**2)
274 if (tmp.gt.max_a(imonpgm,8)) max_a(imonpgm,8) = tmp
275
276 tmp = dsqrt((axm)**2 + (aym)**2 + (azm)**2)
277 if (tmp.gt.max_a(imonpgm,9)) max_a(imonpgm,9) = tmp
278
279
280
281
282
283 do k = 1,nn
284 do j = 1,nn
285 do i = 1,nn
286 is = nn*nn*(k -1) +nn*(j -1) +i
287 in = cs_loc(cs_loc(ie -1) + is)
288
289 iaz = 3*(in -1) +1
290 variable1_el(i,j,k) = omega(iaz)
291 iaz = 3*(in -1) +2
292 variable2_el(i,j,k) = omega(iaz)
293 iaz = 3*(in -1) +3
294 variable3_el(i,j,k) = omega(iaz)
295 enddo
296 enddo
297 enddo
298
299
301 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm),variable1m
302
304 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm),variable2m
305
307 xr_mpgm(imonpgm),yr_mpgm(imonpgm),zr_mpgm(imonpgm),variable3m
308
309 if (dabs(variable1m).lt.1.0e-30) variable1m = 0.0e+00
310 if (dabs(variable2m).lt.1.0e-30) variable2m = 0.0e+00
311 if (dabs(variable3m).lt.1.0e-30) variable3m = 0.0e+00
312 if (dabs(variable1m).gt.max_o(imonpgm,1)) max_o(imonpgm,1) = dabs
313 if (dabs(variable2m).gt.max_o(imonpgm,2)) max_o(imonpgm,2) = dabs
314 if (dabs(variable3m).gt.max_o(imonpgm,3)) max_o(imonpgm,3) = dabs
315
316
317 deallocate(ct,ww,dd)
318 deallocate(ux_el,uy_el,uz_el)
319 deallocate(variable1_el,variable2_el,variable3_el)
320 deallocate(variable4_el,variable5_el,variable6_el)
321
322
323 endif
324 enddo
325
subroutine get_indloc_from_indglo(local_el, nel_loc, ie, ic)
Returns local id from global id.
subroutine get_monitor_value(nb_nod, xq, val, xref, yref, zref, re
Computes the mean value.
subroutine make_lgl_nw(nb_pnt, xq, wq, dd)
Makes Gauss-Legendre-Lobatto nodes, weigths and spectral derivatives.