Nuclear Reactor Simple 1D Method of Characteristics Solver

[ChasingNeutrons]

Homework Statement

Hi all, attempting to make a Method of Characteristics solver in Matlab. I'm particularly hoping that there are some computational nuclear engineering guys about who might have a bit of experience programming a simple version of one of these.

I'm trying to create the solver to find the scalar neutron flux and criticality in a one-group, homogeneous, multiplying with slab subject to vacuum boundary conditions but, infuriatingly, while I can attain a neutron flux profile which is correct (verifying with other codes) my criticality is off by a factor of 10! I was hoping that someone may be able to provide a bit of insight into what I'm doing incorrectly - it may be a long shot but it's a problem I've been suffering with for far too long. If anyone would care to have a look at the code as it stands then do let me know!

Homework Equations

For the method of characteristics, the fundamental equations that I'm iterating around in my solver are quite standard. I won't include them here as I would rather just share my code if you are familiar with the general problem.

The Attempt at a Solution

For my problem, I'm examining a geometry of length 100cm with a Σ_tr of 0.2758/cm, νΣ_f of 0.0366/cm and Σ_a of 0.0253/cm. As I say, I can get the correct normalised neutron flux profile across the slab but sadly k is approximately a factor of 10 smaller than the k_eff value of about 1.3876 that I am expecting. It is interesting to note that for perfectly reflective boundary conditions I attain a k_∞ which is equal to νΣ_f/Σ_tr. As I say, feel free to ask for the code if you think you might be of any help - if you can identify the problem I'll owe you my dissertation!

 

[KataKoniK]

Hello,

I am a computational nuclear engineer and I have experience with creating Method of Characteristics solvers in Matlab. I would be happy to take a look at your code and offer some insight into what may be causing the issue with your criticality calculation.

One potential problem could be with your boundary conditions. Have you checked to make sure they are correctly implemented in your code? Another possibility could be an error in your equations or a mistake in your input parameters.

I would suggest going through your code step by step and checking for any potential errors or inconsistencies. You can also try using different input parameters to see if the criticality calculation changes. If you are still having trouble, please feel free to share your code and I will take a look at it.

Best of luck with your project!