Fix Field3D::setBoundaryTo parallel boundary conditions #2962

bendudson · 2024-08-17T19:14:31Z

Didn't handle parallel boundaries correctly. Now handles input arguments with and without parallelslices. If the input doesn't have parallel slices then they should be cleared in the target field, and to/fromFieldAligned used to set Y boundaries in field-aligned coordinates.

See discussion here of issue caused by previous implementation: boutproject/hermes-3@97a21fd

Didn't handle parallel boundaries correctly. Now handles input arguments with and without parallelslices. If the input doesn't have parallel slices then they should be cleared in the target field, and to/fromFieldAligned used to set Y boundaries in field-aligned coordinates.

Must clear P parallel slices or Grad_par will use them. Alternatively, setBoundaryTo should clear the parallel slices if they are not set in the argument.

github-actions

clang-tidy made some suggestions

src/field/field3d.cxx

github-actions · 2024-08-17T19:16:52Z

src/field/field3d.cxx

+    ASSERT1(hasParallelSlices());
+
+    for (auto& bndry : fieldmesh->getBoundariesPar()) {
+      const Field3D& f3d_next = f3d.ynext(bndry->dir);


warning: no member named 'x' in 'BoundaryRegionPar' [clang-diagnostic-error]

const int x = bndry->x; ^

github-actions · 2024-08-17T19:16:52Z

src/field/field3d.cxx

+
+    for (auto& bndry : fieldmesh->getBoundariesPar()) {
+      const Field3D& f3d_next = f3d.ynext(bndry->dir);
+      Field3D& this_next = ynext(bndry->dir);


warning: no member named 'y' in 'BoundaryRegionPar' [clang-diagnostic-error]

const int y = bndry->y; ^

github-actions · 2024-08-17T19:16:52Z

src/field/field3d.cxx

+    for (auto& bndry : fieldmesh->getBoundariesPar()) {
+      const Field3D& f3d_next = f3d.ynext(bndry->dir);
+      Field3D& this_next = ynext(bndry->dir);
+


warning: no member named 'z' in 'BoundaryRegionPar' [clang-diagnostic-error]

const int z = bndry->z; ^

github-actions · 2024-08-17T19:16:52Z

src/field/field3d.cxx

+    clearParallelSlices();
+
+    // Shift into field-aligned coordinates to apply Y boundary conditions.
+


warning: loop variable name 'r' is too short, expected at least 2 characters [readability-identifier-length]

for (RangeIterator r = fieldmesh->iterateBndryLowerY(); !r.isDone(); r++) { ^

src/field/field3d.cxx

github-actions · 2024-08-17T19:16:52Z

src/field/field3d.cxx

+    for (RangeIterator r = fieldmesh->iterateBndryLowerY(); !r.isDone(); r++) {
+      for (int jz = 0; jz < fieldmesh->LocalNz; jz++) {
+        BoutReal val = 0.5
+                       * (f3d_fa(r.ind, fieldmesh->ystart, jz)


warning: loop variable name 'r' is too short, expected at least 2 characters [readability-identifier-length]

for (RangeIterator r = fieldmesh->iterateBndryUpperY(); !r.isDone(); r++) { ^

github-actions · 2024-08-17T19:16:53Z

src/field/field3d.cxx

+        BoutReal val = 0.5
+                       * (f3d_fa(r.ind, fieldmesh->ystart, jz)
+                          + f3d_fa(r.ind, fieldmesh->ystart - 1, jz));
+        (*this)(r.ind, fieldmesh->ystart - 1, jz) =


warning: variable 'val' of type 'BoutReal' (aka 'double') can be declared 'const' [misc-const-correctness]

Suggested change

(*this)(r.ind, fieldmesh->ystart - 1, jz) =

BoutReal const val = 0.5 * (f3d_fa(r.ind, fieldmesh->yend, jz) + f3d_fa(r.ind, fieldmesh->yend + 1, jz));

dschwoerer

I thought the yboundary iterator had changed, but if it compiles it seems that has not yet been merged.

On a more general question, what is the purpose of this function? I do understand what it does, but why would one want to use it? Why not simply copy the boundaries? In what cases is this useful?

src/field/field3d.cxx

dschwoerer · 2024-08-19T09:56:03Z

src/field/field3d.cxx

+    for (auto& bndry : fieldmesh->getBoundariesPar()) {
+      const Field3D& f3d_next = f3d.ynext(bndry->dir);
+      Field3D& this_next = ynext(bndry->dir);
+
+      for (bndry->first(); !bndry->isDone(); bndry->next()) {
+        const int x = bndry->x;
+        const int y = bndry->y;
+        const int z = bndry->z;
+        BoutReal val = 0.5 * (f3d(x, y, z) + f3d_next(x, y, z));
+        this_next(x, y, z) = 2. * val - (*this)(x, y, z);
+      }
+    }


I think this should be roughly the new syntax:

Suggested change

for (auto& bndry : fieldmesh->getBoundariesPar()) {

const Field3D& f3d_next = f3d.ynext(bndry->dir);

Field3D& this_next = ynext(bndry->dir);

for (bndry->first(); !bndry->isDone(); bndry->next()) {

const int x = bndry->x;

const int y = bndry->y;

const int z = bndry->z;

BoutReal val = 0.5 * (f3d(x, y, z) + f3d_next(x, y, z));

this_next(x, y, z) = 2. * val - (*this)(x, y, z);

}

}

for (auto& region : fieldmesh->getBoundariesPar()) {

for (auto& pnt: region) {

pnt.dirichlet_o1(*this, pnt.interpolate_sheath_o1(f3d));

}

}

Co-authored-by: David Bold <dschwoerer@users.noreply.github.com>

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

Co-authored-by: David Bold <dschwoerer@users.noreply.github.com>

github-actions

clang-tidy made some suggestions

github-actions · 2024-08-19T15:40:22Z

src/field/field3d.cxx

+      for (int jy = fieldmesh->ystart; jy <= fieldmesh->yend; ++jy) {
+        for (int jz = 0; jz < fieldmesh->LocalNz; jz++) {
+          BoutReal const val =
+              0.5 * (f3d(fieldmesh->xstart, jy, jz) + f3d(fieldmesh->xstart - 1, jy, jz));


warning: expression result unused [clang-diagnostic-unused-value]

0.5 * (f3d(fieldmesh->xstart, jy, jz) + f3d(fieldmesh->xstart - 1, jy, jz)); ^

github-actions · 2024-08-19T15:40:22Z

src/field/field3d.cxx

+        }
+      }
+    }
+    if (fieldmesh->lastX()) {


warning: use of undeclared identifier 'jy' [clang-diagnostic-error]

BoutReal const val = 0.5 * (f3d(fieldmesh->xend, jy, jz) + f3d(fieldmesh->xend + 1, jy, jz)); ^

github-actions · 2024-08-19T15:40:23Z

src/field/field3d.cxx

+        }
+      }
+    }
+    if (fieldmesh->lastX()) {


warning: use of undeclared identifier 'jz' [clang-diagnostic-error]

BoutReal const val = 0.5 * (f3d(fieldmesh->xend, jy, jz) + f3d(fieldmesh->xend + 1, jy, jz)); ^

github-actions · 2024-08-19T15:40:23Z

src/field/field3d.cxx

+        }
+      }
+    }
+    if (fieldmesh->lastX()) {


warning: use of undeclared identifier 'jy' [clang-diagnostic-error]

BoutReal const val = 0.5 * (f3d(fieldmesh->xend, jy, jz) + f3d(fieldmesh->xend + 1, jy, jz)); ^

github-actions · 2024-08-19T15:40:23Z

src/field/field3d.cxx

+        }
+      }
+    }
+    if (fieldmesh->lastX()) {


warning: use of undeclared identifier 'jz' [clang-diagnostic-error]

BoutReal const val = 0.5 * (f3d(fieldmesh->xend, jy, jz) + f3d(fieldmesh->xend + 1, jy, jz)); ^

github-actions · 2024-08-19T15:40:23Z

src/field/field3d.cxx

+      BoutReal const val =
+          0.5 * (f3d(fieldmesh->xend, jy, jz) + f3d(fieldmesh->xend + 1, jy, jz));
+      for (int jz = 0; jz < fieldmesh->LocalNz; jz++) {
        BoutReal val =


warning: use of undeclared identifier 'jy' [clang-diagnostic-error]

(*this)(fieldmesh->xend + 1, jy, jz) = ^

github-actions · 2024-08-19T15:40:23Z

src/field/field3d.cxx

-        // Set to this value
-        (*this)(reg->x, reg->y, z) =
-            2. * val - (*this)(reg->x - reg->bx, reg->y - reg->by, z);
+            0.5 * (f3d(fieldmesh->xend, jy, jz) + f3d(fieldmesh->xend + 1, jy, jz));


warning: use of undeclared identifier 'jy' [clang-diagnostic-error]

2. * val - (*this)(fieldmesh->xend, jy, jz); ^

github-actions · 2024-08-19T15:40:23Z

src/field/field3d.cxx

  }
 }

+// Y boundaries


warning: expected unqualified-id [clang-diagnostic-error]

if (f3d.hasParallelSlices()) { ^

github-actions · 2024-08-19T15:40:23Z

src/field/field3d.cxx

+      BoutReal val = 0.5 * (f3d(x, y, z) + f3d_next(x, y + bndr->dir, z));
+      this_next(x, y + bndry->dir, z) = 2. * val - (*this)(x, y, z);
+    }
+  }


