statistics_moments.c File Reference

Weighted centered-moment kernels for the field-statistics pipeline. More...

#include "statistics_moments.h"

Include dependency graph for statistics_moments.c:

Go to the source code of this file.

Functions
void	PicurvMomentStateReset (PicurvMomentState *state)
	Implementation of PicurvMomentStateReset().
void	PicurvCoMomentStateReset (PicurvCoMomentState *state)
	Implementation of PicurvCoMomentStateReset().
PetscErrorCode	PicurvMomentStateUpdate (PicurvMomentState *state, PetscReal value, PetscReal weight)
	Implementation of PicurvMomentStateUpdate().
PetscErrorCode	PicurvCoMomentStateUpdate (PicurvCoMomentState *state, PetscReal value_x, PetscReal value_y, PetscReal weight)
	Implementation of PicurvCoMomentStateUpdate().
PetscErrorCode	PicurvMomentStateMerge (PicurvMomentState *result, const PicurvMomentState *a, const PicurvMomentState *b)
	Implementation of PicurvMomentStateMerge().
PetscErrorCode	PicurvCoMomentStateMerge (PicurvCoMomentState *result, const PicurvCoMomentState *a, const PicurvCoMomentState *b)
	Implementation of PicurvCoMomentStateMerge().
PetscReal	PicurvMomentStateVariance (const PicurvMomentState *state)
	Implementation of PicurvMomentStateVariance().
PetscReal	PicurvCoMomentStateCovariance (const PicurvCoMomentState *state)
	Implementation of PicurvCoMomentStateCovariance().
PetscReal	PicurvMomentStateEffectiveCount (const PicurvMomentState *state)
	Implementation of PicurvMomentStateEffectiveCount().

Detailed Description

Weighted centered-moment kernels for the field-statistics pipeline.

Full API contract is documented with the declarations in include/statistics_moments.h.

Definition in file statistics_moments.c.

Function Documentation

PicurvMomentStateReset()

void PicurvMomentStateReset

(

PicurvMomentState *

state

)

Implementation of PicurvMomentStateReset().

Resets a scalar moment accumulator to the empty state.

See also

PicurvMomentStateReset()

Definition at line 15 of file statistics_moments.c.

{
    if (state == NULL) return;
    state->count = 0.0;
    state->weight = 0.0;
    state->weight_sq = 0.0;
    state->mean = 0.0;
    state->m2 = 0.0;
}

Here is the caller graph for this function:

PicurvCoMomentStateReset()

void PicurvCoMomentStateReset

(

PicurvCoMomentState *

state

)

Implementation of PicurvCoMomentStateReset().

Resets a co-moment accumulator to the empty state.

See also

PicurvCoMomentStateReset()

Definition at line 29 of file statistics_moments.c.

{
    if (state == NULL) return;
    state->count = 0.0;
    state->weight = 0.0;
    state->weight_sq = 0.0;
    state->mean_x = 0.0;
    state->mean_y = 0.0;
    state->cm = 0.0;
}

Here is the caller graph for this function:

PicurvMomentStateUpdate()

PetscErrorCode PicurvMomentStateUpdate	(	PicurvMomentState *	state,
		PetscReal	value,
		PetscReal	weight
	)

Implementation of PicurvMomentStateUpdate().

Applies one weighted sample to a scalar moment accumulator.

See also

PicurvMomentStateUpdate()

Definition at line 44 of file statistics_moments.c.

