MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
mod_small_values.t
Go to the documentation of this file.
1 !> Module for handling problematic values in simulations, such as negative
2 !> pressures
4 
5  implicit none
6  private
7 
8  !> How to handle small values
9  character(len=20), public :: small_values_method = "error"
10 
11  !> Average over this many cells in each direction
12  integer, public :: small_values_daverage = 1
13 
14  public :: small_values_error
15  public :: small_values_average
16 
17 contains
18 
19  subroutine small_values_error(w, x, ixI^L, ixO^L, w_flag, subname)
21  integer, intent(in) :: ixI^L, ixO^L
22  double precision, intent(in) :: w(ixi^s, 1:nw)
23  double precision, intent(in) :: x(ixi^s, 1:ndim)
24  integer, intent(in) :: w_flag(ixi^s)
25  integer :: ix_bad(ndim), iw
26  character(len=*), intent(in) :: subname
27 
28  ix_bad = maxloc(w_flag(ixo^s)) + [ ixomin^d-1 ]
29 
30  if (.not.crash) then
31  write(*,*) "Error: small value of ", trim(prim_wnames(maxval(w_flag(ixo^s)))), &
32  " encountered when call ", subname
33  write(*,*) "Iteration: ", it, " Time: ", global_time
34  write(*,*) "Location: ", x({ix_bad(^d)}, :)
35  write(*,*) "Cell number: ", ix_bad(:)
36  do iw = 1, nw
37  write(*, '(A20,A,E14.7)') trim(cons_wnames(iw)), ": ", &
38  w({ix_bad(^d)}, iw)
39  end do
40  write(*,*) "Saving status at the previous time step"
41  crash=.true.
42  end if
43  end subroutine small_values_error
44 
45  subroutine small_values_average(ixI^L, ixO^L, w, x, w_flag)
47  integer, intent(in) :: ixI^L, ixO^L
48  integer, intent(in) :: w_flag(ixi^s)
49  double precision, intent(inout) :: w(ixi^s, 1:nw)
50  double precision, intent(in) :: x(ixi^s, 1:ndim)
51  integer :: iw, kxO^L, ix^D, i
52 
53  {do ix^db= ixo^lim^db\}
54 
55  ! point with local failure identified by w_flag
56  if (w_flag(ix^d) /= 0) then
57  ! verify in cube with border width small_values_daverage the presence of
58  ! cells where all went ok
59  do i = 1, max(small_values_daverage, 1)
60  {kxomin^d= max(ix^d-i, ixomin^d);
61  kxomax^d= min(ix^d+i, ixomax^d);\}
62 
63  ! in case cells are fine within smaller cube than
64  ! the userset small_values_daverage: use that smaller cube
65  if (any(w_flag(kxo^s) == 0)) exit
66  end do
67 
68  if (any(w_flag(kxo^s) == 0)) then
69  ! within surrounding cube, cells without problem were found
70 
71  ! faulty cells are corrected by averaging here
72  ! only average those which were ok and replace faulty cells
73  do iw = 1, nw
74  if (small_values_fix_iw(iw)) then
75  w(ix^d, iw) = sum(w(kxo^s, iw), w_flag(kxo^s) == 0)&
76  / count(w_flag(kxo^s) == 0)
77  end if
78  end do
79  else
80  write(*,*) "no cells without error were found in cube of size", &
82  write(*,*) "at location:", x(ix^d, 1:ndim)
83  write(*,*) "at index:", ix^d
84  write(*,*) "w_flag(ix^D):", w_flag(ix^d)
85  write(*,*) "Saving status at the previous time step"
86  crash=.true.
87  end if
88  end if
89  {enddo^d&\}
90 
91  end subroutine small_values_average
92 
93 end module mod_small_values
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, public small_values_daverage
Average over this many cells in each direction.
subroutine, public small_values_average(ixIL, ixOL, w, x, w_flag)
subroutine, public small_values_error(w, x, ixIL, ixOL, w_flag, subname)
logical, dimension(:), allocatable small_values_fix_iw
Whether to apply small value fixes to certain variables.
integer, dimension(:), allocatable, parameter d
double precision global_time
The global simulation time.
Module for handling problematic values in simulations, such as negative pressures.
integer, parameter ndim
Number of spatial dimensions for grid variables.
logical crash
Save a snapshot before crash a run met unphysical values.
integer it
Number of time steps taken.
character(len=20), public small_values_method
How to handle small values.