Line data Source code
1 : #include "Constant.H"
2 :
3 : namespace BC
4 : {
5 :
6 :
7 0 : Constant::Constant (int a_ncomp,
8 : amrex::Vector<std::string> bc_hi_str,
9 : amrex::Vector<std::string> bc_lo_str,
10 : AMREX_D_DECL(amrex::Vector<amrex::Real> _bc_lo_1,
11 : amrex::Vector<amrex::Real> _bc_lo_2,
12 : amrex::Vector<amrex::Real> _bc_lo_3),
13 : AMREX_D_DECL(amrex::Vector<amrex::Real> _bc_hi_1,
14 : amrex::Vector<amrex::Real> _bc_hi_2,
15 0 : amrex::Vector<amrex::Real> _bc_hi_3))
16 : //:
17 : //AMREX_D_DECL(bc_lo_1(_bc_lo_1),bc_lo_2(_bc_lo_2),bc_lo_3(_bc_lo_3)),
18 : //AMREX_D_DECL(bc_hi_1(_bc_hi_1),bc_hi_2(_bc_hi_2),bc_hi_3(_bc_hi_3))
19 : {
20 0 : Util::Warning(INFO,"This method is going away. Please use pp.queryclass() instead.");
21 :
22 0 : m_ncomp = a_ncomp;
23 :
24 0 : m_bc_type[Face::XLO].resize(m_ncomp,BCUtil::ReadString(bc_lo_str[0]));
25 0 : m_bc_type[Face::XHI].resize(m_ncomp,BCUtil::ReadString(bc_hi_str[0]));
26 0 : m_bc_type[Face::YLO].resize(m_ncomp,BCUtil::ReadString(bc_lo_str[1]));
27 0 : m_bc_type[Face::YHI].resize(m_ncomp,BCUtil::ReadString(bc_hi_str[1]));
28 : #if AMREX_SPACEDIM == 3
29 0 : m_bc_type[Face::ZLO].resize(m_ncomp,BCUtil::ReadString(bc_lo_str[2]));
30 0 : m_bc_type[Face::ZHI].resize(m_ncomp,BCUtil::ReadString(bc_hi_str[2]));
31 : #endif
32 :
33 :
34 0 : m_bc_val[Face::XLO].resize(m_ncomp,NAN);
35 0 : m_bc_val[Face::XHI].resize(m_ncomp,NAN);
36 0 : m_bc_val[Face::YLO].resize(m_ncomp,NAN);
37 0 : m_bc_val[Face::YHI].resize(m_ncomp,NAN);
38 : #if AMREX_SPACEDIM == 3
39 0 : m_bc_val[Face::ZLO].resize(m_ncomp,NAN);
40 0 : m_bc_val[Face::ZHI].resize(m_ncomp,NAN);
41 : #endif
42 :
43 0 : for (unsigned int i=0;i<m_ncomp;i++)
44 : {
45 0 : if (_bc_lo_1.size() > 0) m_bc_val[Face::XLO][i] = _bc_lo_1[i];
46 0 : if (_bc_hi_1.size() > 0) m_bc_val[Face::XHI][i] = _bc_hi_1[i];
47 0 : if (_bc_lo_2.size() > 0) m_bc_val[Face::YLO][i] = _bc_lo_2[i];
48 0 : if (_bc_hi_2.size() > 0) m_bc_val[Face::YHI][i] = _bc_hi_2[i];
49 : #if AMREX_SPACEDIM == 3
50 0 : if (_bc_lo_3.size() > 0) m_bc_val[Face::ZLO][i] = _bc_lo_3[i];
51 0 : if (_bc_hi_3.size() > 0) m_bc_val[Face::ZHI][i] = _bc_hi_3[i];
52 : #endif
53 : }
54 0 : }
55 :
56 :
57 : //amrex::Mask& m
58 : void
59 1545226 : Constant::FillBoundary (amrex::BaseFab<Set::Scalar> &a_in,
60 : const amrex::Box &a_box,
61 : int ngrow, int /*dcomp*/, int /*ncomp*/, amrex::Real time,
62 : Orientation face, const amrex::Mask * /*mask*/)
63 : {
64 1545226 : const amrex::Real* DX = m_geom.CellSize();
65 :
66 3090442 : Util::Assert(INFO,TEST(a_in.nComp() == (int)m_ncomp));
67 :
68 1545226 : amrex::Box box = a_box;
69 1545226 : box.grow(ngrow);
70 4635668 : const amrex::Dim3 lo= amrex::lbound(m_geom.Domain()), hi = amrex::ubound(m_geom.Domain());
71 :
72 1545226 : amrex::Array4<amrex::Real> const& in = a_in.array();
73 :
74 3748252 : for (int n = 0; n < a_in.nComp(); n++)
75 2203026 : amrex::ParallelFor (box,[=] AMREX_GPU_DEVICE(int i, int j, int k)
76 : {
77 942492470 : amrex::IntVect glevel;
78 2835953900 : AMREX_D_TERM(glevel[0] = std::max(std::min(0,i-lo.x),i-hi.x); ,
79 : glevel[1] = std::max(std::min(0,j-lo.y),j-hi.y); ,
80 : glevel[2] = std::max(std::min(0,k-lo.z),k-hi.z); );
81 :
82 942492470 : if (glevel[0]<0 && (face == Orientation::xlo || face == Orientation::All)) // Left boundary
83 : {
84 56210464 : if (BCUtil::IsDirichlet(m_bc_type[Face::XLO][n]))
85 89774230 : in(i,j,k,n) = m_bc_val[Face::XLO][n](time);
86 11323300 : else if(BCUtil::IsNeumann(m_bc_type[Face::XLO][n]))
87 44115900 : in(i,j,k,n) = in(i-glevel[0],j,k,n) - (m_bc_val[Face::XLO].size() > 0 ? m_bc_val[Face::XLO][n](time)*DX[0] : 0);
88 294360 : else if(BCUtil::IsReflectEven(m_bc_type[Face::XLO][n]))
89 0 : in(i,j,k,n) = in(1-glevel[0],j,k,n);
90 294360 : else if(BCUtil::IsReflectOdd(m_bc_type[Face::XLO][n]))
91 0 : in(i,j,k,n) = -in(1-glevel[0],j,k,n);
92 294360 : else if(BCUtil::IsPeriodic(m_bc_type[Face::XLO][n])) {}
93 : else
94 0 : Util::Abort(INFO, "Incorrect boundary conditions");
95 : }
96 886282210 : else if (glevel[0]>0 && (face == Orientation::xhi || face == Orientation::All)) // Right boundary
97 : {
98 1736668 : if (BCUtil::IsDirichlet(m_bc_type[Face::XHI][n]))
99 2505966 : in(i,j,k,n) = m_bc_val[Face::XHI][n](time);
100 483680 : else if(BCUtil::IsNeumann(m_bc_type[Face::XHI][n]))
101 1483520 : in(i,j,k,n) = in(i-glevel[0],j,k,n) - (m_bc_val[Face::XHI].size() > 0 ? m_bc_val[Face::XHI][n](time)*DX[0] : 0);
102 112800 : else if(BCUtil::IsReflectEven(m_bc_type[Face::XHI][n]))
103 0 : in(i,j,k,n) = in(hi.x-glevel[0],j,k,n);
104 112800 : else if(BCUtil::IsReflectOdd(m_bc_type[Face::XHI][n]))
105 0 : in(i,j,k,n) = -in(hi.x-glevel[0],j,k,n);
106 112800 : else if(BCUtil::IsPeriodic(m_bc_type[Face::XHI][n])) {}
107 : else
108 0 : Util::Abort(INFO, "Incorrect boundary conditions");
109 : }
110 :
111 884545640 : else if (glevel[1]<0 && (face == Orientation::ylo || face == Orientation::All)) // Bottom boundary
112 : {
113 11683412 : if (BCUtil::IsDirichlet(m_bc_type[Face::YLO][n]))
114 4344760 : in(i,j,k,n) = m_bc_val[Face::YLO][n](time);
115 9511000 : else if (BCUtil::IsNeumann(m_bc_type[Face::YLO][n]))
116 36152500 : in(i,j,k,n) = in(i,j-glevel[1],k,n) - (m_bc_val[Face::YLO].size() > 0 ? m_bc_val[Face::YLO][n](time)*DX[1] : 0);
117 472862 : else if (BCUtil::IsReflectEven(m_bc_type[Face::YLO][n]))
118 0 : in(i,j,k,n) = in(i,j-glevel[1],k,n);
119 472862 : else if (BCUtil::IsReflectOdd(m_bc_type[Face::YLO][n]))
120 0 : in(i,j,k,n) = -in(i,j-glevel[1],k,n);
121 472862 : else if(BCUtil::IsPeriodic(m_bc_type[Face::YLO][n])) {}
122 : else
123 0 : Util::Abort(INFO, "Incorrect boundary conditions");
124 : }
125 872862730 : else if (glevel[1]>0 && (face == Orientation::yhi || face == Orientation::All)) // Top boundary
126 : {
127 8324760 : if (BCUtil::IsDirichlet(m_bc_type[Face::YHI][n]))
128 2117170 : in(i,j,k,n) = m_bc_val[Face::YHI][n](time);
129 7266180 : else if (BCUtil::IsNeumann(m_bc_type[Face::YHI][n]))
130 28114100 : in(i,j,k,n) = in(i,j-glevel[1],k,n) - (m_bc_val[Face::YHI].size() > 0 ? m_bc_val[Face::YHI][n](time)*DX[1] : 0);
131 237662 : else if (BCUtil::IsReflectEven(m_bc_type[Face::YHI][n]))
132 0 : in(i,j,k,n) = in(i,hi.y-glevel[1],k,n);
133 237662 : else if (BCUtil::IsReflectOdd(m_bc_type[Face::YHI][n]))
134 0 : in(i,j,k,n) = -in(i,hi.y-glevel[1],k,n);
135 237662 : else if(BCUtil::IsPeriodic(m_bc_type[Face::YHI][n])) {}
136 : else
137 0 : Util::Abort(INFO, "Incorrect boundary conditions");
138 : }
139 :
140 : #if AMREX_SPACEDIM>2
141 6658050 : else if (glevel[2]<0 && (face == Orientation::zlo || face == Orientation::All))
142 : {
143 600256 : if (BCUtil::IsDirichlet(m_bc_type[Face::ZLO][n]))
144 1200510 : in(i,j,k,n) = m_bc_val[Face::ZLO][n](time);
145 0 : else if (BCUtil::IsNeumann(m_bc_type[Face::ZLO][n]))
146 0 : in(i,j,k,n) = in(i,j,k-glevel[2],n) - (m_bc_val[Face::ZLO].size() > 0 ? m_bc_val[Face::ZLO][n](time)*DX[2] : 0);
147 0 : else if (BCUtil::IsReflectEven(m_bc_type[Face::ZLO][n]))
148 0 : in(i,j,k,n) = in(i,j,1-glevel[2],n);
149 0 : else if (BCUtil::IsReflectOdd(m_bc_type[Face::ZLO][n]))
150 0 : in(i,j,k,n) = -in(i,j,1-glevel[2],n);
151 0 : else if(BCUtil::IsPeriodic(m_bc_type[Face::ZLO][n])) {}
152 0 : else Util::Abort(INFO, "Incorrect boundary conditions");
153 : }
154 6057790 : else if (glevel[2]>0 && (face == Orientation::zhi || face == Orientation::All))
155 : {
156 102144 : if (BCUtil::IsDirichlet(m_bc_type[Face::ZHI][n]))
157 204288 : in(i,j,k,n) = m_bc_val[Face::ZHI][n](time);
158 0 : else if(BCUtil::IsNeumann(m_bc_type[Face::ZHI][n]))
159 0 : in(i,j,k,n) = in(i,j,k-glevel[2],n) - (m_bc_val[Face::ZHI].size() > 0 ? m_bc_val[Face::ZHI][n](time)*DX[2] : 0);
160 0 : else if(BCUtil::IsReflectEven(m_bc_type[Face::ZHI][n]))
161 0 : in(i,j,k,n) = in(i,j,hi.z-glevel[2],n);
162 0 : else if(BCUtil::IsReflectOdd(m_bc_type[Face::ZHI][n]))
163 0 : in(i,j,k,n) = -in(i,j,hi.z-glevel[2],n);
164 0 : else if(BCUtil::IsPeriodic(m_bc_type[Face::ZHI][n])) {}
165 0 : else Util::Abort(INFO, "Incorrect boundary conditions");
166 : }
167 : #endif
168 :
169 :
170 942492470 : });
171 1545226 : }
172 :
173 : amrex::BCRec
174 159336 : Constant::GetBCRec()
175 : {
176 159336 : int bc_lo[BL_SPACEDIM] = {AMREX_D_DECL(m_bc_type[Face::XLO][0],m_bc_type[Face::YLO][0],m_bc_type[Face::XLO][0])};
177 159336 : int bc_hi[BL_SPACEDIM] = {AMREX_D_DECL(m_bc_type[Face::XHI][0],m_bc_type[Face::YHI][0],m_bc_type[Face::XHI][0])};
178 :
179 159336 : return amrex::BCRec(bc_lo,bc_hi);
180 : }
181 :
182 : amrex::Array<int,AMREX_SPACEDIM>
183 0 : Constant::IsPeriodic()
184 : {
185 0 : return {AMREX_D_DECL(BCUtil::IsPeriodic(m_bc_type[Face::XLO][0]),
186 : BCUtil::IsPeriodic(m_bc_type[Face::YLO][0]),
187 0 : BCUtil::IsPeriodic(m_bc_type[Face::ZLO][0]))};
188 : }
189 0 : amrex::Periodicity Constant::Periodicity () const
190 : {
191 0 : return amrex::Periodicity(amrex::IntVect(AMREX_D_DECL(m_geom.Domain().length(0) * BCUtil::IsPeriodic(m_bc_type[Face::XLO][0]),
192 : m_geom.Domain().length(1) * BCUtil::IsPeriodic(m_bc_type[Face::YLO][0]),
193 0 : m_geom.Domain().length(2) * BCUtil::IsPeriodic(m_bc_type[Face::ZLO][0]))));
194 : }
195 0 : amrex::Periodicity Constant::Periodicity (const amrex::Box& b) {
196 0 : return amrex::Periodicity(amrex::IntVect(AMREX_D_DECL(b.length(0) * BCUtil::IsPeriodic(m_bc_type[Face::XLO][0]),
197 : b.length(1) * BCUtil::IsPeriodic(m_bc_type[Face::YLO][0]),
198 0 : b.length(2) * BCUtil::IsPeriodic(m_bc_type[Face::ZLO][0]))));
199 :
200 : }
201 :
202 :
203 : }
|
