aboutsummaryrefslogtreecommitdiff
path: root/Source/FieldSolver/WarpXPushFieldsEM_K.H
blob: e8e675cd276ca01c988c0cf623564233ade388a8 (plain) (blame)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
/* Copyright 2019-2020
 *
 * This file is part of WarpX.
 *
 * License: BSD-3-Clause-LBNL
 */

#ifndef WarpXPushFieldsEM_K_h
#define WarpXPushFieldsEM_K_h

#include "Utils/WarpXConst.H"

#include <AMReX.H>

/*
 * \brief Return a tilebox that only covers the outer half of the guard cells.
 *        For tileboxes that don't include cells beyond the whole domain,
 *        an empty box is returned.
 *
 * \param[in,out] input_tb tilebox to be modified
 * \param[in] dir direction where the tilebox smallEnd/bigEnd is modified
 * \param[in] n_domain number of valid cells in the whole simulation domain
 * \param[in] tb_smallEnd the lowest index of the tilebox, including guard cells
 * \param[in] tb_bigEnd the highest index of the tilebox, including guard cells
 */
AMREX_GPU_HOST_DEVICE AMREX_INLINE
amrex::Box constrain_tilebox_to_guards(
    const amrex::Box& input_tb,
    const int dir,
    const int iside,
    const int n_domain,
    const int tb_smallEnd,
    const int tb_bigEnd)
{
    using namespace amrex;

    amrex::Box constrained_tb = input_tb;

    // If the input_tb does not overlap either the lower or upper guard,
    // an empty box is returned.

    if (iside == 0)
    {
        // Lower guard
        const int n_guard = -tb_smallEnd;
        const int upper_bound = (tb_smallEnd < 0 ? -n_guard/2 : tb_smallEnd);
        constrained_tb.setBig(dir, upper_bound - 1);
    }
    else if (iside == 1)
    {
        // Upper guard
        const int n_guard = tb_bigEnd - n_domain;
        const int lower_bound = (tb_bigEnd > n_domain ? n_domain + n_guard/2 : tb_bigEnd);
        constrained_tb.setSmall(dir, lower_bound + 1);
    }

    return constrained_tb;
}

/*
 * \brief Damp a given field in the guard cells along a given direction
 *
 * \param[in,out] mf_arr array that contains the field values to be damped
 * \oaram[in] i index along x
 * \oaram[in] j index along y (in 3D) or z (in 2D/RZ)
 * \oaram[in] k index along z (in 3D, \c k = 0 in 2D/RZ)
 * \param[in] icomp index along the fourth component of the array
 * \param]in] dir direction where the field will be damped
 * \param[in] n_domain number of valid cells in the whole simulation domain
 * \param[in] tb_smallEnd the lowest index of the tilebox, including guard cells
 * \param[in] tb_bigEnd the highest index of the tilebox, including guard cells
 */
AMREX_GPU_HOST_DEVICE AMREX_INLINE
void damp_field_in_guards(
    amrex::Array4<amrex::Real> const& mf_arr,
    const int i,
    const int j,
    const int k,
    const int icomp,
    const int dir,
    const int n_domain,
    const int tb_smallEnd,
    const int tb_bigEnd)
{
    using namespace amrex;

    // dir = 0: idx = i (x)
    // dir = 1: idx = j (y in 3D, z in 2D/RZ)
    // dir = 2: idx = k (z in 3D)
    const int idx = ((dir == 0) ? i : ((dir == 1) ? j : k));

    if (idx < 0)
    {
        // Apply damping factor in guards cells below the lower end of the domain
        const int n_guard = -tb_smallEnd;

        const amrex::Real cell = static_cast<amrex::Real>(idx + n_guard);

        const amrex::Real phase = MathConst::pi * cell / n_guard;
        const amrex::Real sin_phase = std::sin(phase);
        const amrex::Real damp_factor = sin_phase * sin_phase;

        mf_arr(i,j,k,icomp) *= damp_factor;
    }
    else if (idx > n_domain)
    {
        // Apply damping factor in guards cells above the upper end of the domain
        const int n_guard = tb_bigEnd - n_domain;

        const amrex::Real cell = static_cast<amrex::Real>(tb_bigEnd - idx);

        const amrex::Real phase = MathConst::pi * cell / n_guard;
        const amrex::Real sin_phase = std::sin(phase);
        const amrex::Real damp_factor = sin_phase * sin_phase;

        mf_arr(i,j,k,icomp) *= damp_factor;
    }
}

#endif