test_boundaries.c File Reference

C unit tests for boundary factories, inlet ownership, and face helpers. More...

#include "test_support.h"
#include "Boundaries.h"
#include "wallfunction.h"

Include dependency graph for test_boundaries.c:

Go to the source code of this file.

Functions
static PetscErrorCode	AppendBCParam (BC_Param **head, const char *key, const char *value)
	Appends one key/value pair to a linked list of boundary-condition parameters.
static PetscErrorCode	DestroyBoundaryHandler (BoundaryCondition **bc_ptr)
	Destroys one boundary-condition handler allocated by a boundary test.
static PetscInt	InteriorSampleIndex (PetscInt nodes)
	Chooses a stable interior node index for face-sample assertions on tiny test grids.
static PetscReal	GetFaceNormalComponent (Cmpnts value, BCFace face)
	Extracts the velocity component aligned with the supplied boundary-face normal.
static PetscReal	GetFaceOrientationSign (BCFace face)
	Returns the sign convention used for face-normal flux expectations.
static const char *	GetInletParamKey (BCFace face)
	Maps an inlet face to the matching configuration key used by the handler parser.
static PetscReal	GetFaceInteriorPointCount (const UserCtx *user, BCFace face)
	Computes the number of face-interior sample points for a given boundary face.
static PetscErrorCode	GetRepresentativeFaceSlots (const UserCtx *user, BCFace face, PetscInt *ucont_k, PetscInt *ucont_j, PetscInt *ucont_i, PetscInt *ubcs_k, PetscInt *ubcs_j, PetscInt *ubcs_i)
	Selects representative Ucont and Ubcs slots for face-matrix assertions.
static PetscErrorCode	GetParabolicSampleSlots (const UserCtx *user, BCFace face, PetscInt *center_k, PetscInt *center_j, PetscInt *center_i, PetscInt *off_k, PetscInt *off_j, PetscInt *off_i, PetscInt *wall_k, PetscInt *wall_j, PetscInt *wall_i, PetscInt *ubcs_center_k, PetscInt *ubcs_center_j, PetscInt *ubcs_center_i)
	Selects center, off-center, and wall-adjacent slots for parabolic-face checks.
static void	ResetBoundaryFaceConfig (BoundaryFaceConfig *cfg)
	Resets one boundary-face configuration entry to a neutral test-local baseline.
static PetscErrorCode	TestCanRankServiceFaceMatchesInletWhenDefined (void)
	Tests that face-service detection matches a defined inlet face.
static PetscErrorCode	TestCanRankServiceInletFaceRequiresDefinition (void)
	Tests that inlet-face service requires an inlet face to be defined.
static PetscErrorCode	TestBoundaryConditionFactoryAssignments (void)
	Tests the boundary-condition factory assignments for representative handlers.
static PetscErrorCode	ApplyInletAcrossFirstSteps (BoundaryCondition *bc, BCContext *ctx)
	Reapplies a static inlet through the first three time-step contexts.
static PetscErrorCode	TestBoundaryConditionFactoryImplementedHandlerMatrix (void)
	Tests the implemented-handler matrix exposed by the factory.
static PetscErrorCode	WritePicSliceForTests (const char *path, PetscInt n1, PetscInt n2, PetscReal base)
	Writes a tiny canonical PICSLICE profile for boundary handler tests.
static PetscErrorCode	TestInletProfileFromFileHandlerBehavior (void)
	Tests prescribed-flow inlet handler loading and applying one Z-face profile.
static PetscErrorCode	TestBoundaryConditionFactoryRejectsUnsupportedHandler (void)
	Tests that unsupported handler types are rejected by the factory.
static PetscErrorCode	TestGetDeterministicFaceGridLocationFaceMatrix (void)
	Tests deterministic inlet-face grid-location helpers across all faces.
static PetscErrorCode	TestGetRandomCellAndLogicalCoordsOnInletFaceMatrix (void)
	Tests random inlet-face cell selection across all supported faces.
static PetscErrorCode	TestInletConstantVelocityHandlerBehavior (void)
	Tests constant inlet handler initialization and face writes on a Z inlet.
static PetscErrorCode	TestInletCommonApplicationCurvilinearMetric (void)
	Verifies curvilinear projection and Cartesian dummy reflection for the shared inlet path.
static PetscErrorCode	TestWallNoSlipHandlerFaceMatrix (void)
	Tests wall no-slip application across the full face matrix.
static PetscErrorCode	TestInletConstantVelocityHandlerFaceMatrix (void)
	Tests constant-inlet initialization, flux accounting, and face writes across all faces.
static PetscErrorCode	TestInletParabolicProfileHandlerBehavior (void)
	Tests parabolic inlet handler shape on a tiny Z-face.
static PetscErrorCode	TestInletParabolicProfileUniformAlongPeriodicAxis (void)
	Tests that the parabolic inlet is uniform along a periodic cross-stream axis.
static PetscErrorCode	TestInletParabolicProfileHandlerFaceMatrix (void)
	Tests parabolic-inlet centerline, wall, and flux behavior across all faces.
static PetscErrorCode	TestOutletConservationHandlerBehavior (void)
	Tests outlet conservation handler correction and post-step flux accounting.
static PetscErrorCode	TestOutletConservationHandlerFaceMatrix (void)
	Tests outlet-conservation correction and flux accounting across all outlet faces.
static PetscErrorCode	TestFinalizePostProjectionCellFieldsNonPeriodic (void)
	Tests non-periodic cell-field finalization without modifying Ucont.
static PetscErrorCode	TestEnforceRHSBoundaryConditionsMatchesRowClassification (void)
	EnforceRHSBoundaryConditions() must zero a row exactly when it is not an unknown.
static PetscErrorCode	SeedWallFrictionVelocityStorage (UserCtx *user)
	Allocates the wall-model storage the way production allocates it.
static PetscErrorCode	SeedWallFunctionFixture (UserCtx *user, PetscReal tangential_speed)
	Seeds the fields ApplyWallFunction() reads on a single wall face.
static PetscErrorCode	TestApplyWallFunctionMatchesTheLogLawItDispatches (void)
	Tests that the wall-function integration reproduces the model it dispatches to.
static PetscErrorCode	TestWallFunctionModelSelectorChoosesTheLaw (void)
	Tests that the model selector actually changes which law is applied.
static PetscErrorCode	TestApplyWallFunctionSkipsNonWallFaces (void)
	Tests that a face which is not a wall is left alone by the wall function.
static PetscErrorCode	TestWallModelFrictionVelocityRespondsToTheFlow (void)
	Tests that every wall law's friction velocity responds to the flow it is given.
static PetscErrorCode	TestWallFrictionVelocityStorageMatchesTheCatalog (void)
	Tests that the friction-velocity storage matches what the field catalog declares.
static PetscErrorCode	TestWallModelDiagnosticsRecordTheCorrectedCells (void)
	Tests that the wall pass records the statistics its diagnostics report.
static PetscErrorCode	TestDummyCellsFollowEachFieldsBoundaryRule (void)
	Checks each cell-centred field's dummy-cell rule, faces then edges.
int	main (int argc, char **argv)
	Runs the unit-boundaries PETSc test binary.

Detailed Description

C unit tests for boundary factories, inlet ownership, and face helpers.

Definition in file test_boundaries.c.

Function Documentation

AppendBCParam()

static PetscErrorCode AppendBCParam ( BC_Param **  head, const char *  key, const char *  value  )

static

Appends one key/value pair to a linked list of boundary-condition parameters.

Definition at line 14 of file test_boundaries.c.

