MPI-AMRVAC  3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
mod_collapse.t
Go to the documentation of this file.
1 !> Collapses D-dimensional output to D-1 view by line-of-sight integration
2 ! for now: only for non-stretched cartesian cases, with LOS along coordinate
3 ! writes out as either csv or vti file (collapse_type)
6  use mod_comm_lib, only: mpistop
7  implicit none
8 
9 contains
10 !=============================================================================
11 subroutine write_collapsed
13 ! Writes a collapsed view of the data integrated over one grid-direction.
14 ! E.g. column density maps.
15 ! Uses flat interpolation throughout.
16 ! by Oliver Porth
17 ! 6.Nov 2013
18 integer :: idir
19 logical, save :: firstcollapse=.true.
20 !-----------------------------------------------------------------------------
21 
22 if (firstcollapse) then
23  firstcollapse=.false.
24 end if
25 
26 do idir=1, ndim
27  if (collapse(idir)) call put_collapse(idir)
28 end do
29 
31 end subroutine write_collapsed
32 !=============================================================================
33 subroutine put_collapse(dir)
37 integer, intent(in) :: dir
38 ! .. local ..
39 integer :: jgrid, igrid, Morton_no
40 double precision,dimension(0:nw+nwauxio) :: normconv
41 !-----------------------------------------------------------------------------
42 if(.not.slab) call mpistop("collapse only for slab cartesian cases")
43 
44 call allocate_collapsed(dir)
45 
46 jgrid=0
47 do morton_no=morton_start(mype),morton_stop(mype)
48  igrid=sfc_to_igrid(morton_no)
49  jgrid=jgrid+1
50  call alloc_subnode(jgrid,dir,nwauxio)
51  call collapse_subnode(igrid,jgrid,dir,normconv)
52  call integrate_subnode(igrid,jgrid,dir)
53 end do
54 
55 ! Reduce to head-node:
56 if (mype==0) then
57  call mpi_reduce(mpi_in_place,collapseddata,size(collapseddata),mpi_double_precision,mpi_sum,0,icomm,ierrmpi)
58 else
59  call mpi_reduce(collapseddata,collapseddata,size(collapseddata),mpi_double_precision,mpi_sum,0,icomm,ierrmpi)
60 end if
61 call mpi_barrier(icomm, ierrmpi)
62 
63 {^nooned
64 select case(collapse_type)
65 case('vti')
66  call output_collapsed_vti(dir,normconv)
67 case('csv')
68  call output_collapsed_csv(dir,normconv)
69 case('default')
70  call mpistop("Unknown filetype for collapsed views")
71 end select
72  }{^ifoned
73 call mpistop("sorry, 1D collapsed output routine not yet implemented (should be easy)...")
74 }
75 
76 ! If we need the subnodes later, remove deallocation here:
77 do jgrid=1,morton_stop(mype)-morton_start(mype)+1
78  call dealloc_subnode(jgrid)
79 end do
80 deallocate(collapseddata)
81 
82 end subroutine put_collapse
83 !=============================================================================
84 subroutine output_collapsed_csv(dir,normconv)
87 integer, intent(in) :: dir
88 double precision,dimension(0:nw+nwauxio),intent(in):: normconv
89 character(len=1024) :: filename, outfilehead, line
90 logical :: fileopen
91 integer :: iw
92 character(len=name_len) :: wnamei(1:nw+nwauxio),xandwnamei(1:ndim+nw+nwauxio)
93 integer, dimension(ndim) :: myshape
94 {^nooned
95 integer :: ix^DM
96 double precision :: dxdim^DM, xdim^LIM^DM
97 }
98 !-----------------------------------------------------------------------------
99 
100 if (mype==0) then
101 
102 ! Get coordinates:
103 {^ifthreed
104 select case(dir)
105 case (1)
106  dxdim1 = dx(2,collapselevel)
107  dxdim2 = dx(3,collapselevel)
108  xdim^lim1=xprob^lim2;
109  xdim^lim2=xprob^lim3;
110 case (2)
111  dxdim1 = dx(1,collapselevel)
112  dxdim2 = dx(3,collapselevel)
113  xdim^lim1=xprob^lim1;
114  xdim^lim2=xprob^lim3;
115 case (3)
116  dxdim1 = dx(1,collapselevel)
117  dxdim2 = dx(2,collapselevel)
118  xdim^lim1=xprob^lim1;
119  xdim^lim2=xprob^lim2;
120 case default
121  dxdim1 = 1
122  dxdim2 = 1
123  xdim^lim1=xprob^lim2;
124  xdim^lim2=xprob^lim3;
125  call mpistop("slice direction not clear in output_collapsed_csv")
126 end select
127 }
128 {^iftwod
129 select case(dir)
130 case (1)
131  dxdim1 = dx(2,collapselevel)
132  xdim^lim1=xprob^lim2;
133 case (2)
134  dxdim1 = dx(1,collapselevel)
135  xdim^lim1=xprob^lim1;
136 case default
137  dxdim1 = 1
138  xdim^lim1=xprob^lim2;
139  call mpistop("slice direction not clear in output_collapsed_csv")
140 end select
141 }
142 
143  inquire(unitcollapse,opened=fileopen)
144  if(.not.fileopen)then
145  ! generate filename:
146  write(filename,"(a,i1.1,a,i1.1,a,i4.4,a)") trim(base_filename) // '_d', &
147  dir,'_l',collapselevel,'_n',collapsenext,'.csv'
148  open(unitcollapse,file=filename,status='unknown',form='formatted')
149  end if
150  ! get and write the header:
151  call getheadernames(wnamei,xandwnamei,outfilehead)
152  line=''
153  do iw=1,ndim+nw+nwauxio-1
154  if (iw .eq. dir) cycle
155  line = trim(line)//trim(xandwnamei(iw))//', '
156  end do
157  line = trim(line)//trim(xandwnamei(ndim+nw+nwauxio))
158  write(unitcollapse,'(a)')trim(line)
159  myshape = shape(collapseddata)
160 {^nooned
161 {^dm& do ix^dmb = 1,myshape(^dmb)\}
162  ! following assumes uniform cartesian grid
163  write(unitcollapse,'(200(1pe20.12))') &
164  {^dm& dxdim^dm*dble(ix^dm)+xdimmin^dm}, &
165  (collapseddata(ix^dm,iw)*normconv(iw),iw=1,nw+nwauxio)
166 {^dm& enddo\}
167 }
168 close(unitcollapse)
169 
170 end if
171 end subroutine output_collapsed_csv
172 !=============================================================================
173 subroutine output_collapsed_vti(dir,normconv)
176 integer, intent(in) :: dir
177 double precision,dimension(0:nw+nwauxio),intent(in):: normconv
178 character(len=1024) :: filename, outfilehead, line
179 logical :: fileopen
180 integer :: iw
181 character(len=name_len) :: wnamei(1:nw+nwauxio),xandwnamei(1:ndim+nw+nwauxio)
182 integer, dimension(ndim) :: myshape
183 {^nooned
184 integer :: ix^DM
185 double precision :: dxdim^DM, xdim^LIM^DM
186 }
187 double precision :: origin(1:3), spacing(1:3)
188 integer :: wholeExtent(1:6), size_single, length, size_length
189 integer*8 :: offset
190 character :: buf
191 !-----------------------------------------------------------------------------
192 if (mype==0) then
193 offset=0
194 size_single=4
195 size_length=4
196 ! Get coordinates:
197 {^ifthreed
198 select case(dir)
199 case (1)
200  dxdim1 = dx(2,1)*2.0d0**(1-collapselevel)
201  dxdim2 = dx(3,1)*2.0d0**(1-collapselevel)
202  xdim^lim1=xprob^lim2;
203  xdim^lim2=xprob^lim3;
204 case (2)
205  dxdim1 = dx(1,1)*2.0d0**(1-collapselevel)
206  dxdim2 = dx(3,1)*2.0d0**(1-collapselevel)
207  xdim^lim1=xprob^lim1;
208  xdim^lim2=xprob^lim3;
209 case (3)
210  dxdim1 = dx(1,1)*2.0d0**(1-collapselevel)
211  dxdim2 = dx(2,1)*2.0d0**(1-collapselevel)
212  xdim^lim1=xprob^lim1;
213  xdim^lim2=xprob^lim2;
214 case default
215  dxdim1 = 1
216  dxdim2 = 1
217  xdim^lim1=xprob^lim2;
218  xdim^lim2=xprob^lim3;
219  call mpistop("slice direction not clear in output_collapsed_vti")
220 end select
221 }
222 {^iftwod
223 select case(dir)
224 case (1)
225  dxdim1 = dx(2,1)*2.0d0**(1-collapselevel)
226  xdim^lim1=xprob^lim2;
227 case (2)
228  dxdim1 = dx(1,1)*2.0d0**(1-collapselevel)
229  xdim^lim1=xprob^lim1;
230 case default
231  dxdim1 = 1
232  xdim^lim1=xprob^lim2;
233  call mpistop("slice direction not clear in output_collapsed_vti")
234 end select
235 }
236 
237 origin=0
238 spacing=zero
239 wholeextent=0
240 myshape = shape(collapseddata)
241 {^ifoned
242 length = size_single
243 }
244 {^nooned
245 length = ^dm&myshape(^dm)*
246 length = length*size_single
247 {^dm&wholeextent(^dm*2)=myshape(^dm); }
248 {^dm&spacing(^dm) = dxdim^dm; }
249 {^dm&origin(^d) = xdimmin^dm; }
250 }
251 
252  inquire(unitcollapse,opened=fileopen)
253  if(.not.fileopen)then
254  ! generate filename:
255  write(filename,"(a,i1.1,a,i1.1,a,i4.4,a)") trim(base_filename)//'_d',dir,&
256  '_l',collapselevel,'_n',collapsenext,'.vti'
257  open(unitcollapse,file=filename,status='unknown',form='formatted')
258  end if
259 ! get the header:
260 call getheadernames(wnamei,xandwnamei,outfilehead)
261 
262 ! generate xml header
263 write(unitcollapse,'(a)')'<?xml version="1.0"?>'
264 write(unitcollapse,'(a)',advance='no') '<VTKFile type="ImageData"'
265 if(type_endian==1)then
266  write(unitcollapse,'(a)')' version="0.1" byte_order="LittleEndian">'
267 else
268  write(unitcollapse,'(a)')' version="0.1" byte_order="BigEndian">'
269 endif
270 ! following corresponds to uniform cartesian grid
271 write(unitcollapse,'(a,3(1pe14.6),a,6(i10),a,3(1pe14.6),a)')' <ImageData Origin="',&
272  origin,'" WholeExtent="',wholeextent,'" Spacing="',spacing,'">'
273 write(unitcollapse,'(a)')'<FieldData>'
274 write(unitcollapse,'(2a)')'<DataArray type="Float32" Name="TIME" ',&
275  'NumberOfTuples="1" format="ascii">'
277 write(unitcollapse,'(a)')'</DataArray>'
278 write(unitcollapse,'(a)')'</FieldData>'
279 
280 ! we write one VTK PIECE
281 write(unitcollapse,'(a,6(i10),a)') &
282  '<Piece Extent="',wholeextent,'">'
283 write(unitcollapse,'(a)')'<CellData>'
284 
285 do iw=1,nw+nwauxio
286  if(iw<=nw) then
287  if(.not.w_write(iw)) cycle
288  endif
289  write(unitcollapse,'(a,a,a,i16,a)')&
290  '<DataArray type="Float32" Name="',trim(wnamei(iw)),&
291  '" format="appended" offset="',offset,'"/>'
292  offset = offset + length + size_length
293 enddo
294 
295 write(unitcollapse,'(a)')'</CellData>'
296 write(unitcollapse,'(a)')'</Piece>'
297 write(unitcollapse,'(a)')'</ImageData>'
298 write(unitcollapse,'(a)')'<AppendedData encoding="raw">'
299 close(unitcollapse)
300 
301 open(unitcollapse,file=filename,access='stream',form='unformatted',position='append')
302 buf='_'
303 write(unitcollapse) trim(buf)
304 
305 do iw=1,nw+nwauxio
306  if(iw<=nw) then
307  if(.not.w_write(iw)) cycle
308  endif
309  write(unitcollapse) length
310  {^ifoned
311  write(unitcollapse) real(collapseddata(iw)*normconv(iw))
312  }
313  {^iftwod
314  write(unitcollapse) (real(collapseddata(ix1,iw)*normconv(iw)),ix1=1,myshape(1))
315  }
316  {^ifthreed
317  write(unitcollapse) ((real(collapseddata(ix1,ix2,iw)*normconv(iw)),ix1=1,&
318  myshape(1)),ix2=1,myshape(2))
319  }
320 enddo
321 
322 close(unitcollapse)
323 open(unitcollapse,file=filename,status='unknown',form='formatted',position='append')
324 write(unitcollapse,'(a)')'</AppendedData>'
325 write(unitcollapse,'(a)')'</VTKFile>'
326 close(unitcollapse)
327 
328 end if
329 
330 end subroutine output_collapsed_vti
331 !=============================================================================
332 subroutine allocate_collapsed(dir)
334 integer, intent(in) :: dir
335 integer :: dim^D
336 integer :: nx^D
337 !-----------------------------------------------------------------------------
338 ! allocate array for the collapsed data:
339 ! number of cells per grid.
340 nx^d=ixmhi^d-ixmlo^d+1;
341 {dim^d=ng^d(1)*2**(collapselevel-1)*nx^d\}
342 {^ifthreed
343 select case(dir)
344 case (1)
345  allocate(collapseddata(&
346  dim2,dim3,nw+nwauxio))
347 case (2)
348  allocate(collapseddata(&
349  dim1,dim3,nw+nwauxio))
350 case (3)
351  allocate(collapseddata(&
352  dim1,dim2,nw+nwauxio))
353 case default
354  call mpistop("slice direction not clear in allocate_collapsed")
355 end select
356 }
357 {^iftwod
358 select case(dir)
359 case (1)
360  allocate(collapseddata(&
361  dim2,nw+nwauxio))
362 case (2)
363  allocate(collapseddata(&
364  dim1,nw+nwauxio))
365 case default
366  call mpistop("slice direction not clear in allocate_collapsed")
367 end select
368 }
369 {^ifoned
370  allocate(collapseddata(nw+nwauxio))
371 }
372 collapseddata = zero
373 
374 end subroutine allocate_collapsed
375 !=============================================================================
376 subroutine integrate_subnode(igrid,jgrid,dir)
379 integer, intent(in) :: igrid, jgrid, dir
380 ! .. local ..
381 type(tree_node_ptr) :: tree
382 integer :: nx^D
383 integer :: ig^D, level, ig^Dtargetmin, ig^Dtargetmax
384 integer :: ix^Dtarget^LIM, idim^Dtarget^LIM
385 integer :: ixMdim^LLIM^D, ix^Dorig, ix^D
386 {^nooned
387 integer :: igdim^DM
388 }
389 !-----------------------------------------------------------------------------
390 tree%node => igrid_to_node(igrid, mype)%node
391 {^d& ig^d = tree%node%ig^d; }
392 level = tree%node%level
393 ! number of cells per grid.
394 nx^d=ixmhi^d-ixmlo^d+1;
395 
396 ! assumes uniform cartesian grid with factor 2 refinement per dimension
397 {^d&
398 if (level > collapselevel) then
399  ig^dtargetmin = int(dble(ig^d-1)/2.0d0**(level-collapselevel))+1
400  ig^dtargetmax = ig^dtargetmin
401 else if (level < collapselevel) then
402  ig^dtargetmin = int(2.0d0**(collapselevel-level))*(ig^d-1)+1
403  ig^dtargetmax = int(2.0d0**(collapselevel-level))*ig^d
404 else
405  ig^dtargetmin = ig^d
406  ig^dtargetmax = ig^d
407 end if
408 ix^dtargetmin = nx^d*(ig^dtargetmin-1)+1
409 ix^dtargetmax = nx^d*ig^dtargetmax
410 \}
411 
412 {^ifthreed
413 select case(dir)
414 case (1)
415  igdim1=(ig2-1)*nx2
416  igdim2=(ig3-1)*nx3
417  idim1target^lim=ix2target^lim;
418  idim2target^lim=ix3target^lim;
419  ixmdim^llim1=ixm^llim2;
420  ixmdim^llim2=ixm^llim3;
421 case (2)
422  igdim1=(ig1-1)*nx1
423  igdim2=(ig3-1)*nx3
424  idim1target^lim=ix1target^lim;
425  idim2target^lim=ix3target^lim;
426  ixmdim^llim1=ixm^llim1;
427  ixmdim^llim2=ixm^llim3;
428 case (3)
429  igdim1=(ig1-1)*nx1
430  igdim2=(ig2-1)*nx2
431  idim1target^lim=ix1target^lim;
432  idim2target^lim=ix2target^lim;
433  ixmdim^llim1=ixm^llim1;
434  ixmdim^llim2=ixm^llim2;
435 case default
436  call mpistop("slice direction not clear in integrate_subnode")
437 end select
438 
439 if (level >= collapselevel) then
440  do ix1orig = ixmdimlo1,ixmdimhi1
441  do ix2orig = ixmdimlo2,ixmdimhi2
442 {^dm& ix^dm = int(dble(ix^dmorig-nghostcells+igdim^dm-1)*2.0d0**(collapselevel-level))+1\}
443  collapseddata(ix1,ix2,1:nw+nwauxio) = collapseddata(ix1,ix2,1:nw+nwauxio) &
444  + ps_sub(jgrid)%w(ix1orig,ix2orig,1:nw+nwauxio) / 2.0d0**(2*(level-collapselevel))
445  end do
446  end do
447 else
448 {^dm&
449  do ix^dm = idim^dmtargetmin,idim^dmtargetmax\}
450  {^dm& ix^dmorig = int(dble(ix^dm-idim^dmtargetmin)/2.0d0**(collapselevel-level))+1+nghostcells \}
451  collapseddata(ix1,ix2,1:nw+nwauxio) = collapseddata(ix1,ix2,1:nw+nwauxio) + ps_sub(jgrid)%w(ix1orig,ix2orig,1:nw+nwauxio)
452  {^dm& enddo\}
453 end if
454 }
455 {^iftwod
456 select case(dir)
457 case (1)
458  igdim1=(ig2-1)*nx2
459  idim1target^lim=ix2target^lim;
460  ixmdim^llim1=ixm^llim2;
461 case (2)
462  igdim1=(ig1-1)*nx1
463  idim1target^lim=ix1target^lim;
464  ixmdim^llim1=ixm^llim1;
465 case default
466  call mpistop("slice direction not clear in integrate_subnode")
467 end select
468 
469 if (level >= collapselevel) then
470  do ix1orig = ixmdimlo1,ixmdimhi1
471 {^dm& ix^dm = int(dble(ix^dmorig-nghostcells+igdim^dm-1)*2.0d0**(collapselevel-level))+1\}
472  collapseddata(ix1,1:nw+nwauxio) = collapseddata(ix1,1:nw+nwauxio) &
473  + ps_sub(jgrid)%w(ix1orig,1:nw+nwauxio) / 2.0d0**(level-collapselevel)
474  end do
475 else
476 {^dm&
477  do ix^dm = idim^dmtargetmin,idim^dmtargetmax\}
478  {^dm& ix^dmorig = int(dble(ix^dm-idim^dmtargetmin)/2.0d0**(collapselevel-level))+1+nghostcells \}
479  collapseddata(ix1,1:nw+nwauxio) = collapseddata(ix1,1:nw+nwauxio) + ps_sub(jgrid)%w(ix1orig,1:nw+nwauxio)
480  {^dm& enddo\}
481 end if
482 }
483 {^ifoned
484 collapseddata(1:nw+nwauxio) = collapseddata(1:nw+nwauxio) + ps_sub(jgrid)%w(1:nw+nwauxio)
485 }
486 
487 end subroutine integrate_subnode
488 !=============================================================================
489 subroutine collapse_subnode(igrid,jgrid,dir,normconv)
492  use mod_physics, only: phys_to_primitive
493  use mod_calculate_xw
494 
495 integer, intent(in) :: igrid, jgrid, dir
496 double precision,dimension(0:nw+nwauxio),intent(out) :: normconv
497 ! .. local ..
498 integer :: ix, iw
499 double precision :: dx^D
500 double precision, dimension(ixMlo^D-1:ixMhi^D,ndim) :: xC_TMP, xC
501 double precision, dimension(ixMlo^D:ixMhi^D,ndim) :: xCC_TMP, xCC
502 double precision, dimension(ixMlo^D-1:ixMhi^D,nw+nwauxio) :: wC_TMP
503 double precision, dimension(ixMlo^D:ixMhi^D,nw+nwauxio) :: wCC_TMP
504 integer :: ixC^L, ixCC^L
505 !-----------------------------------------------------------------------------
506 ps_sub(jgrid)%w=zero
507 dx^d=rnode(rpdx^d_,igrid);
508 
509 call calc_x(igrid,xc,xcc)
510 call calc_grid(unitslice,igrid,xc,xcc,xc_tmp,xcc_tmp,wc_tmp,wcc_tmp,normconv,&
511  ixc^l,ixcc^l,.true.)
512 
513 {^ifthreed
514 select case (dir)
515 case (1)
516  if(slab_uniform) then
517  do ix=ixmlo1,ixmhi1
518  ps_sub(jgrid)%w(ixmlo2:ixmhi2,ixmlo3:ixmhi3,1:nw+nwauxio) = &
519  ps_sub(jgrid)%w(ixmlo2:ixmhi2,ixmlo3:ixmhi3,1:nw+nwauxio) &
520  + wcc_tmp(ix,ixmlo2:ixmhi2,ixmlo3:ixmhi3,1:nw+nwauxio) * dx1
521  end do
522  else
523  do iw=1,nw+nwauxio
524  do ix=ixmlo1,ixmhi1
525  ps_sub(jgrid)%w(ixmlo2:ixmhi2,ixmlo3:ixmhi3,iw) = &
526  ps_sub(jgrid)%w(ixmlo2:ixmhi2,ixmlo3:ixmhi3,iw) &
527  + wcc_tmp(ix,ixmlo2:ixmhi2,ixmlo3:ixmhi3,iw) * &
528  ps(igrid)%dx(ix,ixmlo2:ixmhi2,ixmlo3:ixmhi3,1)
529  end do
530  end do
531  end if
532  ps_sub(jgrid)%x(ixmlo2:ixmhi2,ixmlo3:ixmhi3,1:ndim) = &
533  ps(igrid)%x(ixmlo1,ixmlo2:ixmhi2,ixmlo3:ixmhi3,1:ndim)
534 case (2)
535  if(slab_uniform) then
536  do ix=ixmlo2,ixmhi2
537  ps_sub(jgrid)%w(ixmlo1:ixmhi1,ixmlo3:ixmhi3,1:nw+nwauxio) = &
538  ps_sub(jgrid)%w(ixmlo1:ixmhi1,ixmlo3:ixmhi3,1:nw+nwauxio) &
539  + wcc_tmp(ixmlo1:ixmhi1,ix,ixmlo3:ixmhi3,1:nw+nwauxio) * dx2
540  end do
541  else
542  do iw=1,nw+nwauxio
543  do ix=ixmlo2,ixmhi2
544  ps_sub(jgrid)%w(ixmlo1:ixmhi1,ixmlo3:ixmhi3,iw) = &
545  ps_sub(jgrid)%w(ixmlo1:ixmhi1,ixmlo3:ixmhi3,iw) &
546  + wcc_tmp(ixmlo1:ixmhi1,ix,ixmlo3:ixmhi3,iw) * &
547  ps(igrid)%dx(ixmlo1:ixmhi1,ix,ixmlo3:ixmhi3,2)
548  end do
549  end do
550  endif
551  ps_sub(jgrid)%x(ixmlo1:ixmhi1,ixmlo3:ixmhi3,1:ndim) = &
552  ps(igrid)%x(ixmlo1:ixmhi1,ixmlo2,ixmlo3:ixmhi3,1:ndim)
553 case (3)
554  if(slab_uniform) then
555  do ix=ixmlo3,ixmhi3
556  ps_sub(jgrid)%w(ixmlo1:ixmhi1,ixmlo2:ixmhi2,1:nw+nwauxio) = &
557  ps_sub(jgrid)%w(ixmlo1:ixmhi1,ixmlo2:ixmhi2,1:nw+nwauxio) &
558  + wcc_tmp(ixmlo1:ixmhi1,ixmlo2:ixmhi2,ix,1:nw+nwauxio) * dx3
559  end do
560  else
561  do iw=1,nw+nwauxio
562  do ix=ixmlo3,ixmhi3
563  ps_sub(jgrid)%w(ixmlo1:ixmhi1,ixmlo2:ixmhi2,iw) = &
564  ps_sub(jgrid)%w(ixmlo1:ixmhi1,ixmlo2:ixmhi2,iw) &
565  + wcc_tmp(ixmlo1:ixmhi1,ixmlo2:ixmhi2,ix,iw) * &
566  ps(igrid)%dx(ixmlo1:ixmhi1,ixmlo2:ixmhi2,ix,3)
567  end do
568  end do
569  endif
570  ps_sub(jgrid)%x(ixmlo1:ixmhi1,ixmlo2:ixmhi2,1:ndim) = &
571  ps(igrid)%x(ixmlo1:ixmhi1,ixmlo2:ixmhi2,ixmlo3,1:ndim)
572 case default
573  print*, 'subnode, dir: ', dir
574  call mpistop("slice direction not clear in collapse_subnode")
575 end select
576 }
577 {^iftwod
578 select case (dir)
579 case (1)
580  if(slab_uniform) then
581  do ix=ixmlo1,ixmhi1
582  ps_sub(jgrid)%w(ixmlo2:ixmhi2,1:nw+nwauxio) = &
583  ps_sub(jgrid)%w(ixmlo2:ixmhi2,1:nw+nwauxio) &
584  + wcc_tmp(ix,ixmlo2:ixmhi2,1:nw+nwauxio) * dx1
585  end do
586  else
587  do iw=1,nw+nwauxio
588  do ix=ixmlo1,ixmhi1
589  ps_sub(jgrid)%w(ixmlo2:ixmhi2,iw) = &
590  ps_sub(jgrid)%w(ixmlo2:ixmhi2,iw) &
591  + wcc_tmp(ix,ixmlo2:ixmhi2,iw) * &
592  ps(igrid)%dx(ix,ixmlo2:ixmhi2,1)
593  end do
594  end do
595  end if
596  ps_sub(jgrid)%x(ixmlo2:ixmhi2,1:ndim) = &
597  ps(igrid)%x(ixmlo1,ixmlo2:ixmhi2,1:ndim)
598 case (2)
599  if(slab_uniform) then
600  do ix=ixmlo2,ixmhi2
601  ps_sub(jgrid)%w(ixmlo1:ixmhi1,1:nw+nwauxio) = &
602  ps_sub(jgrid)%w(ixmlo1:ixmhi1,1:nw+nwauxio) &
603  + wcc_tmp(ixmlo1:ixmhi1,ix,1:nw+nwauxio) * dx2
604  end do
605  else
606  do iw=1,nw+nwauxio
607  do ix=ixmlo2,ixmhi2
608  ps_sub(jgrid)%w(ixmlo1:ixmhi1,iw) = &
609  ps_sub(jgrid)%w(ixmlo1:ixmhi1,iw) &
610  + wcc_tmp(ixmlo1:ixmhi1,ix,iw) * &
611  ps(igrid)%dx(ixmlo1:ixmhi1,ix,2)
612  end do
613  end do
614  end if
615  ps_sub(jgrid)%x(ixmlo1:ixmhi1,1:ndim) = &
616  ps(igrid)%x(ixmlo1:ixmhi1,ixmlo2,1:ndim)
617 case default
618  call mpistop("slice direction not clear in collapse_subnode")
619 end select
620 }
621 {^ifoned
622 if(slab_uniform) then
623  do ix=ixmlo1,ixmhi1
624  ps_sub(jgrid)%w(1:nw+nwauxio) = ps_sub(jgrid)%w(1:nw+nwauxio) + wcc_tmp(ix,1:nw+nwauxio) * dx1
625  end do
626 else
627  do iw=1,nw+nwauxio
628  do ix=ixmlo1,ixmhi1
629  ps_sub(jgrid)%w(iw) = ps_sub(jgrid)%w(iw) + wcc_tmp(ix,iw) * ps(igrid)%dx(ix,1)
630  end do
631  end do
632 end if
633 ps_sub(jgrid)%x(1:ndim) = ps(igrid)%x(ixmlo1,1:ndim)
634 }
635 
636 end subroutine collapse_subnode
637 !=============================================================================
638 end module mod_collapse
Handles computations for coordinates and variables in output.
subroutine calc_grid(qunit, igrid, xC, xCC, xC_TMP, xCC_TMP, wC_TMP, wCC_TMP, normconv, ixCL, ixCCL, first)
Compute both corner as well as cell-centered values for output.
subroutine getheadernames(wnamei, xandwnamei, outfilehead)
get all variables names
subroutine calc_x(igrid, xC, xCC)
computes cell corner (xC) and cell center (xCC) coordinates
Collapses D-dimensional output to D-1 view by line-of-sight integration.
Definition: mod_collapse.t:4
subroutine collapse_subnode(igrid, jgrid, dir, normconv)
Definition: mod_collapse.t:490
subroutine output_collapsed_csv(dir, normconv)
Definition: mod_collapse.t:85
subroutine allocate_collapsed(dir)
Definition: mod_collapse.t:333
subroutine put_collapse(dir)
Definition: mod_collapse.t:34
subroutine output_collapsed_vti(dir, normconv)
Definition: mod_collapse.t:174
subroutine write_collapsed
Definition: mod_collapse.t:12
subroutine integrate_subnode(igrid, jgrid, dir)
Definition: mod_collapse.t:377
subroutine, public mpistop(message)
Exit MPI-AMRVAC with an error message.
Definition: mod_comm_lib.t:208
Module with basic grid data structures.
Definition: mod_forest.t:2
integer, dimension(:), allocatable, save sfc_to_igrid
Go from a Morton number to an igrid index (for a single processor)
Definition: mod_forest.t:53
integer, dimension(:), allocatable, save morton_start
First Morton number per processor.
Definition: mod_forest.t:62
integer, dimension(:), allocatable, save morton_stop
Last Morton number per processor.
Definition: mod_forest.t:65
type(tree_node_ptr), dimension(:,:), allocatable, save igrid_to_node
Array to go from an [igrid, ipe] index to a node pointer.
Definition: mod_forest.t:32
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, parameter unitslice
double precision global_time
The global simulation time.
double precision xprob
minimum and maximum domain boundaries for each dimension
integer, parameter ndim
Number of spatial dimensions for grid variables.
double precision time_convert_factor
Conversion factor for time unit.
integer icomm
The MPI communicator.
integer, dimension(:), allocatable ng
number of grid blocks in domain per dimension, in array over levels
integer mype
The rank of the current MPI task.
logical, dimension(ndim) collapse
If collapse(DIM) is true, generate output integrated over DIM.
integer, dimension(:), allocatable, parameter d
integer ixm
the mesh range of a physical block without ghost cells
integer ierrmpi
A global MPI error return code.
character(len=std_len) collapse_type
logical slab
Cartesian geometry or not.
integer nwauxio
Number of auxiliary variables that are only included in output.
logical, dimension(:), allocatable w_write
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
integer, parameter unitcollapse
double precision, dimension(:,:), allocatable dx
integer nghostcells
Number of ghost cells surrounding a grid.
integer collapselevel
The level at which to produce line-integrated / collapsed output.
logical slab_uniform
uniform Cartesian geometry or not (stretched Cartesian)
character(len=std_len) base_filename
Base file name for simulation output, which will be followed by a number.
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition: mod_physics.t:4
procedure(sub_convert), pointer phys_to_primitive
Definition: mod_physics.t:59
Writes D-1 slice, can do so in various formats, depending on slice_type.
Definition: mod_slice.t:2
subroutine alloc_subnode(jgrid, dir, nwexpand)
Definition: mod_slice.t:994
subroutine dealloc_subnode(jgrid)
Definition: mod_slice.t:1042
Module with all the methods that users can customize in AMRVAC.
procedure(aux_output), pointer usr_aux_output
Pointer to a tree_node.
Definition: mod_forest.t:6