crp/OpenMC
1
0
Fork 0
mirror of https://github.com/openmc-dev/openmc.git synced 2026-07-28 06:05:58 -04:00
Code Issues Projects Releases Packages Wiki Activity Actions 5

removed print statements. Works in serial as expected.

Browse source
This commit is contained in:
John Tramm 2019-11-18 22:28:19 +00:00
parent 3ef3ee9fcb
commit 2f9dba6a58
3 changed files with 37 additions and 32 deletions
Download patch file Download diff file Expand all files Collapse all files

4
examples/xml/pincell/settings.xml
View file

@ -3,8 +3,8 @@
<!-- Define how many particles to run and for how many batches -->
<run_mode>eigenvalue</run_mode>
<batches>100</batches>
<inactive>10</inactive>
<batches>10</batches>
<inactive>5</inactive>
<particles>1000</particles>
<!-- The starting source is a uniform distribution over the entire pin

30
src/physics.cpp
View file

@ -79,7 +79,7 @@ void collision(Particle* p)
void sample_neutron_reaction(Particle* p)
{
std::cout << "particle energy e = " << p->E_ << std::endl;
//std::cout << "particle energy e = " << p->E_ << std::endl;
// Sample a nuclide within the material
int i_nuclide = sample_nuclide(p);
@ -128,12 +128,12 @@ void sample_neutron_reaction(Particle* p)
}
if (!p->alive_) return;
std::cout << "before scatter particle energy e = " << p->E_ << " with nuclide = " << i_nuclide << std::endl;
//std::cout << "before scatter particle energy e = " << p->E_ << " with nuclide = " << i_nuclide << std::endl;
// Sample a scattering reaction and determine the secondary energy of the
// exiting neutron
scatter(p, i_nuclide);
std::cout << "after scatter particle energy e = " << p->E_ << std::endl;
//std::cout << "after scatter particle energy e = " << p->E_ << std::endl;
// Advance URR seed stream 'N' times after energy changes
if (p->E_ != p->E_last_) {
@ -602,7 +602,7 @@ void scatter(Particle* p, int i_nuclide)
// Determine temperature
double kT = nuc->multipole_ ? p->sqrtkT_*p->sqrtkT_ : nuc->kTs_[i_temp];
std::cout << "we are doing an elastic scatter" << std::endl;
//std::cout << "we are doing an elastic scatter" << std::endl;
// Perform collision physics for elastic scattering
elastic_scatter(i_nuclide, *nuc->reactions_[0], kT, p);
@ -615,7 +615,7 @@ void scatter(Particle* p, int i_nuclide)
// =======================================================================
// S(A,B) SCATTERING
std::cout << "we are doing an SAB scatter" << std::endl;
//std::cout << "we are doing an SAB scatter" << std::endl;
sab_scatter(i_nuclide, micro.index_sab, p);
p->event_mt_ = ELASTIC;
@ -623,7 +623,7 @@ void scatter(Particle* p, int i_nuclide)
}
if (!sampled) {
std::cout << "we are doing an Inelastic scatter" << std::endl;
//std::cout << "we are doing an Inelastic scatter" << std::endl;
// =======================================================================
// INELASTIC SCATTERING
@ -650,9 +650,9 @@ void scatter(Particle* p, int i_nuclide)
// Perform collision physics for inelastic scattering
const auto& rx {nuc->reactions_[i]};
std::cout << "Energy before Inelastic scatter = " << p->E_ << std::endl;
//std::cout << "Energy before Inelastic scatter = " << p->E_ << std::endl;
inelastic_scatter(nuc.get(), rx.get(), p);
std::cout << "Energy after Inelastic scatter = " << p->E_ << std::endl;
//std::cout << "Energy after Inelastic scatter = " << p->E_ << std::endl;
p->event_mt_ = rx->mt_;
}
@ -1061,35 +1061,35 @@ void inelastic_scatter(const Nuclide* nuc, const Reaction* rx, Particle* p)
double mu;
// This gives a negative sometimes for some unknown reason
rx->products_[0].sample(E_in, E, mu);
std::cout << "E_in = " << E_in << "E = " << E << " mu = " << mu << std::endl;
//std::cout << "E_in = " << E_in << "E = " << E << " mu = " << mu << std::endl;
// if scattering system is in center-of-mass, transfer cosine of scattering
// angle and outgoing energy from CM to LAB
if (rx->scatter_in_cm_) {
double E_cm = E;
std::cout << "E_cm = " << E_cm << std::endl;
//std::cout << "E_cm = " << E_cm << std::endl;
// determine outgoing energy in lab
double A = nuc->awr_;
std::cout << "A = " << nuc->awr_ << std::endl;
// std::cout << "A = " << nuc->awr_ << std::endl;
E = E_cm + (E_in + 2.0*mu*(A + 1.0) * std::sqrt(E_in*E_cm))
/ ((A + 1.0)*(A + 1.0));
// here's where it goes wrong
std::cout << "E = " << E << std::endl;
//std::cout << "E = " << E << std::endl;
// determine outgoing angle in lab
mu = mu*std::sqrt(E_cm/E) + 1.0/(A+1.0) * std::sqrt(E_in/E);
std::cout << "mu = " << mu << std::endl;
// std::cout << "mu = " << mu << std::endl;
}
// Because of floating-point roundoff, it may be possible for mu to be
// outside of the range [-1,1). In these cases, we just set mu to exactly -1
// or 1
if (std::abs(mu) > 1.0) mu = std::copysign(1.0, mu);
std::cout << "mu = " << mu << std::endl;
// std::cout << "mu = " << mu << std::endl;
// Set outgoing energy and scattering angle
std::cout << "FINAL E = " << E << " mu = " << mu << std::endl;
// std::cout << "FINAL E = " << E << " mu = " << mu << std::endl;
p->E_ = E;
p->mu_ = mu;

35
src/simulation.cpp
View file

@ -371,7 +371,7 @@ void process_collision_events()
//for (auto& p : collision_queue) {
for (int i = 0; i < collision_queue_length; i++) {
Particle * p = particles + collision_queue[i];
std::cout << "Beginning collision of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
//std::cout << "Beginning collision of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
// Score collision estimate of keff
if (settings::run_mode == RUN_MODE_EIGENVALUE &&
p->type_ == Particle::Type::neutron) {
@ -386,7 +386,7 @@ void process_collision_events()
if (!model::active_meshsurf_tallies.empty())
score_surface_tally(p, model::active_meshsurf_tallies);
std::cout << "After surface tally of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
//std::cout << "After surface tally of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
// Clear surface component
p->surface_ = 0;
@ -397,7 +397,7 @@ void process_collision_events()
collision_mg(p);
}
std::cout << "After collision() of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
//std::cout << "After collision() of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
// Score collision estimator tallies -- this is done after a collision
// has occurred rather than before because we need information on the
@ -410,7 +410,7 @@ void process_collision_events()
score_analog_tally_mg(p);
}
}
std::cout << "After analog tally of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
//std::cout << "After analog tally of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
// Reset banked weight during collision
p->n_bank_ = 0;
@ -456,7 +456,7 @@ void process_collision_events()
if (p->alive_)
{
dispatch_xs_event(collision_queue[i]);
std::cout << "Ended collision of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
//std::cout << "Ended collision of particle id " << collision_queue[i] << " with energy E = " << p->E_ << std::endl;
assert(std::isfinite(p->E_) );
}
}
@ -561,42 +561,46 @@ void transport()
dispatch_xs_event(i);
}
int event_kernel_executions = 0;
while (true) {
event_kernel_executions++;
/*
std::cout << "Fuel XS Lookups = " << calculate_fuel_xs_queue_length << std::endl;
std::cout << "Non Fuel XS Lookups = " << calculate_nonfuel_xs_queue_length << std::endl;
std::cout << "Advance Particles = " << advance_particle_queue_length << std::endl;
std::cout << "Surface Crossings = " << surface_crossing_queue_length << std::endl;
std::cout << "Collisions = " << collision_queue_length << std::endl;
*/
int max = std::max({calculate_fuel_xs_queue_length, calculate_nonfuel_xs_queue_length, advance_particle_queue_length, surface_crossing_queue_length, collision_queue_length});
check_energies();
if (max == 0) {
break;
} else if (max == calculate_fuel_xs_queue_length) {
std::cout << "pre fuel XS check..." << std::endl;
//std::cout << "pre fuel XS check..." << std::endl;
//check_energies(calculate_fuel_xs_queue, calculate_fuel_xs_queue_length);
std::cout << "Performing Fuel XS Lookups..." << std::endl;
//std::cout << "Performing Fuel XS Lookups..." << std::endl;
process_calculate_xs_events(calculate_fuel_xs_queue, calculate_fuel_xs_queue_length);
calculate_fuel_xs_queue_length = 0;
} else if (max == calculate_nonfuel_xs_queue_length) {
std::cout << "pre non fuel XS check..." << std::endl;
//std::cout << "pre non fuel XS check..." << std::endl;
//check_energies(calculate_nonfuel_xs_queue, calculate_nonfuel_xs_queue_length);
std::cout << "Performing Non Fuel XS Lookups..." << std::endl;
// std::cout << "Performing Non Fuel XS Lookups..." << std::endl;
process_calculate_xs_events(calculate_nonfuel_xs_queue, calculate_nonfuel_xs_queue_length);
calculate_nonfuel_xs_queue_length = 0;
} else if (max == advance_particle_queue_length) {
std::cout << "pre advancing check..." << std::endl;
//std::cout << "pre advancing check..." << std::endl;
//check_energies(advance_particle_queue, advance_particle_queue_length);
std::cout << "Advancing Particles..." << std::endl;
//std::cout << "Advancing Particles..." << std::endl;
process_advance_particle_events();
} else if (max == surface_crossing_queue_length) {
std::cout << "pre surface crossing check..." << std::endl;
//std::cout << "pre surface crossing check..." << std::endl;
//check_energies(surface_crossing_queue, surface_crossing_queue_length);
std::cout << "Surface Crossings..." << std::endl;
//std::cout << "Surface Crossings..." << std::endl;
process_surface_crossing_events();
} else if (max == collision_queue_length) {
std::cout << "pre Colliding check..." << std::endl;
//std::cout << "pre Colliding check..." << std::endl;
//check_energies(collision_queue, collision_queue_length);
std::cout << "Colliding..." << std::endl;
//std::cout << "Colliding..." << std::endl;
process_collision_events();
}
}
@ -615,6 +619,7 @@ void transport()
// TODO: Do this only once per PI
free_event_queues();
std::cout << "Event kernels retired: " << event_kernel_executions << std::endl;
}
}