Variance Reduced Integrated Charge by gyalla · Pull Request #17 · sandialabs/mfpic

gyalla · 2026-02-27T14:26:26Z

No description provided.

libmfpic/DGEulerOperations.cpp

libmfpic/ParticleOperations.cpp

libmfpic/ParticleOperations.hpp

libmfpic/ParticleOperations.test.cpp

libmfpic/Euler.hpp

libmfpic/Euler.cpp

reuterb · 2026-02-27T18:50:25Z

libmfpic/Euler.cpp

 }

+double evaluateMaxwellian(const mfem::Vector& primitive_state, const mfem::Vector velocity,const Species& species,const int dim){
+  //TODO: Do we want a default value if sigma == 0? 


We shouldn't be able to get T<=0 so you add an error here if you'd like. I think we'd probably fail further upstream.

If error out when sigma<=0, several of our existing tests fail that use zero temperature. We can add a warning or perhaps return a NaN or negative probability

gyalla · 2026-03-19T16:22:58Z

libmfpic/ParticleOperations.test.cpp

+  return out;
+}
+
+TEST(ParticleOperations, ErrorBetweenLowFidelityChargeAndPICChargeReducesWithParticleCountForMaxwellian) {


@reuterb Take a look at this unit test and the following

This is checking that the PIC charge and the Euler charge match as N_part -> \infty? I think the test name is maybe a bit confusing in that regard. You could make this more rigorous by comparing to the expected value, using the expected error (related to the variance) as the tolerance for PIC.

reuterb · 2026-03-19T16:53:45Z

libmfpic/ParticleOperations.test.cpp

+}
+
+
+TEST(ParticleOperations, ErrorBetweenVarianceReducedChargeAndPICChargeReducesWithParticleCountForMaxwellian) {


same thing here. I think we could be more direct.

We know the VR operation is unbiased so we'd expected a value of $E[\phi]$ with a tolerance related to the error $\sqrt{Var[\phi(1-f/f)]/N}$.

I may not be fully grasping what we're trying to do here though.

This unit test and the last are simply checking that the error between the variance-reduced/low-fidelity charge and the standard PIC charge reduces with particle count. I wanted to make sure my variance reduced charge implementation was still converging towards the same value as the original PIC implementation.

I was holding off on computing variance and expectation of the estimator in a unit test since that would require multiple trials and that started to feel more like a regression test.

For a Maxwellian, the (1-fNR/f) factor is identical to 0, so the variance reduced charge is equal to the low fidelity charge. So this is just testing convergence of the PIC charge to the low-fidelity charge. If I use a kappa in the second unit test, that may be more meaningful.

This unit test and the last are simply checking that the error between the variance-reduced/low-fidelity charge and the standard PIC charge reduces with particle count. I wanted to make sure my variance reduced charge implementation was still converging towards the same value as the original PIC implementation.

I was holding off on computing variance and expectation of the estimator in a unit test since that would require multiple trials and that started to feel more like a regression test.

Isn't the variance reduced implementation exact in this case then? So basically these are just checking that the standard PIC charge converges?

Yes, the variance reduced implementation is exact, that's right. That is actually the first unit test I have in the file. If we use a kappa-distribution of particles and check that the variance-reduced charge and PIC charge still converge to the same answer, that seems like a good unit test for the new code.

Co-authored-by: reuterb <bwreute@sandia.gov>

…unit test; switch statement on dimension

…density

…ced charge, and PIC charge