{
    BC_Param *node = NULL;
    BC_Param **cursor = head;
 
    PetscFunctionBeginUser;
    PetscCheck(head != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL, "BC param list head cannot be NULL.");
    PetscCall(PetscCalloc1(1, &node));
    PetscCall(PetscStrallocpy(key, &node->key));
    PetscCall(PetscStrallocpy(value, &node->value));
    while (*cursor) {
        cursor = &(*cursor)->next;
    }
    *cursor = node;
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

DestroyBoundaryHandler()

static PetscErrorCode DestroyBoundaryHandler ( BoundaryCondition **  bc_ptr )

static

Destroys one boundary-condition handler allocated by a boundary test.

Definition at line 34 of file test_boundaries.c.

{
    BoundaryCondition *bc = NULL;
    PetscFunctionBeginUser;
    if (!bc_ptr || !*bc_ptr) PetscFunctionReturn(0);
 
    bc = *bc_ptr;
    if (bc->Destroy) {
        PetscCall(bc->Destroy(bc));
    }
    PetscCall(PetscFree(bc));
    *bc_ptr = NULL;
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

InteriorSampleIndex()

static PetscInt InteriorSampleIndex ( PetscInt  nodes )

static

Chooses a stable interior node index for face-sample assertions on tiny test grids.

Definition at line 52 of file test_boundaries.c.

{
    PetscInt index = 3;
 
    if (index >= nodes - 1) {
        index = nodes - 2;
    }
    if (index < 1) {
        index = 1;
    }
    return index;
}

Here is the caller graph for this function:

GetFaceNormalComponent()

static PetscReal GetFaceNormalComponent ( Cmpnts  value, BCFace  face  )

static

Extracts the velocity component aligned with the supplied boundary-face normal.

Definition at line 68 of file test_boundaries.c.

{
    switch (face) {
        case BC_FACE_NEG_X:
        case BC_FACE_POS_X:
            return value.x;
        case BC_FACE_NEG_Y:
        case BC_FACE_POS_Y:
            return value.y;
        case BC_FACE_NEG_Z:
        case BC_FACE_POS_Z:
            return value.z;
    }
    return 0.0;
}

Here is the caller graph for this function:

GetFaceOrientationSign()

static PetscReal GetFaceOrientationSign ( BCFace  face )

static

Returns the sign convention used for face-normal flux expectations.

Definition at line 87 of file test_boundaries.c.

{
    switch (face) {
        case BC_FACE_NEG_X:
        case BC_FACE_NEG_Y:
        case BC_FACE_NEG_Z:
            return 1.0;
        case BC_FACE_POS_X:
        case BC_FACE_POS_Y:
        case BC_FACE_POS_Z:
            return -1.0;
    }
    return 0.0;
}

Here is the caller graph for this function:

GetInletParamKey()

static const char * GetInletParamKey ( BCFace  face )

static

Maps an inlet face to the matching configuration key used by the handler parser.

Definition at line 105 of file test_boundaries.c.

{
    switch (face) {
        case BC_FACE_NEG_X:
        case BC_FACE_POS_X:
            return "vx";
        case BC_FACE_NEG_Y:
        case BC_FACE_POS_Y:
            return "vy";
        case BC_FACE_NEG_Z:
        case BC_FACE_POS_Z:
            return "vz";
    }
    return "";
}

Here is the caller graph for this function:

GetFaceInteriorPointCount()

static PetscReal GetFaceInteriorPointCount ( const UserCtx *  user, BCFace  face  )

static

Computes the number of face-interior sample points for a given boundary face.

Definition at line 124 of file test_boundaries.c.

{
    switch (face) {
        case BC_FACE_NEG_X:
        case BC_FACE_POS_X:
            return (PetscReal)(user->info.my - 2) * (PetscReal)(user->info.mz - 2);
        case BC_FACE_NEG_Y:
        case BC_FACE_POS_Y:
            return (PetscReal)(user->info.mx - 2) * (PetscReal)(user->info.mz - 2);
        case BC_FACE_NEG_Z:
        case BC_FACE_POS_Z:
            return (PetscReal)(user->info.mx - 2) * (PetscReal)(user->info.my - 2);
    }
    return 0.0;
}

Here is the caller graph for this function:

GetRepresentativeFaceSlots()

static PetscErrorCode GetRepresentativeFaceSlots ( const UserCtx *  user, BCFace  face, PetscInt *  ucont_k, PetscInt *  ucont_j, PetscInt *  ucont_i, PetscInt *  ubcs_k, PetscInt *  ubcs_j, PetscInt *  ubcs_i  )

static

Selects representative Ucont and Ubcs slots for face-matrix assertions.

Definition at line 143 of file test_boundaries.c.

{
    const PetscInt sample_i = InteriorSampleIndex(user->info.mx);
    const PetscInt sample_j = InteriorSampleIndex(user->info.my);
    const PetscInt sample_k = InteriorSampleIndex(user->info.mz);
 
    PetscFunctionBeginUser;
    switch (face) {
        case BC_FACE_NEG_X:
            *ucont_k = sample_k; *ucont_j = sample_j; *ucont_i = 0;
            *ubcs_k = sample_k; *ubcs_j = sample_j; *ubcs_i = 0;
            break;
        case BC_FACE_POS_X:
            *ucont_k = sample_k; *ucont_j = sample_j; *ucont_i = user->info.mx - 2;
            *ubcs_k = sample_k; *ubcs_j = sample_j; *ubcs_i = user->info.mx - 1;
            break;
        case BC_FACE_NEG_Y:
            *ucont_k = sample_k; *ucont_j = 0; *ucont_i = sample_i;
            *ubcs_k = sample_k; *ubcs_j = 0; *ubcs_i = sample_i;
            break;
        case BC_FACE_POS_Y:
            *ucont_k = sample_k; *ucont_j = user->info.my - 2; *ucont_i = sample_i;
            *ubcs_k = sample_k; *ubcs_j = user->info.my - 1; *ubcs_i = sample_i;
            break;
        case BC_FACE_NEG_Z:
            *ucont_k = 0; *ucont_j = sample_j; *ucont_i = sample_i;
            *ubcs_k = 0; *ubcs_j = sample_j; *ubcs_i = sample_i;
            break;
        case BC_FACE_POS_Z:
            *ucont_k = user->info.mz - 2; *ucont_j = sample_j; *ucont_i = sample_i;
            *ubcs_k = user->info.mz - 1; *ubcs_j = sample_j; *ubcs_i = sample_i;
            break;
    }
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

GetParabolicSampleSlots()

static PetscErrorCode GetParabolicSampleSlots ( const UserCtx *  user, BCFace  face, PetscInt *  center_k, PetscInt *  center_j, PetscInt *  center_i, PetscInt *  off_k, PetscInt *  off_j, PetscInt *  off_i, PetscInt *  wall_k, PetscInt *  wall_j, PetscInt *  wall_i, PetscInt *  ubcs_center_k, PetscInt *  ubcs_center_j, PetscInt *  ubcs_center_i  )

static

Selects center, off-center, and wall-adjacent slots for parabolic-face checks.

Definition at line 189 of file test_boundaries.c.

{
    PetscFunctionBeginUser;
    switch (face) {
        case BC_FACE_NEG_X:
            *center_k = 3; *center_j = 3; *center_i = 0;
            *off_k = 2; *off_j = 2; *off_i = 0;
            *wall_k = 1; *wall_j = 1; *wall_i = 0;
            *ubcs_center_k = 3; *ubcs_center_j = 3; *ubcs_center_i = 0;
            break;
        case BC_FACE_POS_X:
            *center_k = 3; *center_j = 3; *center_i = user->info.mx - 2;
            *off_k = 2; *off_j = 2; *off_i = user->info.mx - 2;
            *wall_k = 1; *wall_j = 1; *wall_i = user->info.mx - 2;
            *ubcs_center_k = 3; *ubcs_center_j = 3; *ubcs_center_i = user->info.mx - 1;
            break;
        case BC_FACE_NEG_Y:
            *center_k = 3; *center_j = 0; *center_i = 3;
            *off_k = 2; *off_j = 0; *off_i = 2;
            *wall_k = 1; *wall_j = 0; *wall_i = 1;
            *ubcs_center_k = 3; *ubcs_center_j = 0; *ubcs_center_i = 3;
            break;
        case BC_FACE_POS_Y:
            *center_k = 3; *center_j = user->info.my - 2; *center_i = 3;
            *off_k = 2; *off_j = user->info.my - 2; *off_i = 2;
            *wall_k = 1; *wall_j = user->info.my - 2; *wall_i = 1;
            *ubcs_center_k = 3; *ubcs_center_j = user->info.my - 1; *ubcs_center_i = 3;
            break;
        case BC_FACE_NEG_Z:
            *center_k = 0; *center_j = 3; *center_i = 3;
            *off_k = 0; *off_j = 2; *off_i = 2;
            *wall_k = 0; *wall_j = 1; *wall_i = 1;
            *ubcs_center_k = 0; *ubcs_center_j = 3; *ubcs_center_i = 3;
            break;
        case BC_FACE_POS_Z:
            *center_k = user->info.mz - 2; *center_j = 3; *center_i = 3;
            *off_k = user->info.mz - 2; *off_j = 2; *off_i = 2;
            *wall_k = user->info.mz - 2; *wall_j = 1; *wall_i = 1;
            *ubcs_center_k = user->info.mz - 1; *ubcs_center_j = 3; *ubcs_center_i = 3;
            break;
    }
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

ResetBoundaryFaceConfig()

static void ResetBoundaryFaceConfig ( BoundaryFaceConfig *  cfg )

static

Resets one boundary-face configuration entry to a neutral test-local baseline.

Definition at line 249 of file test_boundaries.c.

{
    cfg->mathematical_type = INTERFACE;
    cfg->handler_type = BC_HANDLER_UNDEFINED;
    cfg->face_id = BC_FACE_NEG_X;
    if (cfg->params) {
        FreeBC_ParamList(cfg->params);
        cfg->params = NULL;
    }
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestCanRankServiceFaceMatchesInletWhenDefined()

static PetscErrorCode TestCanRankServiceFaceMatchesInletWhenDefined ( void  )

static

Tests that face-service detection matches a defined inlet face.

Definition at line 263 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    PetscBool generic_face = PETSC_FALSE;
    PetscBool inlet_face = PETSC_FALSE;
    PetscBool any_face_serviceable = PETSC_FALSE;
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
 
    for (BCFace face = BC_FACE_NEG_X; face <= BC_FACE_POS_Z; face++) {
        PetscCall(CanRankServiceFace(&user->info, user->IM, user->JM, user->KM, face, &generic_face));
        if (generic_face) {
            any_face_serviceable = PETSC_TRUE;
        }
 
        user->inletFaceDefined = PETSC_TRUE;
        user->identifiedInletBCFace = face;
        PetscCall(CanRankServiceInletFace(user, &user->info, user->IM, user->JM, user->KM, &inlet_face));
        PetscCall(PicurvAssertBool((PetscBool)(inlet_face == generic_face), "inlet face service check should match generic face check"));
    }
    PetscCall(PicurvAssertBool(any_face_serviceable, "at least one global face should be serviceable in a non-degenerate domain"));
 
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestCanRankServiceInletFaceRequiresDefinition()

static PetscErrorCode TestCanRankServiceInletFaceRequiresDefinition ( void  )

static

Tests that inlet-face service requires an inlet face to be defined.

Definition at line 294 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    PetscBool can_service = PETSC_FALSE;
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
 
    user->inletFaceDefined = PETSC_FALSE;
    for (BCFace face = BC_FACE_NEG_X; face <= BC_FACE_POS_Z; face++) {
        user->identifiedInletBCFace = face;
        PetscCall(CanRankServiceInletFace(user, &user->info, user->IM, user->JM, user->KM, &can_service));
        PetscCall(PicurvAssertBool((PetscBool)!can_service, "undefined inlet face should not be serviceable"));
    }
 
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestBoundaryConditionFactoryAssignments()

static PetscErrorCode TestBoundaryConditionFactoryAssignments ( void  )

static

Tests the boundary-condition factory assignments for representative handlers.

Definition at line 317 of file test_boundaries.c.

{
    BoundaryCondition *wall = NULL;
    BoundaryCondition *inlet = NULL;
    BoundaryCondition *parabolic = NULL;
    BoundaryCondition *profile = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(BoundaryCondition_Create(BC_HANDLER_WALL_NOSLIP, &wall));
    PetscCall(PicurvAssertBool((PetscBool)(wall != NULL), "factory should return wall handler"));
    PetscCall(PicurvAssertIntEqual(BC_PRIORITY_WALL, wall->priority, "wall handler priority"));
    PetscCall(PicurvAssertBool((PetscBool)(wall->Apply != NULL), "wall handler should expose Apply"));
    PetscCall(PicurvAssertBool((PetscBool)(wall->Destroy == NULL), "wall handler should not require destroy hook"));
    PetscCall(DestroyBoundaryHandler(&wall));
 
    PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_CONSTANT_VELOCITY, &inlet));
    PetscCall(PicurvAssertBool((PetscBool)(inlet != NULL), "factory should return inlet handler"));
    PetscCall(PicurvAssertIntEqual(BC_PRIORITY_INLET, inlet->priority, "inlet handler priority"));
    PetscCall(PicurvAssertBool((PetscBool)(inlet->Initialize != NULL), "inlet handler should expose Initialize"));
    PetscCall(PicurvAssertBool((PetscBool)(inlet->Apply != NULL), "inlet handler should expose Apply"));
    PetscCall(PicurvAssertBool((PetscBool)(inlet->Destroy != NULL), "inlet handler should expose its destroy hook"));
    PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_PARABOLIC, &parabolic));
    PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_PROFILE_FROM_FILE, &profile));
    PetscCall(PicurvAssertBool((PetscBool)(inlet->Apply == parabolic->Apply &&
                                           inlet->Apply == profile->Apply),
                               "all profile-producing inlet modes should share one application hook"));
    PetscCall(DestroyBoundaryHandler(&inlet));
    PetscCall(DestroyBoundaryHandler(&parabolic));
    PetscCall(DestroyBoundaryHandler(&profile));
 
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

ApplyInletAcrossFirstSteps()

static PetscErrorCode ApplyInletAcrossFirstSteps ( BoundaryCondition *  bc, BCContext *  ctx  )

static

Reapplies a static inlet through the first three time-step contexts.

Definition at line 353 of file test_boundaries.c.

