MPI-AMRVAC  3.1
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  !> trace small values in the source file using traceback flag of compiler
15  logical, public :: trace_small_values=.false.
16 
17  !> Whether to apply small value fixes to certain variables
18  logical, public, allocatable :: small_values_fix_iw(:)
19 
20  public :: small_values_error
21  public :: small_values_average
22 
23 contains
24 
25  subroutine small_values_error(wprim, x, ixI^L, ixO^L, w_flag, subname)
27  integer, intent(in) :: ixi^l, ixo^l
28  double precision, intent(in) :: wprim(ixi^s, 1:nw)
29  double precision, intent(in) :: x(ixi^s, 1:ndim)
30  logical, intent(in) :: w_flag(ixi^s,1:nw)
31  character(len=*), intent(in) :: subname
32  integer :: iw,iiw,ix^d
33 
34  if (.not.crash) then
35  do iw=1,nw
36  {do ix^db= ixo^lim^db\}
37  if(w_flag(ix^d,iw)) then
38  write(*,*) "Error: small value of ", trim(prim_wnames(iw)),wprim(ix^d,iw),&
39  " encountered when call ", subname
40  write(*,*) "Iteration: ", it, " Time: ", global_time, "Processor: ",mype
41  write(*,*) "Location: ", x(ix^d,:)
42  write(*,*) "Cell number: ", ix^d
43  do iiw=1,nw
44  write(*,*) trim(prim_wnames(iiw)),": ",wprim(ix^d,iiw)
45  end do
46  ! use erroneous arithmetic operation to crash the run
47  if(trace_small_values) write(*,*) sqrt(wprim(ix^d,iw)-bigdouble)
48  write(*,*) "Saving status at the previous time step"
49  crash=.true.
50  end if
51  {enddo^d&\}
52  end do
53  end if
54 
55  end subroutine small_values_error
56 
57  subroutine small_values_average(ixI^L, ixO^L, w, x, w_flag, windex)
59  integer, intent(in) :: ixi^l, ixo^l
60  logical, intent(in) :: w_flag(ixi^s,1:nw)
61  double precision, intent(inout) :: w(ixi^s, 1:nw)
62  double precision, intent(in) :: x(ixi^s, 1:ndim)
63  integer, optional, intent(in) :: windex
64  integer :: iw, kxo^l, ix^d, i, nwstart, nwend
65 
66  if(present(windex)) then
67  nwstart=windex
68  nwend=windex
69  else
70  nwstart=1
71  nwend=nw
72  end if
73 
74  do iw = nwstart, nwend
75  {do ix^db= ixo^lim^db\}
76  ! point with local failure identified by w_flag
77  if (w_flag(ix^d,iw)) then
78  ! verify in cube with border width small_values_daverage the presence of
79  ! cells where all went ok
80  do i = 1, max(small_values_daverage, 1)
81  {kxomin^d= max(ix^d-i, iximin^d);
82  kxomax^d= min(ix^d+i, iximax^d);\}
83  ! in case cells are fine within smaller cube than
84  ! the userset small_values_daverage: use that smaller cube
85  if(any(w_flag(kxo^s,iw) .eqv. .false.)) exit
86  end do
87 
88  if(any(w_flag(kxo^s,iw) .eqv. .false.)) then
89  ! within surrounding cube, cells without problem were found
90 
91  ! faulty cells are corrected by averaging here
92  ! only average those which were ok and replace faulty cells
93  if(small_values_fix_iw(iw)) then
94  w(ix^d, iw) = sum(w(kxo^s, iw), w_flag(kxo^s,iw) .eqv. .false.)&
95  / count(w_flag(kxo^s,iw) .eqv. .false.)
96  end if
97  else
98  write(*,*) "no cells without error were found in cube of size", &
100  write(*,*) "at location:", x(ix^d, 1:ndim)
101  write(*,*) "at index:", ix^d
102  write(*,*) "w numer:", iw
103  !write(*,*) "Saving status at the previous time step"
104  !crash=.true.
105  write(*,*) "replace with small values"
106  if(iw==iw_e) w(ix^d, iw)=small_pressure
107  if(iw==iw_rho) w(ix^d, iw)=small_density
108  end if
109  end if
110  {enddo^d&\}
111  end do
112 
113  end subroutine small_values_average
114 
115 end module mod_small_values
This module contains definitions of global parameters and variables and some generic functions/subrou...
double precision small_pressure
double precision global_time
The global simulation time.
integer it
Number of time steps taken.
integer, parameter ndim
Number of spatial dimensions for grid variables.
integer mype
The rank of the current MPI task.
integer, dimension(:), allocatable, parameter d
logical crash
Save a snapshot before crash a run met unphysical values.
double precision small_density
Module for handling problematic values in simulations, such as negative pressures.
integer, public small_values_daverage
Average over this many cells in each direction.
logical, public trace_small_values
trace small values in the source file using traceback flag of compiler
subroutine, public small_values_average(ixIL, ixOL, w, x, w_flag, windex)
logical, dimension(:), allocatable, public small_values_fix_iw
Whether to apply small value fixes to certain variables.
subroutine, public small_values_error(wprim, x, ixIL, ixOL, w_flag, subname)
character(len=20), public small_values_method
How to handle small values.