warning: expected unqualified-id [clang-diagnostic-error]

} else { ^

github-actions · 2024-08-19T15:40:23Z

src/field/field3d.cxx

+      (*this)(r.ind, fieldmesh->yend + 1, jz) =
+          2. * val - (*this)(r.ind, fieldmesh->yend, jz);
+    }
+  }


warning: extraneous closing brace ('}') [clang-diagnostic-error]

} ^

bendudson · 2024-08-19T15:43:29Z

Thanks @dschwoerer ! What is being done is:

A field (P) comes from the solver
It is modified, removing negative values and handling low density regions
The modified field has a boundary condition applied
The original field now should be modified to have the same boundary condition as the modified field.

One option is to just copy the boundary cells, but then the boundary condition would be slightly different. I'm not sure which would be better behaved numerically.

dschwoerer · 2024-08-21T12:06:53Z

I see, so the assumption is that both fields are very similar.
One thing that could happen with this method, if you want to ensure positive also of the boundary cells, and you modify them in this operation, you lose that property, and there is no way around this. Also note that this only copies Dirichlet boundary conditions, that should probably be documented.

I have pushed a fix, either as an additional commit (https://github.com/boutproject/BOUT-dev/tree/fix-setboundaryto-fixup) or squashed to remove the broken commits from history, for easier bisecting in the future (https://github.com/boutproject/BOUT-dev/tree/fix-setboundaryto-rebase)

dschwoerer · 2024-08-21T12:29:51Z

One test is failing as it triggers the assertion about parallel slices:

26/50 Test #28: test-invertable-operator ............***Failed    4.89 sec

dschwoerer · 2024-11-27T07:49:29Z

This has been replaced by #3030

bendudson added the bugfix label Aug 17, 2024

Apply clang-format changes

1aed8c5

bendudson referenced this pull request in boutproject/hermes-3 Aug 17, 2024

evolve_pressure: Important fix to P boundary

97a21fd

Must clear P parallel slices or Grad_par will use them. Alternatively, setBoundaryTo should clear the parallel slices if they are not set in the argument.

github-actions bot reviewed Aug 17, 2024

View reviewed changes

bendudson added the work in progress Not ready for merging label Aug 17, 2024

dschwoerer requested changes Aug 19, 2024

View reviewed changes

src/field/field3d.cxx Outdated Show resolved Hide resolved

src/field/field3d.cxx Outdated Show resolved Hide resolved

dschwoerer reviewed Aug 19, 2024

View reviewed changes

bendudson and others added 6 commits August 19, 2024 08:36

Update src/field/field3d.cxx

3e5e8dd

Co-authored-by: David Bold <dschwoerer@users.noreply.github.com>

Update src/field/field3d.cxx

c76446d

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

Update src/field/field3d.cxx

51e21b5

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

Update src/field/field3d.cxx

e5f5541

Co-authored-by: github-actions[bot] <41898282+github-actions[bot]@users.noreply.github.com>

Update src/field/field3d.cxx

8726ba9

Co-authored-by: David Bold <dschwoerer@users.noreply.github.com>

Apply clang-format changes

6d22c7d

github-actions bot reviewed Aug 19, 2024

View reviewed changes

dschwoerer mentioned this pull request Nov 26, 2024

Fix Field3D::setBoundaryTo for FCI methods #3030

Merged

dschwoerer closed this Nov 27, 2024

		clearParallelSlices();

		// Shift into field-aligned coordinates to apply Y boundary conditions.

	(*this)(r.ind, fieldmesh->ystart - 1, jz) =
	BoutReal const val = 0.5 * (f3d_fa(r.ind, fieldmesh->yend, jz) + f3d_fa(r.ind, fieldmesh->yend + 1, jz));

Fix Field3D::setBoundaryTo parallel boundary conditions #2962

Fix Field3D::setBoundaryTo parallel boundary conditions #2962

Uh oh!

Conversation

bendudson commented Aug 17, 2024

Uh oh!

github-actions bot left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Uh oh!

github-actions bot Aug 17, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 17, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 17, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 17, 2024

Choose a reason for hiding this comment

Uh oh!

Uh oh!

github-actions bot Aug 17, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 17, 2024

Choose a reason for hiding this comment

Uh oh!

dschwoerer left a comment

Choose a reason for hiding this comment

Uh oh!

Uh oh!

Uh oh!

dschwoerer Aug 19, 2024 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

github-actions bot left a comment

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

github-actions bot Aug 19, 2024

Choose a reason for hiding this comment

Uh oh!

bendudson commented Aug 19, 2024

Uh oh!

dschwoerer commented Aug 21, 2024

Uh oh!

dschwoerer commented Aug 21, 2024

Uh oh!

dschwoerer commented Nov 27, 2024

Uh oh!

dschwoerer Aug 19, 2024 •

edited

Loading