{
    PetscFunctionBeginUser;
    for (PetscInt step = 0; step < 3; step++) {
        PetscReal local_inflow = 0.0, local_outflow = 0.0;
        ctx->user->simCtx->step = step;
        ctx->user->simCtx->ti = 0.1 * (PetscReal)step;
        PetscCall(VecSet(ctx->user->Ucont, -91.0));
        PetscCall(VecSet(ctx->user->Bcs.Ubcs, -92.0));
        PetscCall(bc->PreStep(bc, ctx, &local_inflow, &local_outflow));
        PetscCall(bc->Apply(bc, ctx));
    }
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

TestBoundaryConditionFactoryImplementedHandlerMatrix()

static PetscErrorCode TestBoundaryConditionFactoryImplementedHandlerMatrix ( void  )

static

Tests the implemented-handler matrix exposed by the factory.

Definition at line 371 of file test_boundaries.c.

{
    struct HandlerExpectation {
        BCHandlerType handler;
        BCPriorityType priority;
        PetscBool expect_apply;
        PetscBool expect_initialize;
    };
    const struct HandlerExpectation expectations[] = {
        {BC_HANDLER_WALL_NOSLIP, BC_PRIORITY_WALL, PETSC_TRUE, PETSC_FALSE},
        {BC_HANDLER_INLET_CONSTANT_VELOCITY, BC_PRIORITY_INLET, PETSC_TRUE, PETSC_TRUE},
        {BC_HANDLER_INLET_PARABOLIC, BC_PRIORITY_INLET, PETSC_TRUE, PETSC_TRUE},
        {BC_HANDLER_INLET_PROFILE_FROM_FILE, BC_PRIORITY_INLET, PETSC_TRUE, PETSC_TRUE},
        {BC_HANDLER_OUTLET_CONSERVATION, BC_PRIORITY_OUTLET, PETSC_TRUE, PETSC_FALSE},
    };
 
    PetscFunctionBeginUser;
    for (size_t i = 0; i < sizeof(expectations) / sizeof(expectations[0]); ++i) {
        BoundaryCondition *bc = NULL;
        PetscCall(BoundaryCondition_Create(expectations[i].handler, &bc));
        PetscCall(PicurvAssertBool((PetscBool)(bc != NULL), "factory should allocate a handler object"));
        PetscCall(PicurvAssertIntEqual(expectations[i].priority, bc->priority, "handler priority should match expectation"));
        PetscCall(PicurvAssertBool((PetscBool)((bc->Apply != NULL) == expectations[i].expect_apply), "Apply hook expectation mismatch"));
        PetscCall(PicurvAssertBool((PetscBool)((bc->Initialize != NULL) == expectations[i].expect_initialize), "Initialize hook expectation mismatch"));
        PetscCall(DestroyBoundaryHandler(&bc));
    }
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

WritePicSliceForTests()

static PetscErrorCode WritePicSliceForTests ( const char *  path, PetscInt  n1, PetscInt  n2, PetscReal  base  )

static

Writes a tiny canonical PICSLICE profile for boundary handler tests.

Definition at line 403 of file test_boundaries.c.

{
    FILE *fd = NULL;
 
    PetscFunctionBeginUser;
    fd = fopen(path, "w");
    PetscCheck(fd != NULL, PETSC_COMM_SELF, PETSC_ERR_FILE_OPEN, "Could not open %s for writing.", path);
    PetscCheck(fprintf(fd, "PICSLICE\n1\n%d %d\n", (int)n1, (int)n2) >= 0,
               PETSC_COMM_SELF, PETSC_ERR_FILE_WRITE, "Could not write PICSLICE header.");
    for (PetscInt a = 0; a < n1; a++) {
        for (PetscInt b = 0; b < n2; b++) {
            PetscReal value = base + (PetscReal)(10 * a + b);
            PetscCheck(fprintf(fd, "%.16e\n", (double)value) >= 0,
                       PETSC_COMM_SELF, PETSC_ERR_FILE_WRITE, "Could not write PICSLICE value.");
        }
    }
    fclose(fd);
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

TestInletProfileFromFileHandlerBehavior()

static PetscErrorCode TestInletProfileFromFileHandlerBehavior ( void  )

static

Tests prescribed-flow inlet handler loading and applying one Z-face profile.

Definition at line 426 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BoundaryCondition *bc = NULL;
    BCContext ctx;
    Cmpnts ***ucont = NULL;
    Cmpnts ***ubcs = NULL;
    char profile_path[PETSC_MAX_PATH_LEN];
    PetscReal local_inflow = 0.0;
    PetscReal local_outflow = 0.0;
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 5, 5, 5));
    PetscCall(PicurvPopulateIdentityMetrics(user));
 
    PetscCall(PetscSNPrintf(profile_path, sizeof(profile_path), "/tmp/picurv_unit_profile_%d.picslice", (int)getpid()));
    PetscCall(WritePicSliceForTests(profile_path, user->JM - 1, user->IM - 1, 2.0));
 
    user->boundary_faces[BC_FACE_NEG_Z].face_id = BC_FACE_NEG_Z;
    user->boundary_faces[BC_FACE_NEG_Z].mathematical_type = INLET;
    user->boundary_faces[BC_FACE_NEG_Z].handler_type = BC_HANDLER_INLET_PROFILE_FROM_FILE;
    PetscCall(AppendBCParam(&user->boundary_faces[BC_FACE_NEG_Z].params, "source_file", profile_path));
 
    PetscCall(VecSet(user->Ucont, 0.0));
    PetscCall(VecSet(user->Bcs.Ubcs, 0.0));
 
    ctx.user = user;
    ctx.face_id = BC_FACE_NEG_Z;
    PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_PROFILE_FROM_FILE, &bc));
    PetscCall(bc->PreStep(bc, &ctx, &local_inflow, &local_outflow));
    PetscCall(PicurvAssertRealNear(0.0, local_inflow, 1.0e-12, "profile inlet PreStep should leave inflow unchanged"));
    PetscCall(PicurvAssertRealNear(0.0, local_outflow, 1.0e-12, "profile inlet PreStep should leave outflow unchanged"));
    PetscCall(bc->Initialize(bc, &ctx));
    PetscCall(ApplyInletAcrossFirstSteps(bc, &ctx));
    PetscCall(bc->PostStep(bc, &ctx, &local_inflow, &local_outflow));
 
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
    PetscCall(PicurvAssertRealNear(2.0, ucont[0][1][1].z, 1.0e-12, "profile inlet should apply first scalar speed"));
    PetscCall(PicurvAssertRealNear(13.0, ucont[0][2][2].z, 1.0e-12, "profile inlet should preserve PICSLICE ordering"));
    PetscCall(PicurvAssertRealNear(24.0, ubcs[0][3][3].z, 1.0e-12, "profile inlet should write boundary velocity"));
    PetscCall(PicurvAssertRealNear(35.0, ucont[0][4][4].z, 1.0e-12, "profile inlet should consume the last face-cell slot"));
    PetscCall(PicurvAssertRealNear(35.0, ubcs[0][4][4].z, 1.0e-12, "profile inlet should write the last face-cell boundary value"));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
    PetscCall(PicurvAssertBool((PetscBool)(local_inflow > 0.0), "profile inlet PostStep should accumulate positive negative-face flux"));
    PetscCall(PicurvAssertRealNear(0.0, local_outflow, 1.0e-12, "profile inlet should not add to outflow"));
 
    PetscCall(DestroyBoundaryHandler(&bc));
    FreeBC_ParamList(user->boundary_faces[BC_FACE_NEG_Z].params);
    user->boundary_faces[BC_FACE_NEG_Z].params = NULL;
    remove(profile_path);
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestBoundaryConditionFactoryRejectsUnsupportedHandler()

static PetscErrorCode TestBoundaryConditionFactoryRejectsUnsupportedHandler ( void  )

static

Tests that unsupported handler types are rejected by the factory.

Definition at line 487 of file test_boundaries.c.

{
    BoundaryCondition *bc = NULL;
    PetscErrorCode ierr_create;
 
    PetscFunctionBeginUser;
    PetscCall(PetscPushErrorHandler(PetscIgnoreErrorHandler, NULL));
    ierr_create = BoundaryCondition_Create(BC_HANDLER_OUTLET_PRESSURE, &bc);
    PetscCall(PetscPopErrorHandler());
 
    PetscCall(PicurvAssertBool((PetscBool)(ierr_create != 0), "unsupported handler should return a non-zero error code"));
    if (bc) {
        PetscCall(PetscFree(bc));
    }
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestGetDeterministicFaceGridLocationFaceMatrix()

static PetscErrorCode TestGetDeterministicFaceGridLocationFaceMatrix ( void  )

static

Tests deterministic inlet-face grid-location helpers across all faces.

Definition at line 507 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 8, 8, 8));
    simCtx->np = 128;
 
    for (BCFace face = BC_FACE_NEG_X; face <= BC_FACE_POS_Z; ++face) {
        PetscInt ci = -1, cj = -1, ck = -1;
        PetscReal xi = -1.0, eta = -1.0, zta = -1.0;
        PetscBool placed = PETSC_FALSE;
 
        user->identifiedInletBCFace = face;
        PetscCall(GetDeterministicFaceGridLocation(
            user, &user->info,
            user->info.xs, user->info.ys, user->info.zs,
            user->info.mx, user->info.my, user->info.mz,
            0,
            &ci, &cj, &ck, &xi, &eta, &zta, &placed));
 
        PetscCall(PicurvAssertBool(placed, "single-rank deterministic face placement should succeed"));
        PetscCall(PicurvAssertBool((PetscBool)(ci >= user->info.xs && ci < user->info.xs + user->info.xm), "ci must map to owned node window"));
        PetscCall(PicurvAssertBool((PetscBool)(cj >= user->info.ys && cj < user->info.ys + user->info.ym), "cj must map to owned node window"));
        PetscCall(PicurvAssertBool((PetscBool)(ck >= user->info.zs && ck < user->info.zs + user->info.zm), "ck must map to owned node window"));
        PetscCall(PicurvAssertBool((PetscBool)(xi >= 0.0 && xi < 1.0), "xi should be in [0,1)"));
        PetscCall(PicurvAssertBool((PetscBool)(eta >= 0.0 && eta < 1.0), "eta should be in [0,1)"));
        PetscCall(PicurvAssertBool((PetscBool)(zta >= 0.0 && zta < 1.0), "zta should be in [0,1)"));
 
        if (face == BC_FACE_NEG_X || face == BC_FACE_POS_X) {
            PetscCall(PicurvAssertRealNear(0.5, xi, 1.0e-10, "deterministic x-face placement should sit halfway into boundary-adjacent cell"));
        } else if (face == BC_FACE_NEG_Y || face == BC_FACE_POS_Y) {
            PetscCall(PicurvAssertRealNear(0.5, eta, 1.0e-10, "deterministic y-face placement should sit halfway into boundary-adjacent cell"));
        } else {
            PetscCall(PicurvAssertRealNear(0.5, zta, 1.0e-10, "deterministic z-face placement should sit halfway into boundary-adjacent cell"));
        }
    }
 
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestGetRandomCellAndLogicalCoordsOnInletFaceMatrix()

static PetscErrorCode TestGetRandomCellAndLogicalCoordsOnInletFaceMatrix ( void  )

static

Tests random inlet-face cell selection across all supported faces.

Definition at line 553 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    PetscRandom rand_i = NULL, rand_j = NULL, rand_k = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 8, 8, 8));
 
    PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rand_i));
    PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rand_j));
    PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rand_k));
    PetscCall(PetscRandomSetInterval(rand_i, 0.0, 1.0));
    PetscCall(PetscRandomSetInterval(rand_j, 0.0, 1.0));
    PetscCall(PetscRandomSetInterval(rand_k, 0.0, 1.0));
    PetscCall(PetscRandomSetFromOptions(rand_i));
    PetscCall(PetscRandomSetFromOptions(rand_j));
    PetscCall(PetscRandomSetFromOptions(rand_k));
 
    for (BCFace face = BC_FACE_NEG_X; face <= BC_FACE_POS_Z; ++face) {
        PetscInt ci = -1, cj = -1, ck = -1;
        PetscReal xi = -1.0, eta = -1.0, zta = -1.0;
        const PetscReal boundary_eps = 5.0e-4;
 
        user->identifiedInletBCFace = face;
        PetscCall(GetRandomCellAndLogicalCoordsOnInletFace(
            user, &user->info,
            user->info.xs, user->info.ys, user->info.zs,
            user->info.mx, user->info.my, user->info.mz,
            &rand_i, &rand_j, &rand_k,
            &ci, &cj, &ck,
            &xi, &eta, &zta));
 
        PetscCall(PicurvAssertBool((PetscBool)(ci >= user->info.xs && ci < user->info.xs + user->info.xm), "ci must map to owned node window"));
        PetscCall(PicurvAssertBool((PetscBool)(cj >= user->info.ys && cj < user->info.ys + user->info.ym), "cj must map to owned node window"));
        PetscCall(PicurvAssertBool((PetscBool)(ck >= user->info.zs && ck < user->info.zs + user->info.zm), "ck must map to owned node window"));
        PetscCall(PicurvAssertBool((PetscBool)(xi >= 0.0 && xi <= 1.0), "xi should be in [0,1]"));
        PetscCall(PicurvAssertBool((PetscBool)(eta >= 0.0 && eta <= 1.0), "eta should be in [0,1]"));
        PetscCall(PicurvAssertBool((PetscBool)(zta >= 0.0 && zta <= 1.0), "zta should be in [0,1]"));
 
        switch (face) {
            case BC_FACE_NEG_X:
                PetscCall(PicurvAssertBool((PetscBool)(xi <= boundary_eps), "NEG_X should pin xi near 0"));
                break;
            case BC_FACE_POS_X:
                PetscCall(PicurvAssertBool((PetscBool)(xi >= 1.0 - boundary_eps), "POS_X should pin xi near 1"));
                break;
            case BC_FACE_NEG_Y:
                PetscCall(PicurvAssertBool((PetscBool)(eta <= boundary_eps), "NEG_Y should pin eta near 0"));
                break;
            case BC_FACE_POS_Y:
                PetscCall(PicurvAssertBool((PetscBool)(eta >= 1.0 - boundary_eps), "POS_Y should pin eta near 1"));
                break;
            case BC_FACE_NEG_Z:
                PetscCall(PicurvAssertBool((PetscBool)(zta <= boundary_eps), "NEG_Z should pin zta near 0"));
                break;
            case BC_FACE_POS_Z:
                PetscCall(PicurvAssertBool((PetscBool)(zta >= 1.0 - boundary_eps), "POS_Z should pin zta near 1"));
                break;
        }
    }
 
    PetscCall(PetscRandomDestroy(&rand_i));
    PetscCall(PetscRandomDestroy(&rand_j));
    PetscCall(PetscRandomDestroy(&rand_k));
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestInletConstantVelocityHandlerBehavior()

static PetscErrorCode TestInletConstantVelocityHandlerBehavior ( void  )

static

Tests constant inlet handler initialization and face writes on a Z inlet.

Definition at line 625 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BoundaryCondition *bc = NULL;
    BCContext ctx;
    Cmpnts ***ucont = NULL;
    Cmpnts ***ubcs = NULL;
    PetscReal local_inflow = 0.0;
    PetscReal local_outflow = 0.0;
    PetscReal expected_flux = 0.0;
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(VecSet(user->Ucont, 0.0));
    PetscCall(VecSet(user->Bcs.Ubcs, 0.0));
    PetscCall(PicurvPopulateIdentityMetrics(user));
 
    user->boundary_faces[BC_FACE_NEG_Z].face_id = BC_FACE_NEG_Z;
    user->boundary_faces[BC_FACE_NEG_Z].mathematical_type = INLET;
    user->boundary_faces[BC_FACE_NEG_Z].handler_type = BC_HANDLER_INLET_CONSTANT_VELOCITY;
    PetscCall(AppendBCParam(&user->boundary_faces[BC_FACE_NEG_Z].params, "vz", "2.5"));
 
    ctx.user = user;
    ctx.face_id = BC_FACE_NEG_Z;
    PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_CONSTANT_VELOCITY, &bc));
    PetscCall(bc->Initialize(bc, &ctx));
    PetscCall(ApplyInletAcrossFirstSteps(bc, &ctx));
    PetscCall(bc->PostStep(bc, &ctx, &local_inflow, &local_outflow));
 
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
    PetscCall(PicurvAssertRealNear(2.5, ucont[0][3][3].z, 1.0e-12, "constant inlet should set Ucont normal component"));
    PetscCall(PicurvAssertRealNear(2.5, ubcs[0][3][3].z, 1.0e-12, "constant inlet should set Ubcs normal component"));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
 
    expected_flux = (PetscReal)(user->info.mx - 2) * (PetscReal)(user->info.my - 2) * 2.5;
    PetscCall(PicurvAssertRealNear(expected_flux, local_inflow, 1.0e-12, "constant inlet PostStep should sum the face flux"));
    PetscCall(PicurvAssertRealNear(0.0, local_outflow, 1.0e-12, "constant inlet should not contribute to outflow accumulation"));
 
    PetscCall(DestroyBoundaryHandler(&bc));
    FreeBC_ParamList(user->boundary_faces[BC_FACE_NEG_Z].params);
    user->boundary_faces[BC_FACE_NEG_Z].params = NULL;
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestInletCommonApplicationCurvilinearMetric()

