Skip to content

Commit 4edf615

Browse files
cinaralcinaral
committed
Update matrix_op and matrix_rw, change 2D array indexing
1 parent 4399f2f commit 4edf615

File tree

12 files changed

+81
-70
lines changed

12 files changed

+81
-70
lines changed

‎.gitattributes

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -1,3 +1,4 @@
11
* text=auto
22
*.m linguist-language=MATLAB
3+
*.sh linguist-detectable=false
34
CMakeLists.txt linguist-detectable=false

‎CMakeLists.txt

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -2,15 +2,15 @@ cmake_minimum_required(VERSION 3.13)
22

33
project(
44
rk4_solver
5-
VERSION 1.3.0
5+
VERSION 1.4.0
66
LANGUAGES CXX
77
)
88

99
#* dependencies
1010
include(FetchContent)
1111

12-
set(matrix_op_VERSION 1.0.7) #* matrix_op, used for basic matrix operations
13-
set(matrix_rw_VERSION 1.0.5) #* matrix_rw, used for testing
12+
set(matrix_op_VERSION 1.1.0) #* matrix_op, used for basic matrix operations
13+
set(matrix_rw_VERSION 1.1.0) #* matrix_rw, used for testing
1414

1515
FetchContent_Declare(
1616
matrix_op

‎README.md

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -96,7 +96,7 @@ loop(
9696
const Real_T t_init,
9797
const Real_T (&x0)[X_DIM],
9898
Real_T (&t)[T_DIM],
99-
Real_T (&x)[T_DIM * X_DIM]
99+
Real_T (&x)[T_DIM][X_DIM]
100100
);
101101
```
102102
# 4. Examples
@@ -139,9 +139,9 @@ int main()
139139
{
140140
rk4_solver::Integrator<x_dim, Dynamics> integrator(dynamics, &Dynamics::ode_fun, time_step);
141141
Real_T t_arr[t_dim];
142-
Real_T x_arr[t_dim * x_dim];
142+
Real_T x_arr[t_dim][x_dim];
143143
//* save the intermediate values
144-
rk4_solver::loop<t_dim>(integrator, t_init, x_init, t_arr, x_arr);
144+
rk4_solver::loop(integrator, t_init, x_init, t_arr, x_arr);
145145
//...
146146
}
147147
```

‎benchmarks/loop-benchmark.cpp

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -38,7 +38,7 @@ main()
3838
Real_T x_arr[t_dim * x_dim];
3939
#else
4040
Real_T(&t_arr)[t_dim] = *(Real_T(*)[t_dim]) new Real_T[t_dim];
41-
Real_T(&x_arr)[t_dim * x_dim] = *(Real_T(*)[t_dim * x_dim]) new Real_T[t_dim * x_dim];
41+
Real_T(&x_arr)[t_dim][x_dim] = *(Real_T(*)[t_dim][x_dim]) new Real_T[t_dim][x_dim];
4242
#endif
4343

4444
printf("Cumulatively integrating 3rd order linear ODE for %.3g steps... ",
@@ -53,7 +53,7 @@ main()
5353
const auto since_sample_ns =
5454
std::chrono::duration_cast<std::chrono::nanoseconds>(now_tp - start_tp);
5555

56-
const Real_T(&x_final)[x_dim] = *matrix_op::select_row<t_dim, x_dim>(t_dim - 1, x_arr);
56+
const Real_T(&x_final)[x_dim] = x_arr[t_dim - 1];
5757

5858
printf("Done.\nx at t = %.3g s: [%.3g; %.3g; %.3g]\n", t_arr[t_dim - 1], x_final[0],
5959
x_final[1], x_final[2]);

‎benchmarks/step-benchmark.cpp

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -34,7 +34,9 @@ main()
3434
{
3535
Real_T t = t_init;
3636
Real_T x[x_dim];
37-
matrix_op::replace_row<1>(0, x_init, x); //* initialize x
37+
for (size_t i = 0; i < x_dim; ++i) {
38+
x[i] = x_init[i];
39+
}
3840

3941
printf("Integrating 3rd order linear ODE for %zu ms... ", benchmark_duration_ms);
4042
auto sample_tp = std::chrono::high_resolution_clock::now();

‎examples/loop-example.cpp

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -29,9 +29,9 @@ main()
2929
{
3030
rk4_solver::Integrator<x_dim, Dynamics> integrator(dynamics, &Dynamics::ode_fun, time_step);
3131
Real_T t_arr[t_dim];
32-
Real_T x_arr[t_dim * x_dim];
32+
Real_T x_arr[t_dim][x_dim];
3333
//* save the intermediate values
34-
rk4_solver::loop<t_dim>(integrator, t_init, x_init, t_arr, x_arr);
35-
34+
rk4_solver::loop(integrator, t_init, x_init, t_arr, x_arr);
35+
3636
return 0;
3737
}

‎include/rk4_solver/loop.hpp

Lines changed: 17 additions & 9 deletions
Original file line numberDiff line numberDiff line change
@@ -50,8 +50,11 @@ void
5050
loop(Integrator<X_DIM, T> integrator, const Real_T &t_init, const Real_T (&x_init)[X_DIM],
5151
Real_T &t, Real_T (&x)[X_DIM])
5252
{
53-
t = t_init; //* initialize t
54-
matrix_op::replace_row<1>(0, x_init, x); //* initialize x
53+
t = t_init; //* initialize t
54+
55+
for (size_t i = 0; i < X_DIM; ++i) {
56+
x[i] = x_init[i]; //* initialize x
57+
}
5558

5659
for (size_t i = 0; i < T_DIM - 1; ++i) {
5760
integrator.step(t, x, t, x); //* update t, x to the next t, x
@@ -72,7 +75,7 @@ loop(Integrator<X_DIM, T> integrator, const Real_T &t_init, const Real_T (&x_ini
7275
template <size_t T_DIM, size_t X_DIM, typename T>
7376
void
7477
loop(Integrator<X_DIM, T> integrator, const Real_T &t_init, const Real_T (&x_init)[X_DIM],
75-
Real_T (&t_arr)[T_DIM], Real_T (&x_arr)[T_DIM * X_DIM])
78+
Real_T (&t_arr)[T_DIM], Real_T (&x_arr)[T_DIM][X_DIM])
7679
{
7780
t_arr[0] = t_init; //* initialize t
7881
matrix_op::replace_row<T_DIM>(0, x_init, x_arr); //* initialize x
@@ -82,7 +85,7 @@ loop(Integrator<X_DIM, T> integrator, const Real_T &t_init, const Real_T (&x_ini
8285

8386
for (size_t i = 0; i < T_DIM - 1; ++i) {
8487
const Real_T t = t_arr[i];
85-
const Real_T(&x)[X_DIM] = *matrix_op::select_row<T_DIM, X_DIM>(i, x_arr);
88+
const Real_T(&x)[X_DIM] = x_arr[i];
8689
integrator.step(t, x, t_next, x_next); //* update t, x to the next t, x
8790
t_arr[i + 1] = t_next;
8891
matrix_op::replace_row<T_DIM>(i + 1, x_next, x_arr);
@@ -107,13 +110,18 @@ size_t
107110
loop(Integrator<X_DIM, T> integrator, Event<X_DIM, T> event, const Real_T &t_init,
108111
const Real_T (&x_init)[X_DIM], Real_T &t, Real_T (&x)[X_DIM], bool halt_on_event = false)
109112
{
110-
t = t_init; //* initialize t
111-
matrix_op::replace_row<1>(0, x_init, x); //* initialize x
113+
t = t_init; //* initialize t
114+
115+
for (size_t i = 0; i < X_DIM; ++i) {
116+
x[i] = x_init[i]; //* initialize x
117+
}
112118
Real_T x_plus[X_DIM];
113119

114120
for (size_t i = 0; i < T_DIM - 1; ++i) {
115121
if (event.check(t, x, x_plus)) {
116-
matrix_op::replace_row<1>(0, x_plus, x);
122+
for (size_t j = 0; j < X_DIM; ++j) {
123+
x[j] = x_plus[j];
124+
}
117125

118126
if (halt_on_event) {
119127
break;
@@ -140,7 +148,7 @@ loop(Integrator<X_DIM, T> integrator, Event<X_DIM, T> event, const Real_T &t_ini
140148
template <size_t T_DIM, size_t X_DIM, typename T>
141149
size_t
142150
loop(Integrator<X_DIM, T> integrator, Event<X_DIM, T> event, const Real_T &t_init,
143-
const Real_T (&x_init)[X_DIM], Real_T (&t_arr)[T_DIM], Real_T (&x_arr)[T_DIM * X_DIM],
151+
const Real_T (&x_init)[X_DIM], Real_T (&t_arr)[T_DIM], Real_T (&x_arr)[T_DIM][X_DIM],
144152
bool halt_on_event = false)
145153
{
146154
t_arr[0] = t_init; //* initialize t
@@ -151,7 +159,7 @@ loop(Integrator<X_DIM, T> integrator, Event<X_DIM, T> event, const Real_T &t_ini
151159

152160
for (size_t i = 0; i < T_DIM - 1; ++i) {
153161
const Real_T t = t_arr[i];
154-
const Real_T(&x)[X_DIM] = *matrix_op::select_row<T_DIM, X_DIM>(i, x_arr);
162+
const Real_T(&x)[X_DIM] = x_arr[i];
155163

156164
if (event.check(t, x, x_plus)) {
157165
integrator.step(t, x_plus, t_next, x_next);

‎test/ball-test.cpp

Lines changed: 15 additions & 15 deletions
Original file line numberDiff line numberDiff line change
@@ -58,42 +58,42 @@ int
5858
main()
5959
{
6060
//* 1. read the reference data
61-
Real_T x_arr_ref[t_dim * x_dim];
62-
matrix_rw::read<t_dim, x_dim>(ref_dat_prefix + test_config::x_arr_ref_fname, x_arr_ref);
61+
Real_T x_arr_ref[t_dim][x_dim];
62+
matrix_rw::read(ref_dat_prefix + test_config::x_arr_ref_fname, x_arr_ref);
6363

6464
//* 2. test
6565
Real_T t = 0;
6666
Real_T x[x_dim];
67-
Real_T t_arr[t_dim];
68-
Real_T x_arr[t_dim * x_dim];
67+
Real_T t_arr[1][t_dim];
68+
Real_T x_arr[t_dim][x_dim];
6969

7070
rk4_solver::Integrator<x_dim, Dynamics> integrator(dynamics, &Dynamics::ode_fun, time_step);
7171
rk4_solver::Event<x_dim, Dynamics> event(dynamics, &Dynamics::event_fun);
7272

7373
rk4_solver::loop<t_dim>(integrator, event, t_init, x_init, t, x);
7474

7575
integrator.reset();
76-
rk4_solver::loop<t_dim>(integrator, event, t_init, x_init, t_arr, x_arr);
76+
rk4_solver::loop(integrator, event, t_init, x_init, t_arr[0], x_arr);
7777

7878
//* 3. write the test data
79-
matrix_rw::write<t_dim, 1>(dat_prefix + test_config::t_arr_fname, t_arr);
80-
matrix_rw::write<t_dim, x_dim>(dat_prefix + test_config::x_arr_fname, x_arr);
79+
matrix_rw::write(dat_prefix + test_config::t_arr_fname, t_arr);
80+
matrix_rw::write(dat_prefix + test_config::x_arr_fname, x_arr);
8181

8282
//* 4. verify the results
83-
Real_T x_0_arr[t_dim * 1];
84-
Real_T x_0_arr_ref[t_dim * 1];
83+
Real_T x_0_arr[1][t_dim];
84+
Real_T x_0_arr_ref[1][t_dim];
8585

8686
for (size_t i = 0; i < t_dim; ++i) {
87-
const Real_T(&x)[x_dim] = *matrix_op::select_row<t_dim, x_dim>(i, x_arr);
88-
const Real_T(&x_ref)[x_dim] = *matrix_op::select_row<t_dim, x_dim>(i, x_arr_ref);
89-
x_0_arr[i] = x[verify_idx];
90-
x_0_arr_ref[i] = x_ref[verify_idx];
87+
const Real_T(&x)[x_dim] = x_arr[i];
88+
const Real_T(&x_ref)[x_dim] = x_arr_ref[i];
89+
x_0_arr[0][i] = x[verify_idx];
90+
x_0_arr_ref[0][i] = x_ref[verify_idx];
9191
}
92-
Real_T max_error = test_config::compute_max_error<t_dim, 1>(x_0_arr, x_0_arr_ref);
92+
Real_T max_error = test_config::compute_max_error(x_0_arr, x_0_arr_ref);
9393

9494
//* loop vs cum_loop sanity check
9595
Real_T max_loop_error = 0.;
96-
const Real_T(&x_final)[x_dim] = *matrix_op::select_row<t_dim, x_dim>(t_dim - 1, x_arr);
96+
const Real_T(&x_final)[x_dim] = x_arr[t_dim - 1];
9797

9898
for (size_t i = 0; i < x_dim; ++i) {
9999
const Real_T error = std::abs(x_final[i] - x[i]);

‎test/first_order-test.cpp

Lines changed: 10 additions & 10 deletions
Original file line numberDiff line numberDiff line change
@@ -41,31 +41,31 @@ main()
4141
//* 2. test
4242
Real_T t = 0;
4343
Real_T x[x_dim];
44-
Real_T t_arr[t_dim];
45-
Real_T x_arr[t_dim * x_dim];
44+
Real_T t_arr[1][t_dim];
45+
Real_T x_arr[t_dim][x_dim];
4646

4747
rk4_solver::Integrator<x_dim, Dynamics> integrator(dynamics, &Dynamics::ode_fun, time_step);
4848
rk4_solver::loop<t_dim>(integrator, t_init, x_init, t, x);
4949

5050
integrator.reset();
51-
rk4_solver::loop<t_dim>(integrator, t_init, x_init, t_arr, x_arr);
51+
rk4_solver::loop(integrator, t_init, x_init, t_arr[0], x_arr);
5252

5353
//* 3. write the test data
54-
matrix_rw::write<t_dim, 1>(dat_prefix + test_config::t_arr_fname, t_arr);
55-
matrix_rw::write<t_dim, x_dim>(dat_prefix + test_config::x_arr_fname, x_arr);
54+
matrix_rw::write(dat_prefix + test_config::t_arr_fname, t_arr);
55+
matrix_rw::write(dat_prefix + test_config::x_arr_fname, x_arr);
5656

5757
//* 4. verify the results
58-
Real_T x_arr_ref[t_dim * x_dim];
58+
Real_T x_arr_ref[t_dim][x_dim];
5959

6060
for (size_t i = 0; i < t_dim; ++i) {
61-
Real_T x_ref[x_dim] = {std::exp(a_const * t_arr[i])};
62-
matrix_op::replace_row<t_dim, x_dim>(i, x_ref, x_arr_ref);
61+
Real_T x_ref[x_dim] = {std::exp(a_const * t_arr[0][i])};
62+
matrix_op::replace_row(i, x_ref, x_arr_ref);
6363
}
64-
Real_T max_error = test_config::compute_max_error<t_dim, x_dim>(x_arr, x_arr_ref);
64+
Real_T max_error = test_config::compute_max_error(x_arr, x_arr_ref);
6565

6666
//* cumulative loop sanity check
6767
Real_T max_loop_error = 0.;
68-
const Real_T(&x_final)[x_dim] = *matrix_op::select_row<t_dim, x_dim>(t_dim - 1, x_arr);
68+
const Real_T(&x_final)[x_dim] = x_arr[t_dim - 1];
6969

7070
for (size_t i = 0; i < x_dim; ++i) {
7171
const Real_T error = std::abs(x_final[i] - x[i]);

‎test/motor-test.cpp

Lines changed: 11 additions & 11 deletions
Original file line numberDiff line numberDiff line change
@@ -20,8 +20,8 @@ constexpr Real_T J = 1.29e-4; //* [kg m-2]
2020
constexpr Real_T b = 3.92e-4; //* [N m s]
2121
constexpr Real_T K_t = 6.4e-2; //* [N m A-1]
2222
constexpr Real_T K_b = 6.4e-2; //* [V s]
23-
constexpr Real_T A[x_dim * x_dim] = {0, 1, 0, 0, -b / J, K_t / J, 0, -K_b / L, -R / L};
24-
constexpr Real_T B[x_dim * u_dim] = {0, 0, 1 / L};
23+
constexpr Real_T A[x_dim][x_dim] = {{0, 1, 0}, {0, -b / J, K_t / J}, {0, -K_b / L, -R / L}};
24+
constexpr Real_T B[x_dim][u_dim] = {{0}, {0}, {1 / L}};
2525
constexpr Real_T since_ampl = 10; //* input amplitude
2626
constexpr Real_T sine_freq = 10; //* input frequency
2727

@@ -69,30 +69,30 @@ int
6969
main()
7070
{
7171
//* 1. read the reference data
72-
Real_T x_arr_ref[t_dim * x_dim];
73-
matrix_rw::read<t_dim, x_dim>(ref_dat_prefix + test_config::x_arr_ref_fname, x_arr_ref);
72+
Real_T x_arr_ref[t_dim][x_dim];
73+
matrix_rw::read(ref_dat_prefix + test_config::x_arr_ref_fname, x_arr_ref);
7474

7575
//* 2. test
7676
Real_T t = 0;
7777
Real_T x[x_dim];
78-
Real_T t_arr[t_dim];
79-
Real_T x_arr[t_dim * x_dim];
78+
Real_T t_arr[1][t_dim];
79+
Real_T x_arr[t_dim][x_dim];
8080
rk4_solver::Integrator<x_dim, Dynamics> integrator(dynamics, &Dynamics::ode_fun, time_step);
8181

8282
rk4_solver::loop<t_dim>(integrator, t_init, x_init, t, x);
8383
integrator.reset();
84-
rk4_solver::loop<t_dim>(integrator, t_init, x_init, t_arr, x_arr);
84+
rk4_solver::loop(integrator, t_init, x_init, t_arr[0], x_arr);
8585

8686
//* 3. write the test data
87-
matrix_rw::write<t_dim, 1>(dat_prefix + test_config::t_arr_fname, t_arr);
88-
matrix_rw::write<t_dim, x_dim>(dat_prefix + test_config::x_arr_fname, x_arr);
87+
matrix_rw::write(dat_prefix + test_config::t_arr_fname, t_arr);
88+
matrix_rw::write(dat_prefix + test_config::x_arr_fname, x_arr);
8989

9090
//* 4. verify the results
91-
Real_T max_error = test_config::compute_max_error<t_dim, x_dim>(x_arr, x_arr_ref);
91+
Real_T max_error = test_config::compute_max_error(x_arr, x_arr_ref);
9292

9393
//* loop vs cum_loop sanity check
9494
Real_T max_loop_error = 0.;
95-
const Real_T(&x_final)[x_dim] = *matrix_op::select_row<t_dim, x_dim>(t_dim - 1, x_arr);
95+
const Real_T(&x_final)[x_dim] = x_arr[t_dim - 1];
9696

9797
for (size_t i = 0; i < x_dim; ++i) {
9898
const Real_T error = std::abs(x_final[i] - x[i]);

0 commit comments

Comments
 (0)