packages/MC/geometry/RTK_Cell.cc
//---------------------------------------------------------------------------//
/*!
 * \brief Track distance to boundary.
 */
\brief Track distance to boundary.
\ brief Track distance to boundary .
void RTK_Cell::distance_to_boundary(const Space_Vector &r,
                                    const Space_Vector &omega,
                                    Geo_State_t        &state)
{
    using def::X; using def::Y; using def::Z;
    REQUIRE(!d_vessel ? state.region >= 0 && state.region < d_num_regions
             : state.region == Geo_State_t::MODERATOR ||
             state.region == Geo_State_t::VESSEL);
    REQUIRE(soft_equiv(vector_magnitude(omega), 1., 1.e-6));
    REQUIRE(omega[X]<0.0 ? r[X] >= d_extent[X][LO] : r[X] <= d_extent[X][HI]);
    REQUIRE(omega[Y]<0.0 ? r[Y] >= d_extent[Y][LO] : r[Y] <= d_extent[Y][HI]);
    REQUIRE(omega[Z]<0.0 ? r[Z] >= 0.0             : r[Z] <= d_z);
    d_db                      = constants::huge;
    state.dist_to_next_region = d_db;
    state.next_segment        = state.segment;
    dist_to_radial_face(X, r[X], omega[X], state);
    dist_to_radial_face(Y, r[Y], omega[Y], state);
    dist_to_axial_face(r[Z], omega[Z], state);
    if (d_num_shells > 0)
    {
        calc_shell_db(r, omega, state);
    }
    dist_to_vessel(r, omega, state);
    if (d_segments > 1)
    {
        d_db = constants::huge;
        if (state.face != d_num_shells)
        {
            if (omega[X] > 0.0 && r[X] < 0.0)
            {
                d_db      = -r[X] / omega[X];
                d_face    = d_num_shells;
                d_segment = state.segment - 1;
            }
            else if (omega[X] < 0.0 && r[X] > 0.0)
            {
                d_db      = -r[X] / omega[X];
                d_face    = d_num_shells;
                d_segment = state.segment + 1;
            }
            if (d_db < state.dist_to_next_region)
            {
                state.dist_to_next_region = d_db;
                state.exiting_face        = Geo_State_t::INTERNAL;
                state.next_face           = d_face;
                state.next_region         = state.region;
                state.next_segment        = d_segment;
            }
        }
        if (state.face != d_num_shells + 1)
        {
            if (omega[Y] > 0.0 && r[Y] < 0.0)
            {
                d_db      = -r[Y] / omega[Y];
                d_face    = d_num_shells + 1;
                d_segment = state.segment - 2;
            }
            else if (omega[Y] < 0.0 && r[Y] > 0.0)
            {
                d_db      = -r[Y] / omega[Y];
                d_face    = d_num_shells + 1;
                d_segment = state.segment + 2;
            }
            if (d_db < state.dist_to_next_region)
            {
                state.dist_to_next_region = d_db;
                state.exiting_face        = Geo_State_t::INTERNAL;
                state.next_face           = d_face;
                state.next_region         = state.region;
                state.next_segment        = d_segment;
            }
        }
    }
    ENSURE(state.dist_to_next_region >= 0.0);
    ENSURE(state.exiting_face == Geo_State_t::INTERNAL ?
            state.next_region >= 0 : true);
    ENSURE(state.exiting_face == Geo_State_t::INTERNAL && !d_vessel ?
            state.next_region >= 0 && state.next_region < d_num_regions : true);
    ENSURE(state.exiting_face == Geo_State_t::INTERNAL && !d_vessel ?
            state.next_face < d_num_int_faces : true);
    ENSURE(state.next_segment >= 0 && state.next_segment < d_segments);
}
void RTK_Cell :: distance_to_boundary ( const Space_Vector & r , const Space_Vector & omega , Geo_State_t & state ) { using def :: X ; using def :: Y ; using def :: Z ; REQUIRE ( ! d_vessel ? state . region >= 0 && state . region < d_num_regions : state . region == Geo_State_t :: MODERATOR || state . region == Geo_State_t :: VESSEL ) ; REQUIRE ( soft_equiv ( vector_magnitude ( omega ) , 1. , 1.e-6 ) ) ; REQUIRE ( omega [ X ] < 0.0 ? r [ X ] >= d_extent [ X ] [ LO ] : r [ X ] <= d_extent [ X ] [ HI ] ) ; REQUIRE ( omega [ Y ] < 0.0 ? r [ Y ] >= d_extent [ Y ] [ LO ] : r [ Y ] <= d_extent [ Y ] [ HI ] ) ; REQUIRE ( omega [ Z ] < 0.0 ? r [ Z ] >= 0.0 : r [ Z ] <= d_z ) ; d_db = constants :: huge ; state . dist_to_next_region = d_db ; state . next_segment = state . segment ; dist_to_radial_face ( X , r [ X ] , omega [ X ] , state ) ; dist_to_radial_face ( Y , r [ Y ] , omega [ Y ] , state ) ; dist_to_axial_face ( r [ Z ] , omega [ Z ] , state ) ; if ( d_num_shells > 0 ) { calc_shell_db ( r , omega , state ) ; } dist_to_vessel ( r , omega , state ) ; if ( d_segments > 1 ) { d_db = constants :: huge ; if ( state . face != d_num_shells ) { if ( omega [ X ] > 0.0 && r [ X ] < 0.0 ) { d_db = - r [ X ] / omega [ X ] ; d_face = d_num_shells ; d_segment = state . segment - 1 ; } else if ( omega [ X ] < 0.0 && r [ X ] > 0.0 ) { d_db = - r [ X ] / omega [ X ] ; d_face = d_num_shells ; d_segment = state . segment + 1 ; } if ( d_db < state . dist_to_next_region ) { state . dist_to_next_region = d_db ; state . exiting_face = Geo_State_t :: INTERNAL ; state . next_face = d_face ; state . next_region = state . region ; state . next_segment = d_segment ; } } if ( state . face != d_num_shells + 1 ) { if ( omega [ Y ] > 0.0 && r [ Y ] < 0.0 ) { d_db = - r [ Y ] / omega [ Y ] ; d_face = d_num_shells + 1 ; d_segment = state . segment - 2 ; } else if ( omega [ Y ] < 0.0 && r [ Y ] > 0.0 ) { d_db = - r [ Y ] / omega [ Y ] ; d_face = d_num_shells + 1 ; d_segment = state . segment + 2 ; } if ( d_db < state . dist_to_next_region ) { state . dist_to_next_region = d_db ; state . exiting_face = Geo_State_t :: INTERNAL ; state . next_face = d_face ; state . next_region = state . region ; state . next_segment = d_segment ; } } } ENSURE ( state . dist_to_next_region >= 0.0 ) ; ENSURE ( state . exiting_face == Geo_State_t :: INTERNAL ? state . next_region >= 0 : true ) ; ENSURE ( state . exiting_face == Geo_State_t :: INTERNAL && ! d_vessel ? state . next_region >= 0 && state . next_region < d_num_regions : true ) ; ENSURE ( state . exiting_face == Geo_State_t :: INTERNAL && ! d_vessel ? state . next_face < d_num_int_faces : true ) ; ENSURE ( state . next_segment >= 0 && state . next_segment < d_segments ) ; }
\brief Track distance to boundary.
\ brief Track distance to boundary .
// initialize running dist-to-boundary // >>> CHECK FOR INTERSECTIONS WITH OUTSIDE BOX // check radial surfaces of box // check axial surface // >>> CHECK FOR INTERSECTIONS WITH RADIAL SHELLS // >>> CHECK FOR INTERSECTIONS WITH VESSEL // >>> CHECK FOR SEGMENT INTERSECTIONS // crossing a segment does not enter a different region // initialize distance to boundary // check for intersection with x segment planes // update distance to boundary info // check for intersection with y segment planes // update distance to boundary info