static PetscErrorCode TestInletCommonApplicationCurvilinearMetric ( void  )

static

Verifies curvilinear projection and Cartesian dummy reflection for the shared inlet path.

Definition at line 677 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BoundaryCondition *bc = NULL;
    BCContext ctx;
    Cmpnts ***zet = NULL, ***ucont = NULL, ***ubcs = NULL, ***ucat = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(PicurvPopulateIdentityMetrics(user));
    PetscCall(DMDAVecGetArray(user->fda, user->lZet, &zet));
    zet[0][3][3].x = 0.0;
    zet[0][3][3].y = 3.0;
    zet[0][3][3].z = 4.0;
    PetscCall(DMDAVecRestoreArray(user->fda, user->lZet, &zet));
 
    user->boundary_faces[BC_FACE_NEG_Z].face_id = BC_FACE_NEG_Z;
    user->boundary_faces[BC_FACE_NEG_Z].mathematical_type = INLET;
    user->boundary_faces[BC_FACE_NEG_Z].handler_type = BC_HANDLER_INLET_CONSTANT_VELOCITY;
    PetscCall(AppendBCParam(&user->boundary_faces[BC_FACE_NEG_Z].params, "vz", "2.0"));
    ctx.user = user;
    ctx.face_id = BC_FACE_NEG_Z;
    PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_CONSTANT_VELOCITY, &bc));
    PetscCall(bc->Initialize(bc, &ctx));
    PetscCall(bc->Apply(bc, &ctx));
 
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
    PetscCall(PicurvAssertRealNear(10.0, ucont[0][3][3].z, 1.0e-12,
                                   "Ucont should be the Cartesian inlet velocity dotted with the face metric"));
    PetscCall(PicurvAssertRealNear(0.0, ubcs[0][3][3].x, 1.0e-12,
                                   "curvilinear normal velocity should preserve its Cartesian x component"));
    PetscCall(PicurvAssertRealNear(1.2, ubcs[0][3][3].y, 1.0e-12,
                                   "curvilinear normal velocity should include Cartesian y"));
    PetscCall(PicurvAssertRealNear(1.6, ubcs[0][3][3].z, 1.0e-12,
                                   "curvilinear normal velocity should include Cartesian z"));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
 
    PetscCall(VecSet(user->Ucat, 1.0));
    PetscCall(UpdateDummyCells(user, FIELD_ID_UCAT));
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucat, &ucat));
    PetscCall(PicurvAssertRealNear(-1.0, ucat[0][3][3].x, 1.0e-12,
                                   "dummy reflection should act on Cartesian x"));
    PetscCall(PicurvAssertRealNear(1.4, ucat[0][3][3].y, 1.0e-12,
                                   "dummy reflection should act on Cartesian y"));
    PetscCall(PicurvAssertRealNear(2.2, ucat[0][3][3].z, 1.0e-12,
                                   "dummy reflection should act on Cartesian z"));
    PetscCall(PicurvAssertRealNear(1.0, ucat[1][3][3].y, 1.0e-12,
                                   "dummy reflection should not overwrite the adjacent interior value"));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucat, &ucat));
 
    PetscCall(DestroyBoundaryHandler(&bc));
    ResetBoundaryFaceConfig(&user->boundary_faces[BC_FACE_NEG_Z]);
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestWallNoSlipHandlerFaceMatrix()

static PetscErrorCode TestWallNoSlipHandlerFaceMatrix ( void  )

static

Tests wall no-slip application across the full face matrix.

Definition at line 740 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BoundaryCondition *bc = NULL;
    BCContext ctx;
    Cmpnts ***ucont = NULL;
    Cmpnts ***ubcs = NULL;
    const BCFace faces[] = {
        BC_FACE_NEG_X, BC_FACE_POS_X,
        BC_FACE_NEG_Y, BC_FACE_POS_Y,
        BC_FACE_NEG_Z, BC_FACE_POS_Z
    };
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(BoundaryCondition_Create(BC_HANDLER_WALL_NOSLIP, &bc));
 
    for (size_t n = 0; n < sizeof(faces) / sizeof(faces[0]); ++n) {
        const BCFace face = faces[n];
        PetscInt ucont_k, ucont_j, ucont_i, ubcs_k, ubcs_j, ubcs_i;
 
        PetscCall(VecSet(user->Ucont, 5.0));
        PetscCall(VecSet(user->Bcs.Ubcs, 7.0));
        ctx.user = user;
        ctx.face_id = face;
        PetscCall(bc->Apply(bc, &ctx));
 
        PetscCall(GetRepresentativeFaceSlots(user, face, &ucont_k, &ucont_j, &ucont_i, &ubcs_k, &ubcs_j, &ubcs_i));
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
        PetscCall(PicurvAssertRealNear(0.0, GetFaceNormalComponent(ucont[ucont_k][ucont_j][ucont_i], face), 1.0e-12, "wall no-slip should zero the normal flux component"));
        PetscCall(PicurvAssertRealNear(0.0, ubcs[ubcs_k][ubcs_j][ubcs_i].x, 1.0e-12, "wall no-slip should zero Ubcs.x"));
        PetscCall(PicurvAssertRealNear(0.0, ubcs[ubcs_k][ubcs_j][ubcs_i].y, 1.0e-12, "wall no-slip should zero Ubcs.y"));
        PetscCall(PicurvAssertRealNear(0.0, ubcs[ubcs_k][ubcs_j][ubcs_i].z, 1.0e-12, "wall no-slip should zero Ubcs.z"));
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
    }
 
    PetscCall(DestroyBoundaryHandler(&bc));
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestInletConstantVelocityHandlerFaceMatrix()

static PetscErrorCode TestInletConstantVelocityHandlerFaceMatrix ( void  )

static

Tests constant-inlet initialization, flux accounting, and face writes across all faces.

Definition at line 788 of file test_boundaries.c.

{
    struct ConstantInletCase {
        BCFace face;
        PetscReal velocity;
        const char *value_text;
    };
    const struct ConstantInletCase cases[] = {
        {BC_FACE_NEG_X, 1.5, "1.5"},
        {BC_FACE_POS_X, 2.0, "2.0"},
        {BC_FACE_NEG_Y, 2.5, "2.5"},
        {BC_FACE_POS_Y, 3.0, "3.0"},
        {BC_FACE_NEG_Z, 3.5, "3.5"},
        {BC_FACE_POS_Z, 4.0, "4.0"},
    };
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BCContext ctx;
    Cmpnts ***ucont = NULL;
    Cmpnts ***ubcs = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(PicurvPopulateIdentityMetrics(user));
 
    for (size_t n = 0; n < sizeof(cases) / sizeof(cases[0]); ++n) {
        const struct ConstantInletCase *test_case = &cases[n];
        BoundaryCondition *bc = NULL;
        PetscReal local_inflow = 0.0;
        PetscReal local_outflow = 0.0;
        PetscInt ucont_k, ucont_j, ucont_i, ubcs_k, ubcs_j, ubcs_i;
        const PetscReal expected_value = GetFaceOrientationSign(test_case->face) * test_case->velocity;
        const PetscReal expected_flux = GetFaceInteriorPointCount(user, test_case->face) * expected_value;
 
        PetscCall(VecSet(user->Ucont, 0.0));
        PetscCall(VecSet(user->Bcs.Ubcs, 0.0));
        ResetBoundaryFaceConfig(&user->boundary_faces[test_case->face]);
        user->boundary_faces[test_case->face].face_id = test_case->face;
        user->boundary_faces[test_case->face].mathematical_type = INLET;
        user->boundary_faces[test_case->face].handler_type = BC_HANDLER_INLET_CONSTANT_VELOCITY;
        PetscCall(AppendBCParam(&user->boundary_faces[test_case->face].params,
                                GetInletParamKey(test_case->face),
                                test_case->value_text));
 
        ctx.user = user;
        ctx.face_id = test_case->face;
        PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_CONSTANT_VELOCITY, &bc));
        PetscCall(bc->PreStep(bc, &ctx, &local_inflow, &local_outflow));
        PetscCall(PicurvAssertRealNear(0.0, local_inflow, 1.0e-12, "constant inlet PreStep should leave inflow unchanged"));
        PetscCall(PicurvAssertRealNear(0.0, local_outflow, 1.0e-12, "constant inlet PreStep should leave outflow unchanged"));
        PetscCall(bc->Initialize(bc, &ctx));
        PetscCall(bc->Apply(bc, &ctx));
        PetscCall(bc->PostStep(bc, &ctx, &local_inflow, &local_outflow));
 
        PetscCall(GetRepresentativeFaceSlots(user, test_case->face, &ucont_k, &ucont_j, &ucont_i, &ubcs_k, &ubcs_j, &ubcs_i));
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
        PetscCall(PicurvAssertRealNear(expected_value, GetFaceNormalComponent(ucont[ucont_k][ucont_j][ucont_i], test_case->face), 1.0e-12, "constant inlet should set the face-normal Ucont component"));
        PetscCall(PicurvAssertRealNear(expected_value, GetFaceNormalComponent(ubcs[ubcs_k][ubcs_j][ubcs_i], test_case->face), 1.0e-12, "constant inlet should set the face-normal Ubcs component"));
        PetscCall(PicurvAssertRealNear((test_case->face == BC_FACE_NEG_X || test_case->face == BC_FACE_POS_X) ? expected_value : 0.0, ubcs[ubcs_k][ubcs_j][ubcs_i].x, 1.0e-12, "constant inlet should only set the expected Ubcs axis"));
        PetscCall(PicurvAssertRealNear((test_case->face == BC_FACE_NEG_Y || test_case->face == BC_FACE_POS_Y) ? expected_value : 0.0, ubcs[ubcs_k][ubcs_j][ubcs_i].y, 1.0e-12, "constant inlet should only set the expected Ubcs axis"));
        PetscCall(PicurvAssertRealNear((test_case->face == BC_FACE_NEG_Z || test_case->face == BC_FACE_POS_Z) ? expected_value : 0.0, ubcs[ubcs_k][ubcs_j][ubcs_i].z, 1.0e-12, "constant inlet should only set the expected Ubcs axis"));
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
        PetscCall(PicurvAssertRealNear(expected_flux, local_inflow, 1.0e-12, "constant inlet PostStep should integrate the face flux with orientation"));
        PetscCall(PicurvAssertRealNear(0.0, local_outflow, 1.0e-12, "constant inlet should not add to outflow"));
 
        PetscCall(DestroyBoundaryHandler(&bc));
        ResetBoundaryFaceConfig(&user->boundary_faces[test_case->face]);
    }
 
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestInletParabolicProfileHandlerBehavior()

static PetscErrorCode TestInletParabolicProfileHandlerBehavior ( void  )

static

Tests parabolic inlet handler shape on a tiny Z-face.

Definition at line 867 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BoundaryCondition *bc = NULL;
    BCContext ctx;
    Cmpnts ***ucont = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(PicurvPopulateIdentityMetrics(user));
 
    user->boundary_faces[BC_FACE_NEG_Z].face_id = BC_FACE_NEG_Z;
    user->boundary_faces[BC_FACE_NEG_Z].mathematical_type = INLET;
    user->boundary_faces[BC_FACE_NEG_Z].handler_type = BC_HANDLER_INLET_PARABOLIC;
    PetscCall(AppendBCParam(&user->boundary_faces[BC_FACE_NEG_Z].params, "v_max", "4.0"));
 
    ctx.user = user;
    ctx.face_id = BC_FACE_NEG_Z;
    PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_PARABOLIC, &bc));
    PetscCall(bc->Initialize(bc, &ctx));
    PetscCall(ApplyInletAcrossFirstSteps(bc, &ctx));
 
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(PicurvAssertRealNear(4.0, ucont[0][3][3].z, 1.0e-12, "parabolic inlet should peak at the face centerline"));
    /* Five cells per axis: the wall-adjacent centre sits at s = -0.8, so each factor is 0.36. */
    PetscCall(PicurvAssertRealNear(4.0 * 0.36 * 0.36, ucont[0][1][1].z, 1.0e-12,
                                   "parabolic inlet should vanish on the wall, not at the wall-adjacent centre"));
    PetscCall(PicurvAssertBool((PetscBool)(ucont[0][3][3].z > ucont[0][2][2].z), "parabolic inlet centerline should exceed off-center velocity"));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
 
    PetscCall(DestroyBoundaryHandler(&bc));
    FreeBC_ParamList(user->boundary_faces[BC_FACE_NEG_Z].params);
    user->boundary_faces[BC_FACE_NEG_Z].params = NULL;
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestInletParabolicProfileUniformAlongPeriodicAxis()