{
    PetscReal new_weight = 0.0;
    PetscReal delta = 0.0;
 
    PetscFunctionBeginUser;
    PetscCheck(state != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL,
               "Moment accumulator is required.");
    PetscCheck(weight > 0.0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE,
               "Sample weight must be strictly positive, got %g.", (double)weight);
 
    new_weight = state->weight + weight;
    delta = value - state->mean;
    state->mean += (weight / new_weight) * delta;
    /* Uses the post-update mean deliberately: this is the pairing that keeps the
     * centered sum stable when the mean is large relative to the fluctuation. */
    state->m2 += weight * delta * (value - state->mean);
    state->weight = new_weight;
    state->weight_sq += weight * weight;
    state->count += 1.0;
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

PicurvCoMomentStateUpdate()

PetscErrorCode PicurvCoMomentStateUpdate	(	PicurvCoMomentState *	state,
		PetscReal	value_x,
		PetscReal	value_y,
		PetscReal	weight
	)

Implementation of PicurvCoMomentStateUpdate().

Applies one weighted paired sample to a co-moment accumulator.

See also

PicurvCoMomentStateUpdate()

Definition at line 71 of file statistics_moments.c.

{
    PetscReal new_weight = 0.0;
    PetscReal delta_x = 0.0;
 
    PetscFunctionBeginUser;
    PetscCheck(state != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL,
               "Co-moment accumulator is required.");
    PetscCheck(weight > 0.0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE,
               "Sample weight must be strictly positive, got %g.", (double)weight);
 
    new_weight = state->weight + weight;
    delta_x = value_x - state->mean_x;
    state->mean_x += (weight / new_weight) * delta_x;
    state->mean_y += (weight / new_weight) * (value_y - state->mean_y);
    /* Pre-update mean of x against post-update mean of y. With x == y this
     * reduces exactly to the M2 update, which the unit tests assert. */
    state->cm += weight * delta_x * (value_y - state->mean_y);
    state->weight = new_weight;
    state->weight_sq += weight * weight;
    state->count += 1.0;
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

PicurvMomentStateMerge()

PetscErrorCode PicurvMomentStateMerge	(	PicurvMomentState *	result,
		const PicurvMomentState *	a,
		const PicurvMomentState *	b
	)

Implementation of PicurvMomentStateMerge().

Merges two independently accumulated scalar states.

See also

PicurvMomentStateMerge()

Definition at line 100 of file statistics_moments.c.

{
    PetscReal total_weight = 0.0;
    PetscReal delta = 0.0;
    PicurvMomentState merged;
 
    PetscFunctionBeginUser;
    PetscCheck(result != NULL && a != NULL && b != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL,
               "Merge requires a destination and two source accumulators.");
 
    if (b->weight == 0.0) { *result = *a; PetscFunctionReturn(0); }
    if (a->weight == 0.0) { *result = *b; PetscFunctionReturn(0); }
 
    total_weight = a->weight + b->weight;
    delta = b->mean - a->mean;
    merged.count = a->count + b->count;
    merged.weight = total_weight;
    merged.weight_sq = a->weight_sq + b->weight_sq;
    merged.mean = a->mean + (b->weight / total_weight) * delta;
    merged.m2 = a->m2 + b->m2 + (a->weight * b->weight / total_weight) * delta * delta;
    *result = merged;
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

PicurvCoMomentStateMerge()

PetscErrorCode PicurvCoMomentStateMerge	(	PicurvCoMomentState *	result,
		const PicurvCoMomentState *	a,
		const PicurvCoMomentState *	b
	)

Implementation of PicurvCoMomentStateMerge().

Merges two independently accumulated co-moment states.

See also

PicurvCoMomentStateMerge()

Definition at line 129 of file statistics_moments.c.

{
    PetscReal total_weight = 0.0;
    PetscReal delta_x = 0.0;
    PetscReal delta_y = 0.0;
    PicurvCoMomentState merged;
 
    PetscFunctionBeginUser;
    PetscCheck(result != NULL && a != NULL && b != NULL, PETSC_COMM_SELF, PETSC_ERR_ARG_NULL,
               "Merge requires a destination and two source accumulators.");
 
    if (b->weight == 0.0) { *result = *a; PetscFunctionReturn(0); }
    if (a->weight == 0.0) { *result = *b; PetscFunctionReturn(0); }
 
    total_weight = a->weight + b->weight;
    delta_x = b->mean_x - a->mean_x;
    delta_y = b->mean_y - a->mean_y;
    merged.count = a->count + b->count;
    merged.weight = total_weight;
    merged.weight_sq = a->weight_sq + b->weight_sq;
    merged.mean_x = a->mean_x + (b->weight / total_weight) * delta_x;
    merged.mean_y = a->mean_y + (b->weight / total_weight) * delta_y;
    merged.cm = a->cm + b->cm + (a->weight * b->weight / total_weight) * delta_x * delta_y;
    *result = merged;
    PetscFunctionReturn(0);
}

Here is the caller graph for this function:

PicurvMomentStateVariance()

PetscReal PicurvMomentStateVariance

(

const PicurvMomentState *

state

)

Implementation of PicurvMomentStateVariance().

Returns the weighted variance M2/W, or zero when no weight accumulated.

See also

PicurvMomentStateVariance()

Definition at line 161 of file statistics_moments.c.

{
    if (state == NULL || state->weight == 0.0) return 0.0;
    return state->m2 / state->weight;
}

Here is the caller graph for this function:

PicurvCoMomentStateCovariance()

PetscReal PicurvCoMomentStateCovariance

(

const PicurvCoMomentState *

state

)

Implementation of PicurvCoMomentStateCovariance().

Returns the weighted covariance C/W, or zero when no weight accumulated.

See also

PicurvCoMomentStateCovariance()

Definition at line 171 of file statistics_moments.c.

{
    if (state == NULL || state->weight == 0.0) return 0.0;
    return state->cm / state->weight;
}

Here is the caller graph for this function:

PicurvMomentStateEffectiveCount()

PetscReal PicurvMomentStateEffectiveCount

(

const PicurvMomentState *

state

)

Implementation of PicurvMomentStateEffectiveCount().

Returns Kish effective sample size W^2/W2, or zero when no weight accumulated.

See also

PicurvMomentStateEffectiveCount()

Definition at line 181 of file statistics_moments.c.

{
    if (state == NULL || state->weight_sq == 0.0) return 0.0;
    return (state->weight * state->weight) / state->weight_sq;
}

Here is the caller graph for this function: