1 #ifndef AMREX_MLEB_LEASTSQUARES_K_H_
2 #define AMREX_MLEB_LEASTSQUARES_K_H_
3 #include <AMReX_Config.H>
18 for (
int ii = 0; ii < neq; ii++)
22 for (
int jj = ii; jj < neq; jj++)
26 for (
int kk = ii-1; kk >= 0; kk--)
28 sum1 = sum1 - aa(ii,kk)*aa(jj,kk);
42 aa(jj,ii) = sum1 / p[ii];
50 for (
int ii = 0; ii < neq; ii++)
52 for (
int jj = ii+1; jj < neq; jj++)
55 aa(ii,jj) = aa(jj,ii);
70 for (
int irow = 0; irow < neq; irow++) {
71 for (
int icol = 0; icol < neq; icol++) {
76 for (
int irow = 0; irow < 12; irow++)
78 AtA(0,0) += Amatrix(irow,0)*Amatrix(irow,0);
79 AtA(0,1) += Amatrix(irow,0)*Amatrix(irow,1);
80 AtA(0,2) += Amatrix(irow,0)*Amatrix(irow,2);
81 AtA(0,3) += Amatrix(irow,0)*Amatrix(irow,3);
82 AtA(0,4) += Amatrix(irow,0)*Amatrix(irow,4);
83 AtA(0,5) += Amatrix(irow,0)*Amatrix(irow,5);
84 AtA(1,1) += Amatrix(irow,1)*Amatrix(irow,1);
85 AtA(1,2) += Amatrix(irow,1)*Amatrix(irow,2);
86 AtA(1,3) += Amatrix(irow,1)*Amatrix(irow,3);
87 AtA(1,4) += Amatrix(irow,1)*Amatrix(irow,4);
88 AtA(1,5) += Amatrix(irow,1)*Amatrix(irow,5);
89 AtA(2,2) += Amatrix(irow,2)*Amatrix(irow,2);
90 AtA(2,3) += Amatrix(irow,2)*Amatrix(irow,3);
91 AtA(2,4) += Amatrix(irow,2)*Amatrix(irow,4);
92 AtA(2,5) += Amatrix(irow,2)*Amatrix(irow,5);
93 AtA(3,3) += Amatrix(irow,3)*Amatrix(irow,3);
94 AtA(3,4) += Amatrix(irow,3)*Amatrix(irow,4);
95 AtA(3,5) += Amatrix(irow,3)*Amatrix(irow,5);
96 AtA(4,4) += Amatrix(irow,4)*Amatrix(irow,4);
97 AtA(4,5) += Amatrix(irow,4)*Amatrix(irow,5);
98 AtA(5,5) += Amatrix(irow,5)*Amatrix(irow,5);
101 for (
int irow = 0; irow < neq-1; irow++) {
102 for (
int icol = irow+1; icol < neq; icol++) {
103 AtA(icol,irow) = AtA(irow,icol);
110 b(0) =
b(0) / AtA(0,0);
115 for (
int ii = 1; ii < neq; ii++)
119 for (
int jj = 0; jj < ii; jj++) {
120 b(ii) =
b(ii) - AtA(ii,jj)*
b(jj);
123 b(ii) =
b(ii) / AtA(ii,ii);
131 if (AtA(neq-1,neq-1) > 0.) {
132 b(neq-1) =
b(neq-1) / AtA(neq-1,neq-1);
137 for (
int ii = neq-2; ii >= 0; ii--)
141 for (
int jj = ii+1; jj < neq; jj++) {
142 b(ii) =
b(ii) - AtA(ii,jj)*
b(jj);
145 b(ii) =
b(ii) / AtA(ii,ii);
161 for (
int irow = 0; irow < neq; irow++) {
162 for (
int icol = 0; icol < neq; icol++) {
163 AtA(irow,icol) = 0.0;
167 for (
int irow = 0; irow < 18; irow++)
169 AtA(0,0) += Amatrix(irow,0)*Amatrix(irow,0);
170 AtA(0,1) += Amatrix(irow,0)*Amatrix(irow,1);
171 AtA(0,2) += Amatrix(irow,0)*Amatrix(irow,2);
172 AtA(0,3) += Amatrix(irow,0)*Amatrix(irow,3);
173 AtA(0,4) += Amatrix(irow,0)*Amatrix(irow,4);
174 AtA(0,5) += Amatrix(irow,0)*Amatrix(irow,5);
175 AtA(1,1) += Amatrix(irow,1)*Amatrix(irow,1);
176 AtA(1,2) += Amatrix(irow,1)*Amatrix(irow,2);
177 AtA(1,3) += Amatrix(irow,1)*Amatrix(irow,3);
178 AtA(1,4) += Amatrix(irow,1)*Amatrix(irow,4);
179 AtA(1,5) += Amatrix(irow,1)*Amatrix(irow,5);
180 AtA(2,2) += Amatrix(irow,2)*Amatrix(irow,2);
181 AtA(2,3) += Amatrix(irow,2)*Amatrix(irow,3);
182 AtA(2,4) += Amatrix(irow,2)*Amatrix(irow,4);
183 AtA(2,5) += Amatrix(irow,2)*Amatrix(irow,5);
184 AtA(3,3) += Amatrix(irow,3)*Amatrix(irow,3);
185 AtA(3,4) += Amatrix(irow,3)*Amatrix(irow,4);
186 AtA(3,5) += Amatrix(irow,3)*Amatrix(irow,5);
187 AtA(4,4) += Amatrix(irow,4)*Amatrix(irow,4);
188 AtA(4,5) += Amatrix(irow,4)*Amatrix(irow,5);
189 AtA(5,5) += Amatrix(irow,5)*Amatrix(irow,5);
192 for (
int irow = 0; irow < neq-1; irow++) {
193 for (
int icol = irow+1; icol < neq; icol++) {
194 AtA(icol,irow) = AtA(irow,icol);
201 b(0) =
b(0) / AtA(0,0);
206 for (
int ii = 1; ii < neq; ii++)
210 for (
int jj = 0; jj < ii; jj++) {
211 b(ii) =
b(ii) - AtA(ii,jj)*
b(jj);
214 b(ii) =
b(ii) / AtA(ii,ii);
222 if (AtA(neq-1,neq-1) > 0.) {
223 b(neq-1) =
b(neq-1) / AtA(neq-1,neq-1);
228 for (
int ii = neq-2; ii >= 0; ii--)
232 for (
int jj = ii+1; jj < neq; jj++) {
233 b(ii) =
b(ii) - AtA(ii,jj)*
b(jj);
236 b(ii) =
b(ii) / AtA(ii,ii);
253 bool is_eb_dirichlet,
bool is_eb_inhomog)
262 for (
int irow = 0; irow < 12; irow++) {
263 for (
int icol = 0; icol < 6; icol++) {
264 Amatrix(irow,icol) = 0.0;
269 for (
int ii = i-1; ii <= i; ii++) {
270 for (
int jj = j-1; jj <= j+1; jj++)
272 if (!flag(ii,jj,k).isCovered())
274 int a_ind = (jj-(j-1)) + 3*(ii-(i-1));
278 Real x_off =
static_cast<Real
>(ii-i) + 0.5;
279 Real y_off =
static_cast<Real
>(jj-j);
281 Amatrix(a_ind,0) = 1.0;
282 Amatrix(a_ind,1) = x_off + ccent(ii,jj,k,0);
283 Amatrix(a_ind,2) = y_off + ccent(ii,jj,k,1) - yloc_on_xface;
284 Amatrix(a_ind,3) = Amatrix(a_ind,1) * Amatrix(a_ind,1);
285 Amatrix(a_ind,4) = Amatrix(a_ind,1) * Amatrix(a_ind,2);
286 Amatrix(a_ind,5) = Amatrix(a_ind,2) * Amatrix(a_ind,2);
288 if (!flag(ii,jj,k).isRegular())
290 int a_ind_eb = a_ind + 6;
292 Amatrix(a_ind_eb,0) = 1.0;
293 Amatrix(a_ind_eb,1) = x_off + bcent(ii,jj,k,0);
294 Amatrix(a_ind_eb,2) = y_off + bcent(ii,jj,k,1) - yloc_on_xface;
295 Amatrix(a_ind_eb,3) = Amatrix(a_ind_eb,1) * Amatrix(a_ind_eb,1);
296 Amatrix(a_ind_eb,4) = Amatrix(a_ind_eb,1) * Amatrix(a_ind_eb,2);
297 Amatrix(a_ind_eb,5) = Amatrix(a_ind_eb,2) * Amatrix(a_ind_eb,2);
304 for (
int irow = 0; irow < 6; irow++)
308 for (
int ii = i-1; ii <= i; ii++) {
309 for (
int jj = j-1; jj <= j+1; jj++)
311 if (!flag(ii,jj,k).isCovered())
313 int a_ind = (jj-(j-1)) + 3*(ii-(i-1));
314 rhs(irow) += Amatrix(a_ind ,irow)* phi(ii,jj,k,n);
316 if (flag(ii,jj,k).isSingleValued() &&
317 is_eb_dirichlet && is_eb_inhomog) {
318 rhs(irow) += Amatrix(a_ind+6,irow)*phieb(ii,jj,k,n);
338 bool is_eb_dirichlet,
bool is_eb_inhomog)
347 for (
int irow = 0; irow < 12; irow++) {
348 for (
int icol = 0; icol < 6; icol++) {
349 Amatrix(irow,icol) = 0.0;
354 for (
int jj = j-1; jj <= j; jj++) {
355 for (
int ii = i-1; ii <= i+1; ii++)
357 if (!flag(ii,jj,k).isCovered())
359 int a_ind = (ii-(i-1)) + 3*(jj-(j-1));
361 Real x_off =
static_cast<Real
>(ii-i);
362 Real y_off =
static_cast<Real
>(jj-j) + 0.5;
364 Amatrix(a_ind,0) = 1.0;
365 Amatrix(a_ind,1) = x_off + ccent(ii,jj,k,0) - xloc_on_yface;
366 Amatrix(a_ind,2) = y_off + ccent(ii,jj,k,1);
367 Amatrix(a_ind,3) = Amatrix(a_ind,1) * Amatrix(a_ind,1);
368 Amatrix(a_ind,4) = Amatrix(a_ind,1) * Amatrix(a_ind,2);
369 Amatrix(a_ind,5) = Amatrix(a_ind,2) * Amatrix(a_ind,2);
371 if (!flag(ii,jj,k).isRegular())
373 int a_ind_eb = a_ind + 6;
374 Amatrix(a_ind_eb,0) = 1.0;
375 Amatrix(a_ind_eb,1) = x_off + bcent(ii,jj,k,0) - xloc_on_yface;
376 Amatrix(a_ind_eb,2) = y_off + bcent(ii,jj,k,1);
377 Amatrix(a_ind_eb,3) = Amatrix(a_ind_eb,1) * Amatrix(a_ind_eb,1);
378 Amatrix(a_ind_eb,4) = Amatrix(a_ind_eb,1) * Amatrix(a_ind_eb,2);
379 Amatrix(a_ind_eb,5) = Amatrix(a_ind_eb,2) * Amatrix(a_ind_eb,2);
386 for (
int irow = 0; irow < 6; irow++)
390 for (
int jj = j-1; jj <= j; jj++) {
391 for (
int ii = i-1; ii <= i+1; ii++) {
392 if (!flag(ii,jj,k).isCovered())
394 int a_ind = (ii-(i-1)) + 3*(jj-(j-1));
395 rhs(irow) += Amatrix(a_ind ,irow)* phi(ii,jj,k,n);
397 if (flag(ii,jj,k).isSingleValued() &&
398 is_eb_dirichlet && is_eb_inhomog) {
399 rhs(irow) += Amatrix(a_ind+6,irow)*phieb(ii,jj,k,n);
418 Real& nrmx, Real& nrmy,
426 for (
int irow = 0; irow < 18; irow++) {
427 for (
int icol = 0; icol < 6; icol++) {
428 Amatrix(irow,icol) = 0.0;
433 for (
int ii = i-1; ii <= i+1; ii++) {
434 for (
int jj = j-1; jj <= j+1; jj++)
436 if (!flag(ii,jj,k).isCovered())
438 int a_ind = (jj-(j-1)) + 3*(ii-(i-1));
440 Real x_off =
static_cast<Real
>(ii-i);
441 Real y_off =
static_cast<Real
>(jj-j);
443 if (flag(i,j,k).isConnected((ii-i),(jj-j),0))
445 Amatrix(a_ind,0) = 1.0;
446 Amatrix(a_ind,1) = x_off + ccent(ii,jj,k,0) - bcent(i,j,k,0);
447 Amatrix(a_ind,2) = y_off + ccent(ii,jj,k,1) - bcent(i,j,k,1);
450 if (flag(i,j,k).isConnected((ii-i),(jj-j),0) && !flag(ii,jj,k).isRegular())
452 Amatrix(a_ind+9,0) = 1.0;
453 Amatrix(a_ind+9,1) = x_off + bcent(ii,jj,k,0) - bcent(i,j,k,0);
454 Amatrix(a_ind+9,2) = y_off + bcent(ii,jj,k,1) - bcent(i,j,k,1);
462 for (
int irow = 0; irow < 18; irow++)
464 Amatrix(irow,3) = Amatrix(irow,1) * Amatrix(irow,1);
465 Amatrix(irow,4) = Amatrix(irow,1) * Amatrix(irow,2);
466 Amatrix(irow,5) = Amatrix(irow,2) * Amatrix(irow,2);
470 for (
int irow = 0; irow < 6; irow++)
474 for (
int ii = i-1; ii <= i+1; ii++) {
475 for (
int jj = j-1; jj <= j+1; jj++) {
476 if (!flag(ii,jj,k).isCovered())
478 int a_ind = (jj-(j-1)) + 3*(ii-(i-1));
479 rhs(irow) += Amatrix(a_ind,irow) * phi(ii,jj,k,n);
480 if (flag(ii,jj,k).isSingleValued() && is_eb_inhomog) {
481 rhs(irow) += Amatrix(a_ind+9,irow)*phieb(ii,jj,k,n);
490 Real dphidn = rhs(1)*nrmx + rhs(2)*nrmy;
503 bool is_eb_dirichlet,
bool is_eb_inhomog,
504 const bool on_x_face,
const int domlo_x,
const int domhi_x,
505 const bool on_y_face,
const int domlo_y,
const int domhi_y)
514 for (
int irow = 0; irow < 18; irow++) {
515 for (
int icol = 0; icol < 6; icol++) {
516 Amatrix(irow,icol) = 0.0;
520 const int im = (i > domhi_x) ? 2 : 1;
521 const int ip = 2 - im;
524 for (
int ii = i-im; ii <= i+ip; ii++) {
525 for (
int jj = j-1; jj <= j+1; jj++)
530 if (((on_x_face && ii < domlo_x) && (on_y_face && jj < domlo_y)) ||
531 ((on_x_face && ii > domhi_x) && (on_y_face && jj > domhi_y)) ||
532 ((on_x_face && ii < domlo_x) && (on_y_face && jj > domhi_y)) ||
533 ((on_x_face && ii > domhi_x) && (on_y_face && jj < domlo_y))) {
541 if (!flag(ii,jj,k).isCovered())
544 int a_ind = (jj-(j-1)) + 3*(ii-(i-im));
549 Real x_off =
static_cast<Real
>(ii-i) + 0.5;
550 Real y_off =
static_cast<Real
>(jj-j);
553 if (ii < domlo_x && (vfrac(ii+1,jj,k) != 1.0 || vfrac(ii+2,jj,k) != 1.0) ) {
556 if (ii > domhi_x && (vfrac(ii-1,jj,k) != 1.0 || vfrac(ii-2,jj,k) != 1.0)) {
562 if (jj < domlo_y && (vfrac(ii,jj+1,k) != 1.0 || vfrac(ii,jj+2,k) != 1.0) ) {
565 if (jj > domhi_y && (vfrac(ii,jj-1,k) != 1.0 || vfrac(ii,jj-2,k) != 1.0) ) {
570 Amatrix(a_ind,0) = 1.0;
571 Amatrix(a_ind,1) = x_off + ccent(ii,jj,k,0);
572 Amatrix(a_ind,2) = y_off + ccent(ii,jj,k,1) - yloc_on_xface;
573 Amatrix(a_ind,3) = Amatrix(a_ind,1) * Amatrix(a_ind,1);
574 Amatrix(a_ind,4) = Amatrix(a_ind,1) * Amatrix(a_ind,2);
575 Amatrix(a_ind,5) = Amatrix(a_ind,2) * Amatrix(a_ind,2);
579 if (flag(ii,jj,k).isSingleValued() &&
580 domlo_x <= ii && ii <= domhi_x &&
581 domlo_y <= jj && jj <= domhi_y)
583 int a_ind_eb = a_ind + 9;
586 Amatrix(a_ind_eb,0) = 1.0;
587 Amatrix(a_ind_eb,1) = x_off + bcent(ii,jj,k,0);
588 Amatrix(a_ind_eb,2) = y_off + bcent(ii,jj,k,1) - yloc_on_xface;
589 Amatrix(a_ind_eb,3) = Amatrix(a_ind_eb,1) * Amatrix(a_ind_eb,1);
590 Amatrix(a_ind_eb,4) = Amatrix(a_ind_eb,1) * Amatrix(a_ind_eb,2);
591 Amatrix(a_ind_eb,5) = Amatrix(a_ind_eb,2) * Amatrix(a_ind_eb,2);
597 for (
int irow = 0; irow < 6; irow++)
601 for (
int ii = i-im; ii <= i+ip; ii++) {
602 for (
int jj = j-1; jj <= j+1; jj++)
608 if (!flag(ii,jj,k).isCovered())
610 int a_ind = (jj-(j-1)) + 3*(ii-(i-im));
613 Real phi_val = Amatrix(a_ind,0) * phi(ii,jj,k,n);
615 rhs(irow) += Amatrix(a_ind ,irow)*phi_val;
617 if (flag(ii,jj,k).isSingleValued() &&
618 is_eb_dirichlet && is_eb_inhomog && Amatrix(a_ind+9,irow) != 0.0) {
619 rhs(irow) += Amatrix(a_ind+9,irow)*phieb(ii,jj,k,n);
640 bool is_eb_dirichlet,
bool is_eb_inhomog,
641 const bool on_x_face,
const int domlo_x,
const int domhi_x,
642 const bool on_y_face,
const int domlo_y,
const int domhi_y)
651 for (
int irow = 0; irow < 18; irow++) {
652 for (
int icol = 0; icol < 6; icol++) {
653 Amatrix(irow,icol) = 0.0;
657 const int jm = (j > domhi_y) ? 2 : 1;
658 const int jp = 2 - jm;
661 for (
int jj = j-jm; jj <= j+jp; jj++) {
662 for (
int ii = i-1; ii <= i+1; ii++)
666 if (((on_x_face && ii < domlo_x) && (on_y_face && jj < domlo_y)) ||
667 ((on_x_face && ii > domhi_x) && (on_y_face && jj > domhi_y)) ||
668 ((on_x_face && ii < domlo_x) && (on_y_face && jj > domhi_y)) ||
669 ((on_x_face && ii > domhi_x) && (on_y_face && jj < domlo_y))) {
677 if (!flag(ii,jj,k).isCovered())
679 int a_ind = (ii-(i-1)) + 3*(jj-(j-jm));
682 Real x_off =
static_cast<Real
>(ii-i);
683 Real y_off =
static_cast<Real
>(jj-j) + 0.5;
686 if (ii < domlo_x && (vfrac(ii+1,jj,k) != 1.0 || vfrac(ii+2,jj,k) != 1.0) ) {
689 if (ii > domhi_x && (vfrac(ii-1,jj,k) != 1.0 || vfrac(ii-2,jj,k) != 1.0)) {
695 if (jj < domlo_y && (vfrac(ii,jj+1,k) != 1.0 || vfrac(ii,jj+2,k) != 1.0) ) {
698 if (jj > domhi_y && (vfrac(ii,jj-1,k) != 1.0 || vfrac(ii,jj-2,k) != 1.0) ) {
703 Amatrix(a_ind,0) = 1.0;
704 Amatrix(a_ind,1) = x_off + ccent(ii,jj,k,0) - xloc_on_yface;
705 Amatrix(a_ind,2) = y_off + ccent(ii,jj,k,1);
706 Amatrix(a_ind,3) = Amatrix(a_ind,1) * Amatrix(a_ind,1);
707 Amatrix(a_ind,4) = Amatrix(a_ind,1) * Amatrix(a_ind,2);
708 Amatrix(a_ind,5) = Amatrix(a_ind,2) * Amatrix(a_ind,2);
710 if (flag(ii,jj,k).isSingleValued() &&
711 domlo_x <= ii && ii <= domhi_x &&
712 domlo_y <= jj && jj <= domhi_y)
715 int a_ind_eb = a_ind + 9;
718 Amatrix(a_ind_eb,0) = 1.0;
719 Amatrix(a_ind_eb,1) = x_off + bcent(ii,jj,k,0) - xloc_on_yface;
720 Amatrix(a_ind_eb,2) = y_off + bcent(ii,jj,k,1);
721 Amatrix(a_ind_eb,3) = Amatrix(a_ind_eb,1) * Amatrix(a_ind_eb,1);
722 Amatrix(a_ind_eb,4) = Amatrix(a_ind_eb,1) * Amatrix(a_ind_eb,2);
723 Amatrix(a_ind_eb,5) = Amatrix(a_ind_eb,2) * Amatrix(a_ind_eb,2);
731 for (
int irow = 0; irow < 6; irow++)
735 for (
int jj = j-jm; jj <= j+jp; jj++) {
736 for (
int ii = i-1; ii <= i+1; ii++) {
740 if (!flag(ii,jj,k).isCovered())
742 int a_ind = (ii-(i-1)) + 3*(jj-(j-jm));
745 Real phi_val = Amatrix(a_ind,0)*phi(ii,jj,k,n);
747 rhs(irow) += Amatrix(a_ind,irow)*phi_val;
749 if (flag(ii,jj,k).isSingleValued() &&
750 is_eb_dirichlet && is_eb_inhomog && Amatrix(a_ind+9,irow) != Real(0.0)){
751 rhs(irow) += Amatrix(a_ind+9,irow)*phieb(ii,jj,k,n);
772 Real& nrmx, Real& nrmy,
774 const bool on_x_face,
const int domlo_x,
const int domhi_x,
775 const bool on_y_face,
const int domlo_y,
const int domhi_y)
782 for (
int irow = 0; irow < 18; irow++) {
783 for (
int icol = 0; icol < 6; icol++) {
784 Amatrix(irow,icol) = 0.0;
789 for (
int ii = i-1; ii <= i+1; ii++) {
790 for (
int jj = j-1; jj <= j+1; jj++)
796 if (((on_x_face && ii < domlo_x) && (on_y_face && jj < domlo_y)) ||
797 ((on_x_face && ii > domhi_x) && (on_y_face && jj > domhi_y))){
802 if (!flag(ii,jj,k).isCovered())
804 int a_ind = (jj-(j-1)) + 3*(ii-(i-1));
806 Real x_off =
static_cast<Real
>(ii-i);
807 Real y_off =
static_cast<Real
>(jj-j);
810 if (ii < domlo_x && (vfrac(ii+1,jj,k) != 1.0 || vfrac(ii+2,jj,k) != 1.0) ) {
813 if (ii > domhi_x && (vfrac(ii-1,jj,k) != 1.0 || vfrac(ii-2,jj,k) != 1.0)) {
819 if (jj < domlo_y && (vfrac(ii,jj+1,k) != 1.0 || vfrac(ii,jj+2,k) != 1.0) ) {
822 if (jj > domhi_y && (vfrac(ii,jj-1,k) != 1.0 || vfrac(ii,jj-2,k) != 1.0) ) {
828 if (flag(i,j,k).isConnected((ii-i),(jj-j),0))
830 Amatrix(a_ind,0) = 1.0;
831 Amatrix(a_ind,1) = x_off + ccent(ii,jj,k,0) - bcent(i,j,k,0);
832 Amatrix(a_ind,2) = y_off + ccent(ii,jj,k,1) - bcent(i,j,k,1);
836 if (flag(i,j,k).isConnected((ii-i),(jj-j),0) && flag(ii,jj,k).isSingleValued() &&
837 domlo_x <= ii && ii <= domhi_x && domlo_y <= jj && jj <= domhi_y)
839 Amatrix(a_ind+9,0) = 1.0;
840 Amatrix(a_ind+9,1) = x_off + bcent(ii,jj,k,0) - bcent(i,j,k,0);
841 Amatrix(a_ind+9,2) = y_off + bcent(ii,jj,k,1) - bcent(i,j,k,1);
849 for (
int irow = 0; irow < 18; irow++)
851 Amatrix(irow,3) = Amatrix(irow,1) * Amatrix(irow,1);
852 Amatrix(irow,4) = Amatrix(irow,1) * Amatrix(irow,2);
853 Amatrix(irow,5) = Amatrix(irow,2) * Amatrix(irow,2);
857 for (
int irow = 0; irow < 6; irow++)
861 for (
int ii = i-1; ii <= i+1; ii++) {
862 for (
int jj = j-1; jj <= j+1; jj++) {
866 if (!flag(ii,jj,k).isCovered())
868 int a_ind = (jj-(j-1)) + 3*(ii-(i-1));
872 Real phi_val = Amatrix(a_ind,0)*phi(ii,jj,k,n);
874 rhs(irow) += Amatrix(a_ind,irow) * phi_val;
876 if (flag(ii,jj,k).isSingleValued() && is_eb_inhomog && Amatrix(a_ind+9,irow) != 0.0) {
877 rhs(irow) += Amatrix(a_ind+9,irow)*phieb(ii,jj,k,n);
886 Real dphidn = rhs(1)*nrmx + rhs(2)*nrmy;
#define AMREX_ASSERT(EX)
Definition: AMReX_BLassert.H:38
#define AMREX_ALWAYS_ASSERT(EX)
Definition: AMReX_BLassert.H:50
#define AMREX_FORCE_INLINE
Definition: AMReX_Extension.H:119
#define AMREX_GPU_HOST_DEVICE
Definition: AMReX_GpuQualifiers.H:20
Definition: AMReX_Amr.cpp:49
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_y_of_phi_on_centroids(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Real &xloc_on_yface, bool is_eb_dirichlet, bool is_eb_inhomog)
Definition: AMReX_EB_LeastSquares_2D_K.H:331
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void cholsol_for_eb(Array2D< Real, 0, 17, 0, 5 > &Amatrix, Array1D< Real, 0, 5 > &b)
Definition: AMReX_EB_LeastSquares_2D_K.H:155
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_eb_of_phi_on_centroids(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Real &nrmx, Real &nrmy, bool is_eb_inhomog)
Definition: AMReX_EB_LeastSquares_2D_K.H:412
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void decomp_chol_np6(Array2D< Real, 0, 5, 0, 5 > &aa)
Definition: AMReX_EB_LeastSquares_2D_K.H:10
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_x_of_phi_on_centroids_extdir(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Array4< Real const > const &vfrac, Real &yloc_on_xface, bool is_eb_dirichlet, bool is_eb_inhomog, const bool on_x_face, const int domlo_x, const int domhi_x, const bool on_y_face, const int domlo_y, const int domhi_y)
Definition: AMReX_EB_LeastSquares_2D_K.H:495
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_eb_of_phi_on_centroids_extdir(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Array4< Real const > const &vfrac, Real &nrmx, Real &nrmy, bool is_eb_inhomog, const bool on_x_face, const int domlo_x, const int domhi_x, const bool on_y_face, const int domlo_y, const int domhi_y)
Definition: AMReX_EB_LeastSquares_2D_K.H:765
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_y_of_phi_on_centroids_extdir(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Array4< Real const > const &vfrac, Real &xloc_on_yface, bool is_eb_dirichlet, bool is_eb_inhomog, const bool on_x_face, const int domlo_x, const int domhi_x, const bool on_y_face, const int domlo_y, const int domhi_y)
Definition: AMReX_EB_LeastSquares_2D_K.H:632
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE void cholsol_np6(Array2D< Real, 0, 11, 0, 5 > &Amatrix, Array1D< Real, 0, 5 > &b)
Definition: AMReX_EB_LeastSquares_2D_K.H:64
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE Real grad_x_of_phi_on_centroids(int i, int j, int k, int n, Array4< Real const > const &phi, Array4< Real const > const &phieb, Array4< EBCellFlag const > const &flag, Array4< Real const > const &ccent, Array4< Real const > const &bcent, Real &yloc_on_xface, bool is_eb_dirichlet, bool is_eb_inhomog)
Definition: AMReX_EB_LeastSquares_2D_K.H:246
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE GpuComplex< T > sqrt(const GpuComplex< T > &a_z) noexcept
Return the square root of a complex number.
Definition: AMReX_GpuComplex.H:373
Definition: AMReX_Array.H:161
Definition: AMReX_Array.H:282
Definition: AMReX_Array4.H:61
AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE bool contains(int i, int j, int k) const noexcept
Definition: AMReX_Array4.H:251