static PetscErrorCode TestInletParabolicProfileUniformAlongPeriodicAxis ( void  )

static

Tests that the parabolic inlet is uniform along a periodic cross-stream axis.

A channel periodic in i has walls only in j, so its inlet is the exact one-dimensional parabola in j.

Definition at line 911 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BoundaryCondition *bc = NULL;
    BCContext ctx;
    Cmpnts ***ucont = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(PicurvPopulateIdentityMetrics(user));
    simCtx->i_periodic = 1;
 
    user->boundary_faces[BC_FACE_NEG_Z].face_id = BC_FACE_NEG_Z;
    user->boundary_faces[BC_FACE_NEG_Z].mathematical_type = INLET;
    user->boundary_faces[BC_FACE_NEG_Z].handler_type = BC_HANDLER_INLET_PARABOLIC;
    PetscCall(AppendBCParam(&user->boundary_faces[BC_FACE_NEG_Z].params, "v_max", "4.0"));
    ctx.user = user;
    ctx.face_id = BC_FACE_NEG_Z;
    PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_PARABOLIC, &bc));
    PetscCall(bc->Initialize(bc, &ctx));
 
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(PicurvAssertRealNear(4.0, ucont[0][3][1].z, 1.0e-12, "periodic axis carries no parabola factor"));
    PetscCall(PicurvAssertRealNear(ucont[0][3][1].z, ucont[0][3][4].z, 1.0e-12, "profile is uniform along the periodic axis"));
    PetscCall(PicurvAssertRealNear(4.0 * 0.36, ucont[0][1][3].z, 1.0e-12, "wall-bounded axis keeps its parabola"));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
 
    PetscCall(DestroyBoundaryHandler(&bc));
    FreeBC_ParamList(user->boundary_faces[BC_FACE_NEG_Z].params);
    user->boundary_faces[BC_FACE_NEG_Z].params = NULL;
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestInletParabolicProfileHandlerFaceMatrix()

static PetscErrorCode TestInletParabolicProfileHandlerFaceMatrix ( void  )

static

Tests parabolic-inlet centerline, wall, and flux behavior across all faces.

Definition at line 950 of file test_boundaries.c.

{
    const BCFace faces[] = {
        BC_FACE_NEG_X, BC_FACE_POS_X,
        BC_FACE_NEG_Y, BC_FACE_POS_Y,
        BC_FACE_NEG_Z, BC_FACE_POS_Z
    };
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BCContext ctx;
    Cmpnts ***ucont = NULL;
    Cmpnts ***ubcs = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(PicurvPopulateIdentityMetrics(user));
 
    for (size_t n = 0; n < sizeof(faces) / sizeof(faces[0]); ++n) {
        const BCFace face = faces[n];
        BoundaryCondition *bc = NULL;
        PetscReal local_inflow = 0.0;
        PetscReal local_outflow = 0.0;
        PetscInt center_k, center_j, center_i;
        PetscInt off_k, off_j, off_i;
        PetscInt wall_k, wall_j, wall_i;
        PetscInt ubcs_center_k, ubcs_center_j, ubcs_center_i;
        const PetscReal expected_center = GetFaceOrientationSign(face) * 4.0;
 
        PetscCall(VecSet(user->Ucont, 0.0));
        PetscCall(VecSet(user->Bcs.Ubcs, 0.0));
        ResetBoundaryFaceConfig(&user->boundary_faces[face]);
        user->boundary_faces[face].face_id = face;
        user->boundary_faces[face].mathematical_type = INLET;
        user->boundary_faces[face].handler_type = BC_HANDLER_INLET_PARABOLIC;
        PetscCall(AppendBCParam(&user->boundary_faces[face].params, "v_max", "4.0"));
 
        ctx.user = user;
        ctx.face_id = face;
        PetscCall(BoundaryCondition_Create(BC_HANDLER_INLET_PARABOLIC, &bc));
        PetscCall(bc->PreStep(bc, &ctx, &local_inflow, &local_outflow));
        PetscCall(PicurvAssertRealNear(0.0, local_inflow, 1.0e-12, "parabolic inlet PreStep should leave inflow unchanged"));
        PetscCall(PicurvAssertRealNear(0.0, local_outflow, 1.0e-12, "parabolic inlet PreStep should leave outflow unchanged"));
        PetscCall(bc->Initialize(bc, &ctx));
        PetscCall(bc->Apply(bc, &ctx));
        PetscCall(bc->PostStep(bc, &ctx, &local_inflow, &local_outflow));
 
        PetscCall(GetParabolicSampleSlots(user, face,
                                          &center_k, &center_j, &center_i,
                                          &off_k, &off_j, &off_i,
                                          &wall_k, &wall_j, &wall_i,
                                          &ubcs_center_k, &ubcs_center_j, &ubcs_center_i));
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
        PetscCall(PicurvAssertRealNear(expected_center, GetFaceNormalComponent(ucont[center_k][center_j][center_i], face), 1.0e-12, "parabolic inlet should peak at the face centerline"));
        PetscCall(PicurvAssertRealNear(expected_center, GetFaceNormalComponent(ubcs[ubcs_center_k][ubcs_center_j][ubcs_center_i], face), 1.0e-12, "parabolic inlet should write the centerline boundary velocity"));
        PetscCall(PicurvAssertRealNear(expected_center * 0.36 * 0.36, GetFaceNormalComponent(ucont[wall_k][wall_j][wall_i], face), 1.0e-12,
                                       "parabolic inlet should take its wall-adjacent value, vanishing only on the wall"));
        PetscCall(PicurvAssertBool((PetscBool)(PetscAbsReal(GetFaceNormalComponent(ucont[center_k][center_j][center_i], face)) >
                                               PetscAbsReal(GetFaceNormalComponent(ucont[off_k][off_j][off_i], face))),
                                   "parabolic inlet centerline magnitude should exceed the off-center magnitude"));
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
        PetscCall(PicurvAssertBool((PetscBool)(PetscAbsReal(local_inflow) > 0.0), "parabolic inlet PostStep should accumulate a non-zero flux"));
        PetscCall(PicurvAssertBool((PetscBool)((local_inflow > 0.0) == (GetFaceOrientationSign(face) > 0.0)),
                                   "parabolic inlet PostStep should preserve the face orientation"));
        PetscCall(PicurvAssertRealNear(0.0, local_outflow, 1.0e-12, "parabolic inlet should not add to outflow"));
 
        PetscCall(DestroyBoundaryHandler(&bc));
        ResetBoundaryFaceConfig(&user->boundary_faces[face]);
    }
 
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestOutletConservationHandlerBehavior()

static PetscErrorCode TestOutletConservationHandlerBehavior ( void  )

static

Tests outlet conservation handler correction and post-step flux accounting.

Definition at line 1029 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BoundaryCondition *bc = NULL;
    BCContext ctx;
    PetscReal global_inflow = 30.0;
    PetscReal global_farfield_in = 0.0;
    PetscReal global_farfield_out = 0.0;
    PetscReal global_outflow = 25.0;
    PetscReal local_inflow = 0.0;
    PetscReal local_outflow = 0.0;
    Cmpnts ***ucont = NULL;
    Cmpnts ***ucat = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(PicurvPopulateIdentityMetrics(user));
    PetscCall(VecSet(user->Ucat, 1.0));
    PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
    PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
    simCtx->AreaOutSum = (PetscReal)(user->info.mx - 2) * (PetscReal)(user->info.my - 2);
 
    user->boundary_faces[BC_FACE_POS_Z].face_id = BC_FACE_POS_Z;
    user->boundary_faces[BC_FACE_POS_Z].mathematical_type = OUTLET;
    user->boundary_faces[BC_FACE_POS_Z].handler_type = BC_HANDLER_OUTLET_CONSERVATION;
 
    ctx.user = user;
    ctx.face_id = BC_FACE_POS_Z;
    ctx.global_inflow_sum = &global_inflow;
    ctx.global_farfield_inflow_sum = &global_farfield_in;
    ctx.global_farfield_outflow_sum = &global_farfield_out;
    ctx.global_outflow_sum = &global_outflow;
 
    PetscCall(BoundaryCondition_Create(BC_HANDLER_OUTLET_CONSERVATION, &bc));
    PetscCall(bc->Apply(bc, &ctx));
    PetscCall(bc->PostStep(bc, &ctx, &local_inflow, &local_outflow));
 
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(DMDAVecGetArrayRead(user->fda, user->lUcat, &ucat));
    PetscCall(PicurvAssertRealNear(1.2, ucont[5][3][3].z, 1.0e-12, "outlet correction should add the expected flux trim"));
    PetscCall(PicurvAssertRealNear(1.0, ucat[5][3][3].z, 1.0e-12, "outlet handler should preserve the interior Ucat reference"));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->lUcat, &ucat));
 
    PetscCall(PicurvAssertRealNear(30.0, local_outflow, 1.0e-12, "outlet PostStep should report corrected flux"));
    PetscCall(PicurvAssertRealNear(0.0, local_inflow, 1.0e-12, "outlet handler should not accumulate inflow"));
 
    PetscCall(DestroyBoundaryHandler(&bc));
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestOutletConservationHandlerFaceMatrix()

static PetscErrorCode TestOutletConservationHandlerFaceMatrix ( void  )

static

Tests outlet-conservation correction and flux accounting across all outlet faces.

Definition at line 1086 of file test_boundaries.c.

{
    const BCFace faces[] = {
        BC_FACE_NEG_X, BC_FACE_POS_X,
        BC_FACE_NEG_Y, BC_FACE_POS_Y,
        BC_FACE_NEG_Z, BC_FACE_POS_Z
    };
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    BCContext ctx;
    Cmpnts ***ucont = NULL;
    Cmpnts ***ubcs = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PetscMemzero(&ctx, sizeof(ctx)));
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(PicurvPopulateIdentityMetrics(user));
    PetscCall(VecSet(user->Ucat, 1.0));
    PetscCall(DMGlobalToLocalBegin(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
    PetscCall(DMGlobalToLocalEnd(user->fda, user->Ucat, INSERT_VALUES, user->lUcat));
 
    for (size_t n = 0; n < sizeof(faces) / sizeof(faces[0]); ++n) {
        const BCFace face = faces[n];
        BoundaryCondition *bc = NULL;
        PetscReal global_inflow = 0.0;
        PetscReal global_farfield_in = 0.0;
        PetscReal global_farfield_out = 0.0;
        PetscReal global_outflow = 0.0;
        PetscReal local_inflow = 0.0;
        PetscReal local_outflow = 0.0;
        PetscReal area = GetFaceInteriorPointCount(user, face);
        PetscInt ucont_k, ucont_j, ucont_i, ubcs_k, ubcs_j, ubcs_i;
 
        PetscCall(VecSet(user->Ucont, 0.0));
        PetscCall(VecSet(user->Bcs.Ubcs, 0.0));
        ResetBoundaryFaceConfig(&user->boundary_faces[face]);
        user->boundary_faces[face].face_id = face;
        user->boundary_faces[face].mathematical_type = OUTLET;
        user->boundary_faces[face].handler_type = BC_HANDLER_OUTLET_CONSERVATION;
        user->simCtx->AreaOutSum = area;
 
        ctx.user = user;
        ctx.face_id = face;
        ctx.global_inflow_sum = &global_inflow;
        ctx.global_farfield_inflow_sum = &global_farfield_in;
        ctx.global_farfield_outflow_sum = &global_farfield_out;
        ctx.global_outflow_sum = &global_outflow;
 
        PetscCall(BoundaryCondition_Create(BC_HANDLER_OUTLET_CONSERVATION, &bc));
        PetscCall(bc->PreStep(bc, &ctx, &local_inflow, &local_outflow));
        PetscCall(PicurvAssertRealNear(area, local_outflow, 1.0e-12, "outlet PreStep should measure the uncorrected face flux"));
 
        global_inflow = area + 0.2 * area;
        global_outflow = local_outflow;
        local_outflow = 0.0;
        PetscCall(bc->Apply(bc, &ctx));
        PetscCall(bc->PostStep(bc, &ctx, &local_inflow, &local_outflow));
 
        PetscCall(GetRepresentativeFaceSlots(user, face, &ucont_k, &ucont_j, &ucont_i, &ubcs_k, &ubcs_j, &ubcs_i));
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucont, &ucont));
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
        PetscCall(PicurvAssertRealNear(1.2, GetFaceNormalComponent(ucont[ucont_k][ucont_j][ucont_i], face), 1.0e-12, "outlet correction should add the configured flux trim"));
        PetscCall(PicurvAssertRealNear(1.0, ubcs[ubcs_k][ubcs_j][ubcs_i].x, 1.0e-12, "outlet handler should copy Ucat into Ubcs.x"));
        PetscCall(PicurvAssertRealNear(1.0, ubcs[ubcs_k][ubcs_j][ubcs_i].y, 1.0e-12, "outlet handler should copy Ucat into Ubcs.y"));
        PetscCall(PicurvAssertRealNear(1.0, ubcs[ubcs_k][ubcs_j][ubcs_i].z, 1.0e-12, "outlet handler should copy Ucat into Ubcs.z"));
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucont, &ucont));
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Bcs.Ubcs, &ubcs));
        PetscCall(PicurvAssertRealNear(0.0, local_inflow, 1.0e-12, "outlet handler should not add to inflow"));
        PetscCall(PicurvAssertRealNear(1.2 * area, local_outflow, 1.0e-12, "outlet PostStep should report the corrected flux"));
 
        PetscCall(DestroyBoundaryHandler(&bc));
        ResetBoundaryFaceConfig(&user->boundary_faces[face]);
    }
 
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestFinalizePostProjectionCellFieldsNonPeriodic()

static PetscErrorCode TestFinalizePostProjectionCellFieldsNonPeriodic ( void  )

static

Tests non-periodic cell-field finalization without modifying Ucont.

Definition at line 1167 of file test_boundaries.c.

{
    SimCtx *simCtx = NULL;
    UserCtx *user = NULL;
    Cmpnts ***ucat = NULL;
    Cmpnts ***lucat = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 4, 4, 4));
    PetscCall(VecSet(user->Ucat, 4.0));
    PetscCall(VecSet(user->Bcs.Ubcs, 1.0));
    PetscCall(VecSet(user->P, 3.0));
    PetscCall(VecSet(user->Ucont, 7.0));
 
    PetscCall(FinalizePostProjectionCellFields(user));
 
    PetscCall(PicurvAssertVecConstant(user->Ucont, 7.0, 1.0e-12,
                                      "post-projection cell finalization must preserve non-periodic Ucont"));
    PetscCall(PicurvAssertVecConstant(user->P, 3.0, 1.0e-12,
                                      "constant pressure should remain unchanged during cell finalization"));
 
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucat, &ucat));
    PetscCall(DMDAVecGetArrayRead(user->fda, user->lUcat, &lucat));
    PetscCall(PicurvAssertRealNear(-2.0, ucat[2][2][0].x, 1.0e-12,
                                    "non-periodic Ucat dummy face should use Ubcs extrapolation"));
    PetscCall(PicurvAssertRealNear(-2.0, lucat[2][2][0].x, 1.0e-12,
                                    "non-periodic local Ucat should contain the finalized dummy face"));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucat, &ucat));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->lUcat, &lucat));
 
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestEnforceRHSBoundaryConditionsMatchesRowClassification()

static PetscErrorCode TestEnforceRHSBoundaryConditionsMatchesRowClassification ( void  )

static

EnforceRHSBoundaryConditions() must zero a row exactly when it is not an unknown.

This is the invariant that ties the residual path to ClassifyMomentumRow(). It is checked across a non-periodic, a singly-periodic and a fully-periodic configuration because the defect it guards against was visible only when an axis was periodic: the duplicate column at index 0 was left unzeroed, so ComputeTotalResidual()'s BDF term accumulated there on every call and the residual norm stopped describing the state. A sentinel fill makes an untouched row detectable rather than accidentally correct.

Definition at line 1211 of file test_boundaries.c.

{
    const struct {
        const char *label;
        PetscBool   x_periodic, y_periodic, z_periodic;
    } configs[] = {
        {"non-periodic",     PETSC_FALSE, PETSC_FALSE, PETSC_FALSE},
        {"xi-periodic",      PETSC_TRUE,  PETSC_FALSE, PETSC_FALSE},
        {"fully-periodic",   PETSC_TRUE,  PETSC_TRUE,  PETSC_TRUE },
    };
    const PetscScalar sentinel = 7.0;
 
    PetscFunctionBeginUser;
    for (size_t c = 0; c < sizeof(configs) / sizeof(configs[0]); ++c) {
        SimCtx  *simCtx = NULL;
        UserCtx *user = NULL;
        Cmpnts ***rhs = NULL;
        const PetscBool periodic[3] = {configs[c].x_periodic, configs[c].y_periodic, configs[c].z_periodic};
        const BCFace neg_face[3] = {BC_FACE_NEG_X, BC_FACE_NEG_Y, BC_FACE_NEG_Z};
        const BCFace pos_face[3] = {BC_FACE_POS_X, BC_FACE_POS_Y, BC_FACE_POS_Z};
 
        PetscCall(PicurvCreateMinimalContextsWithPeriodicity(&simCtx, &user, 6, 6, 6,
                                                             configs[c].x_periodic,
                                                             configs[c].y_periodic,
                                                             configs[c].z_periodic));
        for (PetscInt axis = 0; axis < 3; ++axis) {
            const BCType type = periodic[axis] ? PERIODIC : WALL;
            user->boundary_faces[neg_face[axis]].mathematical_type = type;
            user->boundary_faces[pos_face[axis]].mathematical_type = type;
        }
 
        PetscCall(VecDuplicate(user->Ucont, &user->Rhs));
        PetscCall(VecSet(user->Rhs, sentinel));
        PetscCall(EnforceRHSBoundaryConditions(user));
 
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Rhs, &rhs));
        for (PetscInt k = user->info.zs; k < user->info.zs + user->info.zm; ++k) {
            for (PetscInt j = user->info.ys; j < user->info.ys + user->info.ym; ++j) {
                for (PetscInt i = user->info.xs; i < user->info.xs + user->info.xm; ++i) {
                    const PetscScalar *row = &rhs[k][j][i].x;
                    for (PetscInt component = 0; component < 3; ++component) {
                        PetscInt ri, rj, rk;
                        const MomentumRowType type =
                            ClassifyMomentumRow(user, i, j, k, component, &ri, &rj, &rk);
                        const PetscScalar expected = (type == MOM_ROW_PHYSICAL) ? sentinel : 0.0;
                        char context[192];
 
                        PetscCall(PetscSNPrintf(context, sizeof(context),
                            "%s: row (i=%d,j=%d,k=%d,c=%d) classified %d must be %s",
                            configs[c].label, (int)i, (int)j, (int)k, (int)component, (int)type,
                            (type == MOM_ROW_PHYSICAL) ? "left untouched" : "zeroed"));
                        PetscCall(PicurvAssertRealNear((PetscReal)expected, (PetscReal)row[component],
                                                       1.0e-14, context));
                    }
                }
            }
        }
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Rhs, &rhs));
 
        PetscCall(VecDestroy(&user->Rhs));
        PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    }
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

SeedWallFrictionVelocityStorage()

static PetscErrorCode SeedWallFrictionVelocityStorage ( UserCtx *  user )

static

Allocates the wall-model storage the way production allocates it.

Production creates each global on da and ghosts it into the local view. A fixture that allocates differently - a different DM, only one of a pair, or only one of the fields the same block creates - tests a layout the solver never uses, which is how a DM/dof mismatch survived in ApplyWallFunction while these cases passed. The wall eddy viscosity is allocated in that same block and cleared at the top of every pass, so a fixture omitting it dereferences a null vector before reaching any assertion. Every wall case routes through here so the two cannot drift.

Definition at line 1287 of file test_boundaries.c.

{
    PetscFunctionBeginUser;
    PetscCall(DMCreateGlobalVector(user->da, &user->Friction_Velocity));
    PetscCall(VecSet(user->Friction_Velocity, 0.0));
    PetscCall(DMCreateLocalVector(user->da, &user->lFriction_Velocity));
    PetscCall(VecSet(user->lFriction_Velocity, 0.0));
    PetscCall(DMCreateGlobalVector(user->da, &user->Nu_Wall));
    PetscCall(VecSet(user->Nu_Wall, 0.0));
    PetscCall(DMCreateLocalVector(user->da, &user->lNu_Wall));
    PetscCall(VecSet(user->lNu_Wall, 0.0));
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

SeedWallFunctionFixture()

static PetscErrorCode SeedWallFunctionFixture ( UserCtx *  user, PetscReal  tangential_speed  )

static

Seeds the fields ApplyWallFunction() reads on a single wall face.

A uniform tangential velocity and unit metrics make the wall distances and the face normal analytically known, so the test can reproduce the model's own arithmetic without copying the integration's index choices.

Definition at line 1308 of file test_boundaries.c.

{
    Cmpnts ***ucat = NULL, ***csi = NULL;
    PetscReal ***aj = NULL;
    DMDALocalInfo info = user->info;
 
    PetscFunctionBeginUser;
    PetscCall(VecSet(user->Aj, 1.0));
    PetscCall(DMDAVecGetArray(user->fda, user->Ucat, &ucat));
    PetscCall(DMDAVecGetArray(user->fda, user->Csi, &csi));
    PetscCall(DMDAVecGetArray(user->da, user->Aj, &aj));
    for (PetscInt k = info.zs; k < info.zs + info.zm; ++k)
        for (PetscInt j = info.ys; j < info.ys + info.ym; ++j)
            for (PetscInt i = info.xs; i < info.xs + info.xm; ++i) {
                /* Flow along eta, so the whole velocity is tangential to a xi wall. */
                ucat[k][j][i] = (Cmpnts){.x = 0.0, .y = tangential_speed, .z = 0.0};
                csi[k][j][i] = (Cmpnts){.x = 1.0, .y = 0.0, .z = 0.0};
                aj[k][j][i] = 1.0;
            }
    PetscCall(DMDAVecRestoreArray(user->da, user->Aj, &aj));
    PetscCall(DMDAVecRestoreArray(user->fda, user->Csi, &csi));
    PetscCall(DMDAVecRestoreArray(user->fda, user->Ucat, &ucat));
    PetscCall(DMGlobalToLocalBegin(user->da, user->Aj, INSERT_VALUES, user->lAj));
    PetscCall(DMGlobalToLocalEnd(user->da, user->Aj, INSERT_VALUES, user->lAj));
    PetscCall(DMGlobalToLocalBegin(user->fda, user->Csi, INSERT_VALUES, user->lCsi));
    PetscCall(DMGlobalToLocalEnd(user->fda, user->Csi, INSERT_VALUES, user->lCsi));
    PetscFunctionReturn(PETSC_SUCCESS);
}

Here is the caller graph for this function:

TestApplyWallFunctionMatchesTheLogLawItDispatches()

static PetscErrorCode TestApplyWallFunctionMatchesTheLogLawItDispatches ( void  )

static

Tests that the wall-function integration reproduces the model it dispatches to.

The kernels are covered separately; what is unverified without this is the wiring around them, where the plausible failures are silent. Feeding the reference and boundary wall distances in the wrong order, taking the reference velocity from the wrong cell, or normalising the face normal incorrectly all still produce a plausible-looking corrected velocity. The check here is that the corrected first-interior velocity is exactly what the log-law model predicts at that cell's distance, for the friction velocity the integration itself stored.

Definition at line 1348 of file test_boundaries.c.

{
    SimCtx      *simCtx = NULL;
    UserCtx     *user = NULL;
    Cmpnts    ***ucat = NULL, ***ubcs = NULL;
    PetscReal ***utau = NULL;
    const PetscReal speed = 1.0;
    const PetscReal roughness = 1.0e-16;   /* smooth wall */
    PetscReal    kinematic_viscosity, distance_first, distance_second;
    PetscReal    stored_utau, corrected_speed, expected_speed, reference_check;
    PetscReal    untouched_before = 0.0, untouched_after = 0.0;
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    simCtx->ren = 1000.0;
    simCtx->wall_roughness_height = roughness;
    kinematic_viscosity = 1.0 / simCtx->ren;
 
    PetscCall(SeedWallFrictionVelocityStorage(user));
    user->boundary_faces[BC_FACE_NEG_X].mathematical_type = WALL;
    PetscCall(SeedWallFunctionFixture(user, speed));
 
    /* With unit Jacobians and a unit face area the model's distances are exact:
       the first cell centre sits half a cell from the wall, the second one and a
       half. Reproducing them here rather than reading them back is what makes a
       swapped-distance wiring error detectable. */
    distance_first = 0.5;
    distance_second = 2.0 * distance_first + 0.5;
 
    /* Disabled: the routine must not touch the field at all. */
    simCtx->wallfunction = 0;
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
    untouched_before = ucat[2][2][1].y;
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
    PetscCall(ApplyWallFunction(user));
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
    untouched_after = ucat[2][2][1].y;
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
    PetscCall(PicurvAssertRealNear(untouched_before, untouched_after, 0.0,
                                   "a disabled wall function must leave the velocity untouched"));
 
    /* Enabled: the first interior cell is corrected onto the modelled profile. */
    simCtx->wallfunction = 1;
    PetscCall(ApplyWallFunction(user));
 
    PetscCall(DMDAVecGetArrayRead(user->da, user->Friction_Velocity, (const PetscReal ***)&utau));
    stored_utau = utau[2][2][1];
    PetscCall(DMDAVecRestoreArrayRead(user->da, user->Friction_Velocity, (const PetscReal ***)&utau));
    PetscCall(PicurvAssertBool((PetscBool)(stored_utau > 0.0),
                               "the integration must store a positive friction velocity"));
 
    /* The stored friction velocity must be the one that reproduces the *reference*
       cell's speed at the *reference* distance. Swapping the two distances passes the
       next check but fails this one. */
    reference_check = u_hydset_roughness(kinematic_viscosity, distance_second,
                                         stored_utau, roughness);
    PetscCall(PicurvAssertRealNear(speed, reference_check, 1.0e-6,
                                   "the friction velocity must invert the log law at the "
                                   "reference cell's distance"));
 
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
    corrected_speed = ucat[2][2][1].y;
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
 
    expected_speed = u_hydset_roughness(kinematic_viscosity, distance_first,
                                        stored_utau, roughness);
    PetscCall(PicurvAssertRealNear(expected_speed, corrected_speed, 1.0e-9,
                                   "the corrected velocity must be the modelled profile at "
                                   "the first interior cell's distance"));
    PetscCall(PicurvAssertBool((PetscBool)(corrected_speed < speed),
                               "a log-law wall correction must slow the first interior cell "
                               "relative to the reference cell"));
 
    /* The ghost row carries no boundary velocity of its own. */
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Bcs.Ubcs, (const Cmpnts ***)&ubcs));
    PetscCall(PicurvAssertRealNear(0.0, ubcs[2][2][0].y, 0.0,
                                   "the wall ghost row must carry a zero boundary velocity"));
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Bcs.Ubcs, (const Cmpnts ***)&ubcs));
 
    PetscCall(VecDestroy(&user->Friction_Velocity));
    PetscCall(VecDestroy(&user->lFriction_Velocity));
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestWallFunctionModelSelectorChoosesTheLaw()

