Advection¶
\[-
\dtempadv{1}\]
92static int advection_x (
93 const domain_t * domain,
94 const double * restrict t,
95 const double * restrict ux,
96 double * restrict src
97) {
98 const int isize = domain->mysizes[0];
99 const int jsize = domain->mysizes[1];
100#if NDIMS == 3
101 const int ksize = domain->mysizes[2];
102#endif
103 const double * restrict hxxf = domain->hxxf;
104 const double * restrict jdxf = domain->jdxf;
105 const double * restrict jdxc = domain->jdxc;
106 BEGIN
107 const double hx_xm = HXXF(i );
108 const double hx_xp = HXXF(i+1);
109 const double jd_xm = JDXF(i );
110 const double jd_x0 = JDXC(i );
111 const double jd_xp = JDXF(i+1);
112#if NDIMS == 2
113 const double ux_xm = jd_xm / hx_xm * UX(i , j );
114 const double ux_xp = jd_xp / hx_xp * UX(i+1, j );
115#else
116 const double ux_xm = jd_xm / hx_xm * UX(i , j , k );
117 const double ux_xp = jd_xp / hx_xp * UX(i+1, j , k );
118#endif
119 const double l = - 0.5 * ux_xm;
120 const double u = + 0.5 * ux_xp;
121 const double c = - l - u;
122 src[cnt] -= 1. / jd_x0 * (
123#if NDIMS == 2
124 + l * T(i-1, j )
125 + c * T(i , j )
126 + u * T(i+1, j )
127#else
128 + l * T(i-1, j , k )
129 + c * T(i , j , k )
130 + u * T(i+1, j , k )
131#endif
132 );
133 END
134 return 0;
135}
\[-
\dtempadv{2}\]
138static int advection_y (
139 const domain_t * domain,
140 const double * restrict t,
141 const double * restrict uy,
142 double * restrict src
143) {
144 const int isize = domain->mysizes[0];
145 const int jsize = domain->mysizes[1];
146#if NDIMS == 3
147 const int ksize = domain->mysizes[2];
148#endif
149 const double hy = domain->hy;
150 const double * restrict jdxc = domain->jdxc;
151 BEGIN
152 const double jd = JDXC(i );
153#if NDIMS == 2
154 const double uy_ym = jd / hy * UY(i , j );
155 const double uy_yp = jd / hy * UY(i , j+1);
156#else
157 const double uy_ym = jd / hy * UY(i , j , k );
158 const double uy_yp = jd / hy * UY(i , j+1, k );
159#endif
160 const double l = - 0.5 * uy_ym;
161 const double u = + 0.5 * uy_yp;
162 const double c = - l - u;
163 src[cnt] -= 1. / jd * (
164#if NDIMS == 2
165 + l * T(i , j-1)
166 + c * T(i , j )
167 + u * T(i , j+1)
168#else
169 + l * T(i , j-1, k )
170 + c * T(i , j , k )
171 + u * T(i , j+1, k )
172#endif
173 );
174 END
175 return 0;
176}
\[-
\dtempadv{3}\]
180static int advection_z (
181 const domain_t * domain,
182 const double * restrict t,
183 const double * restrict uz,
184 double * restrict src
185) {
186 const int isize = domain->mysizes[0];
187 const int jsize = domain->mysizes[1];
188 const int ksize = domain->mysizes[2];
189 const double hz = domain->hz;
190 const double * restrict jdxc = domain->jdxc;
191 BEGIN
192 const double jd = JDXC(i );
193 const double uz_zm = jd / hz * UZ(i , j , k );
194 const double uz_zp = jd / hz * UZ(i , j , k+1);
195 const double l = - 0.5 * uz_zm;
196 const double u = + 0.5 * uz_zp;
197 const double c = - l - u;
198 src[cnt] -= 1. / jd * (
199 + l * T(i , j , k-1)
200 + c * T(i , j , k )
201 + u * T(i , j , k+1)
202 );
203 END
204 return 0;
205}