Figure 2.2:

Estimated reaction rates from six production-rate measurements subject to measurement noise.

Code for Figure 2.2

Text of the GNU GPL.

main.py

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
# Converted from rate_estimation.m
import numpy as np
import os

stoi = np.array([
    [0, 1, 0, -1, -1, 1],
    [-1, 1, 1, -1, 0, 0],
    [1, 0, -1, 0, -1, 1]
], dtype=float)

r = np.array([1., 2.])
A = stoi[:2, :]
R = A.T @ r

# Set seed for reproducible random numbers
rng = np.random.RandomState(1)

npoints = 6
Rmeas = np.zeros((6, npoints))
for i in range(npoints):
    Rmeas[:, i] = 0.05 * rng.randn(6) + R

rest = (np.linalg.inv(A @ A.T) @ A @ Rmeas).T

with open('rate_estimation.dat', 'w') as myfile:
    for row in rest:
        myfile.write(f'{row[0]:8.4f} {row[1]:8.4f} \n')