static PetscErrorCode TestWallFunctionModelSelectorChoosesTheLaw ( void  )

static

Tests that the model selector actually changes which law is applied.

The selector was previously validated in the configuration layer and then dropped, so every case ran the log law whatever it asked for. Distinct laws must produce distinct corrections from identical inputs, or the selector is decorative again.

Definition at line 1440 of file test_boundaries.c.

{
    SimCtx      *simCtx = NULL;
    UserCtx     *user = NULL;
    Cmpnts    ***ucat = NULL;
    PetscReal    corrected[4] = {0.0, 0.0, 0.0, 0.0};
    const PetscInt models[3] = {WALL_FUNCTION_LOG_LAW, WALL_FUNCTION_WERNER, WALL_FUNCTION_CABOT};
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    simCtx->ren = 1000.0;
    simCtx->wall_roughness_height = 1.0e-16;
    PetscCall(SeedWallFrictionVelocityStorage(user));
    user->boundary_faces[BC_FACE_NEG_X].mathematical_type = WALL;
 
    for (PetscInt m = 0; m < 3; ++m) {
        PetscCall(SeedWallFunctionFixture(user, 1.0));
        simCtx->wallfunction = models[m];
        PetscCall(ApplyWallFunction(user));
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
        corrected[m] = ucat[2][2][1].y;
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
        PetscCall(PicurvAssertBool((PetscBool)(corrected[m] > 0.0 && corrected[m] < 1.0),
                                   "every wall model must slow the first interior cell"));
    }
 
    /* Each law has its own profile, so no two may coincide. */
    PetscCall(PicurvAssertBool((PetscBool)(PetscAbsReal(corrected[0] - corrected[1]) > 1.0e-9),
                               "the Werner-Wengle model must differ from the log law"));
    PetscCall(PicurvAssertBool((PetscBool)(PetscAbsReal(corrected[0] - corrected[2]) > 1.0e-9),
                               "the Cabot model must differ from the log law"));
 
    PetscCall(VecDestroy(&user->Friction_Velocity));
    PetscCall(VecDestroy(&user->lFriction_Velocity));
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestApplyWallFunctionSkipsNonWallFaces()

static PetscErrorCode TestApplyWallFunctionSkipsNonWallFaces ( void  )

static

Tests that a face which is not a wall is left alone by the wall function.

Definition at line 1481 of file test_boundaries.c.

{
    SimCtx      *simCtx = NULL;
    UserCtx     *user = NULL;
    Cmpnts    ***ucat = NULL;
    PetscReal    before = 0.0, after = 0.0;
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    simCtx->ren = 1000.0;
    simCtx->wallfunction = 1;
    simCtx->wall_roughness_height = 1.0e-16;
    PetscCall(SeedWallFrictionVelocityStorage(user));
    /* Every face left at its default, so none is a wall. */
    PetscCall(SeedWallFunctionFixture(user, 1.0));
 
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
    before = ucat[2][2][1].y;
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
    PetscCall(ApplyWallFunction(user));
    PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
    after = ucat[2][2][1].y;
    PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucat, (const Cmpnts ***)&ucat));
 
    PetscCall(PicurvAssertRealNear(before, after, 0.0,
                                   "a face that is not a wall must not be corrected"));
 
    PetscCall(VecDestroy(&user->Friction_Velocity));
    PetscCall(VecDestroy(&user->lFriction_Velocity));
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestWallModelFrictionVelocityRespondsToTheFlow()

static PetscErrorCode TestWallModelFrictionVelocityRespondsToTheFlow ( void  )

static

Tests that every wall law's friction velocity responds to the flow it is given.

A law that ignores the resolved velocity and reports a fixed number still slows the first cell and still differs from the other laws, so the selector test above passes over it. Driving a law at two different speeds and requiring the friction velocity to follow is what distinguishes a model from a constant.

All three laws are covered. Werner-Wengle is the reason this case exists: it used to report its Newton initial guess rather than solving, so every wall cell in a run carried the same u_tau while still slowing the cell and still differing from the other laws, which is everything the selector test above checks.

Definition at line 1527 of file test_boundaries.c.

