layers/craters/ages/steady_state.cpp

Line	Branch	Exec	Source
1			#include "terraingen/layers/craters/ages/steady_state.h"
2
3			#include <cmath>
4			#include <stdexcept>
5			#include <utility>
6
7			namespace terraingen {
8
9			namespace {
10
11			/// An age in [0, `maxAgeGyr`] drawn from the formation rate of `chronology`, by bisection on
12			/// N1(T).
13		782009	double SampleFormationAge(RandomStream& random, const Chronology& chronology, double maxAgeGyr) {
14		782009	const double target = random.NextDouble() * chronology(maxAgeGyr);
15		782009	double low = 0.0;
16		782009	double high = maxAgeGyr;
17	3/4 ✓ Branch 0 taken 39887935 times. ✗ Branch 1 not taken. ✓ Branch 2 taken 39105926 times. ✓ Branch 3 taken 782009 times.	39887935	for (int i = 0; i < 100 && high - low > 1e-15 * maxAgeGyr; ++i) {
18		39105926	const double middle = 0.5 * (low + high);
19	2/2 ✓ Branch 1 taken 19612450 times. ✓ Branch 2 taken 19493476 times.	39105926	if (chronology(middle) < target) {
20		19612450	low = middle;
21			} else {
22		19493476	high = middle;
23			}
24			}
25		782009	return 0.5 * (low + high);
26			}
27
28			} // namespace
29
30		427	SteadyStateAges::SteadyStateAges(std::shared_ptr<const SizeFrequencyDistribution> productionRate,
31			std::shared_ptr<const SizeFrequencyDistribution> equilibrium,
32		427	double ageGyr, Chronology chronology)
33		427	: productionRate_(std::move(productionRate)), equilibrium_(std::move(equilibrium)),
34		854	ageGyr_(ageGyr), chronology_(std::move(chronology)) {
35	6/6 ✓ Branch 1 taken 426 times. ✓ Branch 2 taken 1 time. ✓ Branch 5 taken 1 time. ✓ Branch 6 taken 425 times. ✓ Branch 7 taken 2 times. ✓ Branch 8 taken 425 times.	427	if (!productionRate_ \|\| productionRate_->Kind() != SfdKind::ProductionRate) {
36	1/1 ✓ Branch 2 taken 2 times.	2	throw std::invalid_argument("Steady-state ages need a production rate distribution");
37			}
38	6/6 ✓ Branch 1 taken 424 times. ✓ Branch 2 taken 1 time. ✓ Branch 5 taken 1 time. ✓ Branch 6 taken 423 times. ✓ Branch 7 taken 2 times. ✓ Branch 8 taken 423 times.	425	if (!equilibrium_ \|\| equilibrium_->Kind() != SfdKind::Density) {
39	1/1 ✓ Branch 2 taken 2 times.	2	throw std::invalid_argument("Steady-state ages need an equilibrium density");
40			}
41	2/2 ✓ Branch 1 taken 1 time. ✓ Branch 2 taken 422 times.	423	if (!chronology_) {
42	1/1 ✓ Branch 2 taken 1 time.	1	throw std::invalid_argument("Steady-state ages need a chronology");
43			}
44	1/1 ✓ Branch 1 taken 418 times.	422	chronology_(ageGyr_);
45		454	}
46
47		1017188	double SteadyStateAges::EquilibriumLifetimeGyr(double diameter) const {
48		1017188	return equilibrium_->Incremental(diameter) / productionRate_->Incremental(diameter);
49			}
50
51		1017183	double SteadyStateAges::LifetimeAt(double diameter) const {
52		1017183	return std::fmin(EquilibriumLifetimeGyr(diameter), ageGyr_);
53			}
54
55		782010	double SteadyStateAges::SampleAgeAt(RandomStream& random, double /diameter/,
56			double lifetimeGyr) const {
57	2/2 ✓ Branch 0 taken 1 time. ✓ Branch 1 taken 782009 times.	782010	if (lifetimeGyr == 0.0) {
58		1	return 0.0;
59			}
60		782009	return SampleFormationAge(random, chronology_, lifetimeGyr);
61			}
62
63			} // namespace terraingen
64

TerrainGen coverage