{
    SimCtx        *simCtx = NULL;
    UserCtx       *user = NULL;
    PetscReal   ***utau = NULL;
    const PetscInt models[3] = {WALL_FUNCTION_LOG_LAW, WALL_FUNCTION_WERNER,
                                WALL_FUNCTION_CABOT};
    const char    *names[3] = {"the log law must respond to the flow it is given",
                               "Werner-Wengle must respond to the flow it is given",
                               "Cabot must respond to the flow it is given"};
    const PetscReal speeds[2] = {0.5, 2.0};
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    simCtx->ren = 1000.0;
    simCtx->wall_roughness_height = 1.0e-16;
    PetscCall(SeedWallFrictionVelocityStorage(user));
    user->boundary_faces[BC_FACE_NEG_X].mathematical_type = WALL;
 
    for (PetscInt m = 0; m < 3; ++m) {
        PetscReal stored[2] = {0.0, 0.0};
 
        simCtx->wallfunction = models[m];
        for (PetscInt s = 0; s < 2; ++s) {
            PetscCall(SeedWallFunctionFixture(user, speeds[s]));
            PetscCall(ApplyWallFunction(user));
            PetscCall(DMDAVecGetArrayRead(user->da, user->Friction_Velocity,
                                          (const PetscReal ***)&utau));
            stored[s] = utau[2][2][1];
            PetscCall(DMDAVecRestoreArrayRead(user->da, user->Friction_Velocity,
                                              (const PetscReal ***)&utau));
            PetscCall(PicurvAssertBool((PetscBool)(stored[s] > 0.0),
                                       "a wall law must report a positive friction velocity"));
        }
 
        /* Faster flow, more wall stress. Any law that reports its initial guess, or any
           other constant, fails here whatever that constant is. */
        PetscCall(PicurvAssertBool((PetscBool)(stored[1] > stored[0] * 1.05),
                                   names[m]));
    }
 
    PetscCall(VecDestroy(&user->Friction_Velocity));
    PetscCall(VecDestroy(&user->lFriction_Velocity));
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestWallFrictionVelocityStorageMatchesTheCatalog()

static PetscErrorCode TestWallFrictionVelocityStorageMatchesTheCatalog ( void  )

static

Tests that the friction-velocity storage matches what the field catalog declares.

The wall model reaches its storage through a DM, and PETSc only rejects a mismatched DM/vector pairing at the moment of access - inside a branch that runs only when wall functions are enabled on a case that has a WALL face. That is how storage created from a dof-3 DM and read through a dof-1 one stayed latent. Comparing the allocation against the catalog descriptor turns that into a failure here.

Definition at line 1583 of file test_boundaries.c.

{
    SimCtx                *simCtx = NULL;
    UserCtx               *user = NULL;
    const FieldDescriptor *descriptor = NULL;
    FieldView              view;
    PetscInt               global_size = 0, local_size = 0;
    PetscInt               xm = 0, ym = 0, zm = 0, gxm = 0, gym = 0, gzm = 0;
    PetscReal           ***utau = NULL;
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    PetscCall(SeedWallFrictionVelocityStorage(user));
 
    PetscCall(FieldGetDescriptor(FIELD_ID_U_TAU, &descriptor));
    PetscCall(PicurvAssertBool((PetscBool)(descriptor->dof == 1),
                               "the friction velocity is one scalar per cell"));
    PetscCall(PicurvAssertBool((PetscBool)(descriptor->dm_kind == FIELD_DM_DA),
                               "a scalar cell field belongs on the dof-1 DM"));
 
    PetscCall(FieldGetView(user, FIELD_ID_U_TAU, &view));
    PetscCall(PicurvAssertBool((PetscBool)(view.global_vec == user->Friction_Velocity &&
                                           view.local_vec == user->lFriction_Velocity),
                               "the catalog must resolve to the vectors setup allocated"));
 
    /* The sizes the DM implies, against the sizes the vectors actually have. */
    PetscCall(DMDAGetCorners(user->da, NULL, NULL, NULL, &xm, &ym, &zm));
    PetscCall(DMDAGetGhostCorners(user->da, NULL, NULL, NULL, &gxm, &gym, &gzm));
    PetscCall(VecGetLocalSize(user->Friction_Velocity, &global_size));
    PetscCall(VecGetLocalSize(user->lFriction_Velocity, &local_size));
    PetscCall(PicurvAssertBool((PetscBool)(global_size == xm * ym * zm * descriptor->dof),
                               "the global friction velocity must match the owned extent"));
    PetscCall(PicurvAssertBool((PetscBool)(local_size == gxm * gym * gzm * descriptor->dof),
                               "the local friction velocity must match the ghosted extent"));
 
    /* The access the wall model itself performs. A mismatched pairing fails right here. */
    PetscCall(DMDAVecGetArray(user->da, user->Friction_Velocity, &utau));
    PetscCall(DMDAVecRestoreArray(user->da, user->Friction_Velocity, &utau));
 
    PetscCall(VecDestroy(&user->Friction_Velocity));
    PetscCall(VecDestroy(&user->lFriction_Velocity));
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestWallModelDiagnosticsRecordTheCorrectedCells()

static PetscErrorCode TestWallModelDiagnosticsRecordTheCorrectedCells ( void  )

static

Tests that the wall pass records the statistics its diagnostics report.

y+ cannot be recovered downstream, so it is formed where the wall distance is still in hand. This checks that what is accumulated there is the same quantity the reported column claims to be, and that a pass over no wall face records nothing rather than recording zeros.

Definition at line 1636 of file test_boundaries.c.

{
    SimCtx      *simCtx = NULL;
    UserCtx     *user = NULL;
    PetscReal    kinematic_viscosity = 0.0;
    PetscReal    expected_y_plus = 0.0;
    const PetscReal distance_first = 0.5;
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 6, 6, 6));
    simCtx->ren = 1000.0;
    simCtx->wall_roughness_height = 1.0e-16;
    kinematic_viscosity = 1.0 / simCtx->ren;
    PetscCall(SeedWallFrictionVelocityStorage(user));
    user->boundary_faces[BC_FACE_NEG_X].mathematical_type = WALL;
    PetscCall(SeedWallFunctionFixture(user, 1.0));
 
    /* Disabled: nothing is corrected, so nothing is recorded. */
    simCtx->wallfunction = 0;
    PetscCall(ApplyWallFunction(user));
    PetscCall(PicurvAssertBool((PetscBool)(user->wall_diagnostics.cells == 0),
                               "a disabled wall model must record no corrected cells"));
 
    simCtx->wallfunction = 1;
    PetscCall(ApplyWallFunction(user));
    PetscCall(PicurvAssertBool((PetscBool)(user->wall_diagnostics.cells > 0),
                               "an enabled wall model must record the cells it corrected"));
 
    /* One face, so the sums are one cell's values and y+ is exactly u_tau * y / nu.
       Reproducing it from the stored friction velocity is what makes a wrong distance
       or a wrong viscosity in the accumulation visible. */
    expected_y_plus = user->wall_diagnostics.friction_velocity_sum * distance_first /
                      kinematic_viscosity / (PetscReal)user->wall_diagnostics.cells;
    PetscCall(PicurvAssertRealNear(user->wall_diagnostics.y_plus_sum /
                                       (PetscReal)user->wall_diagnostics.cells,
                                   expected_y_plus, 1.0e-9,
                                   "y+ must be formed from the same distance the "
                                   "correction used"));
    PetscCall(PicurvAssertRealNear(user->wall_diagnostics.wall_distance_sum /
                                       (PetscReal)user->wall_diagnostics.cells,
                                   distance_first, 1.0e-9,
                                   "the recorded wall distance must be the first cell's"));
    PetscCall(PicurvAssertBool((PetscBool)(user->wall_diagnostics.friction_velocity_min <=
                                           user->wall_diagnostics.friction_velocity_max),
                               "the recorded extrema must bracket each other"));
 
    /* A second pass must describe that pass alone, not both. */
    {
        const PetscInt first_pass_cells = user->wall_diagnostics.cells;
 
        PetscCall(ApplyWallFunction(user));
        PetscCall(PicurvAssertBool((PetscBool)(user->wall_diagnostics.cells == first_pass_cells),
                                   "each pass must replace the previous pass's sample set"));
    }
 
    PetscCall(VecDestroy(&user->Friction_Velocity));
    PetscCall(VecDestroy(&user->lFriction_Velocity));
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

TestDummyCellsFollowEachFieldsBoundaryRule()

static PetscErrorCode TestDummyCellsFollowEachFieldsBoundaryRule ( void  )

static

Checks each cell-centred field's dummy-cell rule, faces then edges.

Interior cells carry index-unique sentinels and every dummy cell starts poisoned, so a wrong offset or a missed face shows. Velocity dummies must put the face average on Ubcs; pressure and the particle scalar copy the adjacent cell (zero normal gradient); each edge dummy is the mean of its two face neighbours; and a field with no rule, or one that is not cell-centred, is refused.

Definition at line 1705 of file test_boundaries.c.

{
    SimCtx       *simCtx = NULL;
    UserCtx      *user = NULL;
    const FieldId scalars[] = {FIELD_ID_P, FIELD_ID_PSI};
 
    PetscFunctionBeginUser;
    PetscCall(PicurvCreateMinimalContexts(&simCtx, &user, 5, 6, 7));
    const DMDALocalInfo info = user->info;
    const PetscInt mx = info.mx, my = info.my, mz = info.mz;
    PetscCall(VecSet(user->Bcs.Ubcs, 0.0));
    {
        PetscReal ***p = NULL, ***psi = NULL;
        Cmpnts    ***ucat = NULL, ***ubcs = NULL;
 
        PetscCall(DMDAVecGetArray(user->da, user->P, &p));
        PetscCall(DMDAVecGetArray(user->da, user->Psi, &psi));
        PetscCall(DMDAVecGetArray(user->fda, user->Ucat, &ucat));
        PetscCall(DMDAVecGetArray(user->fda, user->Bcs.Ubcs, &ubcs));
        for (PetscInt k = info.zs; k < info.zs + info.zm; k++)
        for (PetscInt j = info.ys; j < info.ys + info.ym; j++)
        for (PetscInt i = info.xs; i < info.xs + info.xm; i++) {
            const PetscBool interior = (PetscBool)(i > 0 && i < mx - 1 && j > 0 && j < my - 1 && k > 0 && k < mz - 1);
            const PetscReal sentinel = 100.0 * k + 10.0 * j + i;
 
            p[k][j][i] = interior ? sentinel : -999.0;
            psi[k][j][i] = interior ? -sentinel : -999.0;
            ucat[k][j][i].x = interior ? sentinel : -999.0;
            ucat[k][j][i].y = interior ? 2.0 * sentinel : -999.0;
            ucat[k][j][i].z = interior ? 3.0 * sentinel : -999.0;
            ubcs[k][j][i] = (Cmpnts){1.0, 2.0, 3.0};
        }
        PetscCall(DMDAVecRestoreArray(user->fda, user->Bcs.Ubcs, &ubcs));
        PetscCall(DMDAVecRestoreArray(user->fda, user->Ucat, &ucat));
        PetscCall(DMDAVecRestoreArray(user->da, user->Psi, &psi));
        PetscCall(DMDAVecRestoreArray(user->da, user->P, &p));
    }
 
    PetscCall(UpdateDummyCells(user, FIELD_ID_UCAT));
    PetscCall(UpdateCornerNodes(user, FIELD_ID_UCAT));
    for (size_t n = 0; n < sizeof(scalars) / sizeof(scalars[0]); ++n) {
        PetscCall(UpdateDummyCells(user, scalars[n]));
        PetscCall(UpdateCornerNodes(user, scalars[n]));
    }
 
    for (size_t n = 0; n < sizeof(scalars) / sizeof(scalars[0]); ++n) {
        FieldView  view;
        PetscReal ***f = NULL;
 
        PetscCall(FieldGetView(user, scalars[n], &view));
        PetscCall(DMDAVecGetArrayRead(user->da, view.global_vec, &f));
        for (PetscInt k = 1; k < mz - 1; k++) for (PetscInt j = 1; j < my - 1; j++) {
            PetscCall(PicurvAssertRealNear(f[k][j][1], f[k][j][0], 0.0, "-x dummy copies the adjacent cell"));
            PetscCall(PicurvAssertRealNear(f[k][j][mx - 2], f[k][j][mx - 1], 0.0, "+x dummy copies the adjacent cell"));
        }
        for (PetscInt k = 1; k < mz - 1; k++) for (PetscInt i = 1; i < mx - 1; i++) {
            PetscCall(PicurvAssertRealNear(f[k][1][i], f[k][0][i], 0.0, "-y dummy copies the adjacent cell"));
            PetscCall(PicurvAssertRealNear(f[k][my - 2][i], f[k][my - 1][i], 0.0, "+y dummy copies the adjacent cell"));
        }
        for (PetscInt j = 1; j < my - 1; j++) for (PetscInt i = 1; i < mx - 1; i++) {
            PetscCall(PicurvAssertRealNear(f[1][j][i], f[0][j][i], 0.0, "-z dummy copies the adjacent cell"));
            PetscCall(PicurvAssertRealNear(f[mz - 2][j][i], f[mz - 1][j][i], 0.0, "+z dummy copies the adjacent cell"));
        }
        for (PetscInt j = 0; j < my; j++) {
            PetscCall(PicurvAssertRealNear(0.5 * (f[1][j][0] + f[0][j][1]), f[0][j][0], 0.0,
                                           "a (-z,-x) edge dummy is the mean of its face neighbours"));
        }
        PetscCall(PicurvAssertBool((PetscBool)(f[0][0][0] != -999.0 && f[mz - 1][my - 1][mx - 1] != -999.0),
                                   "corner dummies are written"));
        PetscCall(DMDAVecRestoreArrayRead(user->da, view.global_vec, &f));
    }
    {
        const Cmpnts ***ucat = NULL;
 
        PetscCall(DMDAVecGetArrayRead(user->fda, user->Ucat, &ucat));
        for (PetscInt k = 1; k < mz - 1; k++) for (PetscInt j = 1; j < my - 1; j++) {
            PetscCall(PicurvAssertRealNear(1.0, 0.5 * (ucat[k][j][0].x + ucat[k][j][1].x), 1.0e-12,
                                           "a velocity dummy puts the face average on Ubcs"));
            PetscCall(PicurvAssertRealNear(3.0, 0.5 * (ucat[k][j][mx - 1].z + ucat[k][j][mx - 2].z), 1.0e-12,
                                           "a velocity dummy puts the face average on Ubcs"));
        }
        PetscCall(DMDAVecRestoreArrayRead(user->fda, user->Ucat, &ucat));
    }
 
    {
        PetscErrorCode refused;
 
        PetscCall(PetscPushErrorHandler(PetscIgnoreErrorHandler, NULL));
        refused = UpdateDummyCells(user, FIELD_ID_NVERT);
        PetscCall(PetscPopErrorHandler());
        PetscCall(PicurvAssertIntEqual(PETSC_ERR_SUP, refused, "a field with no boundary rule is refused"));
        PetscCall(PetscPushErrorHandler(PetscIgnoreErrorHandler, NULL));
        refused = UpdateDummyCells(user, FIELD_ID_COORDINATES);
        PetscCall(PetscPopErrorHandler());
        PetscCall(PicurvAssertIntEqual(PETSC_ERR_ARG_WRONG, refused, "a field that is not cell-centred is refused"));
    }
 
    PetscCall(PicurvDestroyMinimalContexts(&simCtx, &user));
    PetscFunctionReturn(0);
}

Here is the call graph for this function:

Here is the caller graph for this function:

main()

int main	(	int	argc,
		char **	argv
	)

Runs the unit-boundaries PETSc test binary.

Definition at line 1810 of file test_boundaries.c.

{
    PetscErrorCode ierr;
    const PicurvTestCase cases[] = {
        {"can-rank-service-face-matches-inlet-when-defined", TestCanRankServiceFaceMatchesInletWhenDefined},
        {"can-rank-service-inlet-face-requires-definition", TestCanRankServiceInletFaceRequiresDefinition},
        {"boundary-condition-factory-assignments", TestBoundaryConditionFactoryAssignments},
        {"boundary-condition-factory-implemented-handler-matrix", TestBoundaryConditionFactoryImplementedHandlerMatrix},
        {"boundary-condition-factory-rejects-unsupported-handler", TestBoundaryConditionFactoryRejectsUnsupportedHandler},
        {"deterministic-face-grid-location-matrix", TestGetDeterministicFaceGridLocationFaceMatrix},
        {"random-inlet-face-location-matrix", TestGetRandomCellAndLogicalCoordsOnInletFaceMatrix},
        {"wall-no-slip-handler-face-matrix", TestWallNoSlipHandlerFaceMatrix},
        {"dummy-cells-follow-each-fields-boundary-rule", TestDummyCellsFollowEachFieldsBoundaryRule},
        {"inlet-constant-velocity-handler-behavior", TestInletConstantVelocityHandlerBehavior},
        {"inlet-common-application-curvilinear-metric", TestInletCommonApplicationCurvilinearMetric},
        {"inlet-constant-velocity-handler-face-matrix", TestInletConstantVelocityHandlerFaceMatrix},
        {"inlet-parabolic-profile-handler-behavior", TestInletParabolicProfileHandlerBehavior},
        {"inlet-parabolic-profile-handler-face-matrix", TestInletParabolicProfileHandlerFaceMatrix},
        {"inlet-parabolic-profile-uniform-along-periodic-axis", TestInletParabolicProfileUniformAlongPeriodicAxis},
        {"inlet-profile-from-file-handler-behavior", TestInletProfileFromFileHandlerBehavior},
        {"outlet-conservation-handler-behavior", TestOutletConservationHandlerBehavior},
        {"outlet-conservation-handler-face-matrix", TestOutletConservationHandlerFaceMatrix},
        {"finalize-post-projection-cell-fields-non-periodic", TestFinalizePostProjectionCellFieldsNonPeriodic},
        {"enforce-rhs-bcs-matches-row-classification", TestEnforceRHSBoundaryConditionsMatchesRowClassification},
        {"apply-wall-function-matches-the-log-law", TestApplyWallFunctionMatchesTheLogLawItDispatches},
        {"apply-wall-function-skips-non-wall-faces", TestApplyWallFunctionSkipsNonWallFaces},
        {"wall-function-model-selector-chooses-the-law", TestWallFunctionModelSelectorChoosesTheLaw},
        {"wall-model-friction-velocity-responds-to-the-flow", TestWallModelFrictionVelocityRespondsToTheFlow},
        {"wall-friction-velocity-storage-matches-the-catalog", TestWallFrictionVelocityStorageMatchesTheCatalog},
        {"wall-model-diagnostics-record-the-corrected-cells", TestWallModelDiagnosticsRecordTheCorrectedCells},
    };
 
    ierr = PetscInitialize(&argc, &argv, NULL, "PICurv boundary tests");
    if (ierr) {
        return (int)ierr;
    }
 
    ierr = PicurvRunTests("unit-boundaries", cases, sizeof(cases) / sizeof(cases[0]));
    if (ierr) {
        PetscFinalize();
        return (int)ierr;
    }
 
    ierr = PetscFinalize();
    return (int)ierr;
}

Here is the call graph for this function: