Started on Random Library

I made seed and uniform functions, but they are not omp friendly. They are runnning omp themself, but is not safe for threading/parallizing. Uniform is aproximated and is NOT validaded up against a real uniform distribution.
This commit is contained in:
2026-07-28 17:32:49 +02:00
parent 52684e6b8a
commit 441540a996
7 changed files with 1359 additions and 15 deletions
+110
View File
@@ -0,0 +1,110 @@
/**++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
*
* PANIC
* Portable Algorithms and Numerics In C++
*
* Scientific computing from scratch, with feeling.
*
* Copyright (c) 2026 Michelle Bausager
*
* This file is part of PANIC.
*
* PANIC is free software licensed under the GNU General Public License v3.0 or later.
* You may redistribute and/or modify it under the terms of the GPL.
*
* PANIC is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY;
* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
* See the LICENSE file for the full license text.
*
* SPDX-License-Identifier: GPL-3.0-or-later
*
*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
*
* Project Name: PANIC
* Module Name: random
* File Name: seed.hpp
* Revision: 0.1.0
* Date: 28-07-2026
* Author: Michelle Bausager
*
* Description:
* Defines seed for use in other functions in random/
*
*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
#pragma once
//---------------------------------------------------------------------------------------------------------------------------
// INCLUDE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
#include <config/types.hpp> // panic::uint_t, panic::int_t, and panic::real_t
//--------------------------------------------------------------------------------------------------------------------------
// Function Name : panic::random::seed
//
// Description:
// base for random libary
//--------------------------------------------------------------------------------------------------------------------------
namespace panic{
namespace random{
/**
* @brief struct for seed object used in random/
*
* The struct is used for PANIC random library.
*/
struct seed_t{
// Variable to store the value the random value is based on.
panic::types::uint_t value;
/**
* @brief Empthy constructor
*
*/
seed_t();
/**
* @brief Contructor with seed.
*
* @param seed The seed that is used as the base for the object.
*/
seed_t(panic::types::uint_t seed);
/**
* @brief Function to set the seed..
*
* @param seed The seed that is used as the base for the object.
*/
bool set(panic::types::uint_t seed);
// Creates a deterministic state from the seed and an index.
// This is OMP-friendly because it does not modify shared memory.
/**
* @brief Returns a random number based on seed and index
*
* @param index Value the random number will be drawn from.
*
* @Note this can be used with omp in a loop if the index in the loop
* is used as the input to this function.
*/
panic::types::uint_t state_at(panic::types::uint_t index) const;
/**
* @brief Returns the base seed value.
*
*/
panic::types::uint_t get();
};
} // namespace random
} // namespace panic
//---------------------------------------------------------------------------------------------------------------------------
// VARIABLE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
//---------------------------------------------------------------------------------------------------------------------------
// FUNCTION PROTOTYPE
//---------------------------------------------------------------------------------------------------------------------------
+196
View File
@@ -0,0 +1,196 @@
/**++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
*
* PANIC
* Portable Algorithms and Numerics In C++
*
* Scientific computing from scratch, with feeling.
*
* Copyright (c) 2026 Michelle Bausager
*
* This file is part of PANIC.
*
* PANIC is free software licensed under the GNU General Public License v3.0 or later.
* You may redistribute and/or modify it under the terms of the GPL.
*
* PANIC is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY;
* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
* See the LICENSE file for the full license text.
*
* SPDX-License-Identifier: GPL-3.0-or-later
*
*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
*
* Project Name: PANIC
* Module Name: random
* File Name: uniform.hpp
* Revision: 0.1.0
* Date: 29-06-2026
* Author: Michelle Bausager
*
* Description:
* Defindes the functions that returns values based on a uniform distribution.
*
*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
#pragma once
//---------------------------------------------------------------------------------------------------------------------------
// INCLUDE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
#include <config/types.hpp> // panic::uint_t, panic::int_t, and panic::real_t
#include <tensor/vector.hpp>
#include <tensor/matrix.hpp>
//---------------------------------------------------------------------------------------------------------------------------
// DEFINE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
//---------------------------------------------------------------------------------------------------------------------------
// TYPE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
//--------------------------------------------------------------------------------------------------------------------------
// Function Name : panic::random::seed
//
// Description:
// base for random libary
//--------------------------------------------------------------------------------------------------------------------------
namespace panic{
namespace random{
/**
* @brief Returns a value from a uniform distribution bewteen 0 and 1.
*
* Computes:
* @code
* c = uniform();
* @endcode
*
* @return a value from a uniform distribution between 0 and 1
*
* @note This is not a omp-safe function
*/
panic::types::real_t uniform();
/**
* @brief Returns a value from a uniform distribution with limits.
*
* Computes:
* @code
* c = uniform(0.1f, 42.0f);
* @endcode
*
* @tparam T Numeric element type.
* @param min minimum limit for return value.
* @param max maximum limit for return value.
*
* @return a value from a uniform distribution with limits
*
* @note This is not a omp-safe function
*/
template <typename T>
T uniform(const T min, const T max);
/**
* @brief Filleds a vector with values from a uniform distribution from 0 to 1.
*
* Computes:
* @code
* panic::tensor::vector a(5);
* uniform(a);
* @endcode
*
* @tparam T Numeric element type.
* @param a Output vector filled iwht new values
*
* @return True of vector is filled correctly
*
* @note This is not a omp-safe function
*/
bool uniform(panic::tensor::real_vector& a);
/**
* @brief Filleds a vector with values from a uniform distribution from min to max.
*
* Computes:
* @code
* panic::tensor::real_vector a(5);
* uniform(a, 1.0f, 42.0f);
* @endcode
*
* @tparam T Numeric element type.
* @param a Output vector filled iwht new values
* @param min minimum limit for return value.
* @param max maximum limit for return value.
*
* @return True of vector is filled correctly
*
* @note This is not a omp-safe function
*/
template <typename T>
bool uniform(panic::tensor::vector<T>& a, const T min, const T max);
/**
* @brief Filleds a matrix with values from a uniform distribution from 0 to 1.
*
* Computes:
* @code
* panic::tensor::matrix A(5);
* uniform(A);
* @endcode
*
* @tparam T Numeric element type.
* @param a Output matrix filled with new values
*
* @return True of matrix is filled correctly
*
* @note This is not a omp-safe function
*/
bool uniform(panic::tensor::real_matrix& A);
/**
* @brief Filleds a matrix with values from a uniform distribution from min to max.
*
* Computes:
* @code
* panic::tensor::real_matrix A(5,5);
* uniform(A, 1.0f, 42.0f);
* @endcode
*
* @tparam T Numeric element type.
* @param a Output matrix filled with new values
* @param min minimum limit for return value.
* @param max maximum limit for return value.
*
* @return True of matrix is filled correctly
*
* @note This is not a omp-safe function
*/
template <typename T>
bool uniform(panic::tensor::matrix<T>& A, const T min, const T max);
/*
panic::tensor::real_vector uniform_vector(panic::tensor::real_vector& a);
panic::tensor::real_vector uniform_vector(const panic::types::real_t min, const panic::types::real_t max);
// Create and return vector
// Create and return matrix
*/
} // namespace random
} // namespace panic
//---------------------------------------------------------------------------------------------------------------------------
// VARIABLE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
//---------------------------------------------------------------------------------------------------------------------------
// FUNCTION PROTOTYPE
//---------------------------------------------------------------------------------------------------------------------------
+649
View File
@@ -34,6 +34,7 @@
//---------------------------------------------------------------------------------------------------------------------------
// INCLUDE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
#include <config/omp.hpp>
#include <iostream> // std::cout, std::endl
#include <config/types.hpp> // include types to use
#include <config/constants.hpp> // Math constants
@@ -43,6 +44,21 @@
#include <math/matmul.hpp>
#include <neural_network/layer/layer_dense.hpp>
#include <math/add.hpp>
#include <random/uniform.hpp>
// For omp tesing:
#include <chrono>
#include <cstddef>
#include <cstdint>
#include <iomanip>
#include <iostream>
#include <limits>
#if PANIC_HAS_OPENMP
#include <omp.h>
#endif
//---------------------------------------------------------------------------------------------------------------------------
// DEFINE DESCRIPTION
@@ -51,6 +67,601 @@
// #define TEST_FALG 1
//---------------------------------------------------------------------------------------------------------------------------
// Function Name : benchmark_omp_min_work
//
// Description:
// Measures an operation with one thread and with all available threads.
//
// It then tests several possible omp_min_work values and estimates which
// threshold best matches the measured results.
//---------------------------------------------------------------------------------------------------------------------------
void benchmark_omp_min_work(){
#if !PANIC_HAS_OPENMP
std::cout << "OpenMP is not enabled.\n";
#else
//-----------------------------------------------------------------------------------------------------------------------
// OPENMP SETTINGS
//-----------------------------------------------------------------------------------------------------------------------
omp_set_dynamic(0);
const int parallel_thread_count = omp_get_max_threads();
//-----------------------------------------------------------------------------------------------------------------------
// EDIT 1:
// Sizes to test.
//
// This list is suitable for vector operations.
//-----------------------------------------------------------------------------------------------------------------------
const panic::types::uint_t sizes[] = {
100,
250,
500,
1000,
2000,
4000,
8000,
12000,
16000,
24000,
32000,
50000,
75000,
100000,
150000,
250000,
500000,
1000000,
2000000,
5000000
};
/*
// Example sizes for square matrix multiplication:
const panic::types::uint_t sizes[] = {
4,
6,
8,
10,
12,
16,
20,
24,
32,
48,
64,
96,
128,
192,
256,
384,
512
};
*/
const std::size_t number_of_sizes =
sizeof(sizes) / sizeof(sizes[0]);
//-----------------------------------------------------------------------------------------------------------------------
// EDIT 2:
// Possible omp_min_work values to test.
//
// These are work values, not necessarily vector or matrix sizes.
//-----------------------------------------------------------------------------------------------------------------------
const std::uint64_t omp_min_work_values[] = {
0,
50,
100,
150,
200,
250,
300,
350,
400,
450,
500,
600,
750,
1000,
1500,
2000,
4000,
8000,
16000,
32000
};
/*
// Example omp_min_work values for matrix multiplication:
const std::uint64_t omp_min_work_values[] = {
0,
1000,
5000,
10000,
25000,
50000,
100000,
250000,
500000,
1000000,
2000000,
5000000,
10000000,
25000000,
50000000,
100000000
};
*/
const std::size_t number_of_omp_min_work_values =
sizeof(omp_min_work_values) /
sizeof(omp_min_work_values[0]);
//-----------------------------------------------------------------------------------------------------------------------
// BENCHMARK SETTINGS
//-----------------------------------------------------------------------------------------------------------------------
// The benchmark tries to perform approximately this much work for each size.
const std::uint64_t target_total_work = 50000000;
// Minimum and maximum number of repetitions for each size.
const std::uint64_t minimum_repetitions = 3;
const std::uint64_t maximum_repetitions = 1000;
// Prevent the compiler from treating all calculated results as unused.
volatile panic::types::real_t benchmark_sink = 0.0f;
//-----------------------------------------------------------------------------------------------------------------------
// ARRAYS FOR THE MEASURED RESULTS
//-----------------------------------------------------------------------------------------------------------------------
std::uint64_t measured_work[number_of_sizes];
double one_thread_results[number_of_sizes];
double parallel_results[number_of_sizes];
//-----------------------------------------------------------------------------------------------------------------------
// PRINT BENCHMARK INFORMATION
//-----------------------------------------------------------------------------------------------------------------------
std::cout
<< std::right
<< std::setw(12) << "Size"
<< std::setw(16) << "Work"
<< std::setw(14) << "Repetitions"
<< std::setw(18) << "1 thread (us)"
<< std::setw(18) << "Parallel (us)"
<< std::setw(12) << "Speedup"
<< std::setw(14) << "Fastest"
<< "\n";
std::cout
<< std::string(104, '-')
<< "\n";
//-----------------------------------------------------------------------------------------------------------------------
// TEST EVERY SIZE
//-----------------------------------------------------------------------------------------------------------------------
for (std::size_t size_index = 0;
size_index < number_of_sizes;
++size_index){
const panic::types::uint_t size =
sizes[size_index];
//-------------------------------------------------------------------------------------------------------------------
// EDIT 3:
// Initialize the vectors or matrices used by the operation.
//-------------------------------------------------------------------------------------------------------------------
// Vector-add example:
panic::tensor::real_vector a(size, 1.0f);
panic::tensor::real_vector b(size, 2.0f);
panic::tensor::real_vector c(size);
/*
// Square-matrix multiplication example:
panic::tensor::real_matrix a(
size,
size,
0.01f
);
panic::tensor::real_matrix b(
size,
size,
0.02f
);
panic::tensor::real_matrix c(
size,
size
);
*/
//-------------------------------------------------------------------------------------------------------------------
// EDIT 4:
// Calculate work in the same way as the function being tested.
//-------------------------------------------------------------------------------------------------------------------
// Vector operation:
const std::uint64_t work =
static_cast<std::uint64_t>(size);
/*
// Square matrix add:
const std::uint64_t work =
static_cast<std::uint64_t>(size) *
static_cast<std::uint64_t>(size);
*/
/*
// Square matrix multiplication:
const std::uint64_t work =
static_cast<std::uint64_t>(size) *
static_cast<std::uint64_t>(size) *
static_cast<std::uint64_t>(size);
*/
//-------------------------------------------------------------------------------------------------------------------
// CALCULATE NUMBER OF REPETITIONS
//-------------------------------------------------------------------------------------------------------------------
std::uint64_t repetitions =
target_total_work / work;
if (repetitions < minimum_repetitions){
repetitions = minimum_repetitions;
}
if (repetitions > maximum_repetitions){
repetitions = maximum_repetitions;
}
//-------------------------------------------------------------------------------------------------------------------
// MEASURE WITH ONE THREAD AND ALL THREADS
//
// test == 0: one OpenMP thread
// test == 1: all available OpenMP threads
//-------------------------------------------------------------------------------------------------------------------
double measured_time_us[2] = {
0.0,
0.0
};
for (int test = 0; test < 2; ++test){
if (test == 0){
omp_set_num_threads(1);
}
else{
omp_set_num_threads(
parallel_thread_count
);
}
//----------------------------------------------------------------------------------------------------------------
// WARM-UP
//----------------------------------------------------------------------------------------------------------------
for (std::size_t warmup = 0;
warmup < 2;
++warmup){
//------------------------------------------------------------------------------------------------------------
// EDIT 5:
// Put the operation being tested here.
//------------------------------------------------------------------------------------------------------------
panic::math::add(a, b, c);
// Matrix multiplication:
// panic::math::matmul(a, b, c);
}
//----------------------------------------------------------------------------------------------------------------
// TIMED LOOP
//----------------------------------------------------------------------------------------------------------------
const std::chrono::steady_clock::time_point start =
std::chrono::steady_clock::now();
for (std::uint64_t repetition = 0;
repetition < repetitions;
++repetition){
//------------------------------------------------------------------------------------------------------------
// EDIT 6:
// Put the same operation here.
//------------------------------------------------------------------------------------------------------------
panic::math::add(a, b, c);
// Matrix multiplication:
// panic::math::matmul(a, b, c);
}
const std::chrono::steady_clock::time_point end =
std::chrono::steady_clock::now();
//----------------------------------------------------------------------------------------------------------------
// CALCULATE AVERAGE TIME PER OPERATION
//----------------------------------------------------------------------------------------------------------------
const double total_time_us =
std::chrono::duration<double, std::micro>(
end - start
).count();
measured_time_us[test] =
total_time_us /
static_cast<double>(repetitions);
//----------------------------------------------------------------------------------------------------------------
// READ ONE RESULT
//
// Change this if the output cannot be accessed with c[size / 2].
//----------------------------------------------------------------------------------------------------------------
benchmark_sink += c[size / 2];
/*
// For a matrix:
benchmark_sink += c(
size / 2,
size / 2
);
*/
}
//-------------------------------------------------------------------------------------------------------------------
// STORE RESULTS
//-------------------------------------------------------------------------------------------------------------------
const double one_thread_us =
measured_time_us[0];
const double parallel_us =
measured_time_us[1];
measured_work[size_index] = work;
one_thread_results[size_index] = one_thread_us;
parallel_results[size_index] = parallel_us;
//-------------------------------------------------------------------------------------------------------------------
// PRINT RESULTS FOR THIS SIZE
//-------------------------------------------------------------------------------------------------------------------
const double speedup =
one_thread_us / parallel_us;
const char* fastest;
if (parallel_us < one_thread_us){
fastest = "parallel";
}
else{
fastest = "one_thread";
}
std::cout
<< std::right
<< std::setw(12) << size
<< std::setw(16) << work
<< std::setw(14) << repetitions
<< std::setw(18) << std::fixed << std::setprecision(3)
<< one_thread_us
<< std::setw(18) << parallel_us
<< std::setw(12) << speedup
<< std::setw(14) << fastest
<< "\n";
}
//-----------------------------------------------------------------------------------------------------------------------
// TEST THE POSSIBLE OMP_MIN_WORK VALUES
//-----------------------------------------------------------------------------------------------------------------------
std::uint64_t best_omp_min_work = 0;
double best_score =
std::numeric_limits<double>::max();
std::cout
<< "\n"
<< "omp_min_work"
<< ", average_slowdown"
<< "\n";
for (std::size_t threshold_index = 0;
threshold_index < number_of_omp_min_work_values;
++threshold_index){
const std::uint64_t omp_min_work =
omp_min_work_values[threshold_index];
double score = 0.0;
//-------------------------------------------------------------------------------------------------------------------
// SEE WHICH VERSION THIS THRESHOLD WOULD SELECT
//-------------------------------------------------------------------------------------------------------------------
for (std::size_t size_index = 0;
size_index < number_of_sizes;
++size_index){
double selected_time_us;
//----------------------------------------------------------------------------------------------------------------
// This uses > because the PANIC functions currently use:
//
// work > omp_min_work
//----------------------------------------------------------------------------------------------------------------
if (measured_work[size_index] > omp_min_work){
// This threshold would select OpenMP.
selected_time_us =
parallel_results[size_index];
}
else{
// This threshold would select serial execution.
selected_time_us =
one_thread_results[size_index];
}
//----------------------------------------------------------------------------------------------------------------
// FIND THE FASTEST MEASURED VERSION FOR THIS SIZE
//----------------------------------------------------------------------------------------------------------------
double fastest_time_us =
one_thread_results[size_index];
if (
parallel_results[size_index] <
fastest_time_us
){
fastest_time_us =
parallel_results[size_index];
}
//----------------------------------------------------------------------------------------------------------------
// CALCULATE HOW MUCH SLOWER THE SELECTED VERSION IS
//
// 1.0 means the threshold selected the fastest version.
// 1.1 means it was 10% slower than the fastest version.
//----------------------------------------------------------------------------------------------------------------
score +=
selected_time_us /
fastest_time_us;
}
//-------------------------------------------------------------------------------------------------------------------
// CALCULATE THE AVERAGE SCORE
//-------------------------------------------------------------------------------------------------------------------
score /=
static_cast<double>(number_of_sizes);
//-------------------------------------------------------------------------------------------------------------------
// PRINT THIS OMP_MIN_WORK RESULT
//-------------------------------------------------------------------------------------------------------------------
std::cout
<< omp_min_work
<< ", "
<< std::fixed
<< std::setprecision(4)
<< score
<< "\n";
//-------------------------------------------------------------------------------------------------------------------
// SAVE THE BEST OMP_MIN_WORK
//-------------------------------------------------------------------------------------------------------------------
if (score < best_score){
best_score = score;
best_omp_min_work = omp_min_work;
}
}
//-----------------------------------------------------------------------------------------------------------------------
// PRINT FINAL ESTIMATE
//-----------------------------------------------------------------------------------------------------------------------
std::cout
<< "\n"
<< "Best estimated omp_min_work: "
<< best_omp_min_work
<< "\n";
std::cout
<< "Average slowdown score: "
<< std::fixed
<< std::setprecision(4)
<< best_score
<< "\n";
std::cout
<< "\n"
<< "A score of 1.0000 means the threshold selected the\n"
<< "fastest measured version for every tested size.\n";
//-----------------------------------------------------------------------------------------------------------------------
// PRINT CHECKSUM
//-----------------------------------------------------------------------------------------------------------------------
std::cout
<< "Benchmark checksum: "
<< benchmark_sink
<< "\n";
#endif
}
//---------------------------------------------------------------------------------------------------------------------------
// VARIABLE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
@@ -73,6 +684,10 @@ panic::tensor::int_matrix C(2,2, 3);
int main(void) {
//benchmark_omp_min_work();
// Comment out benchmark_omp_min_work() when it is not needed.
a[2] = 1.2;
b[2] = 3;
@@ -107,6 +722,40 @@ int main(void) {
panic::neural_network::layer_dense dense(100,10);
std::cout << dense.forward(B1) << std::endl;
std::cout << "random" << std::endl;
for (int i = 0; i < 5; ++i)
{
std::cout << panic::random::uniform() << std::endl;
}
std::cout << "random(min, max)" << std::endl;
panic::types::real_t a1 = -15;
panic::types::real_t a2 = 15;
for (int i = 0; i < 5; ++i)
{
std::cout << panic::random::uniform(a1,a2) << std::endl;
}
panic::tensor::real_vector a3(2);
panic::random::uniform(a3);
panic::io::print_vector(a3);
panic::tensor::uint_vector a4(2);
panic::random::uniform(a4, static_cast<panic::types::uint_t>(1), static_cast<panic::types::uint_t>(3));
panic::io::print_vector(a4);
panic::tensor::real_matrix D1(3,3);
panic::random::uniform(D1);
panic::io::print_matrix(D1);
panic::tensor::real_matrix D3(3,3);
panic::random::uniform(D3, 100.f, 200.f);
panic::io::print_matrix(D3);
return 0;
+13 -6
View File
@@ -99,19 +99,26 @@ template bool print_vector<panic::types::real_t>(const panic::tensor::vector<pan
//--------------------------------------------------------------------------------------------------------------------------
template <typename T>
bool print_matrix(const panic::tensor::matrix<T>& A){
std::cout << "[";
for (panic::types::uint_t i = 0; i < A.rows(); ++i){
std::cout << "[";
for (panic::types::uint_t j = 0; j < A.cols(); ++j){
std::cout << A(i,j) << ", ";
}
if (i < A.rows()-1){
std::cout << A(A.rows()-1,A.cols()-1) << "]" << std::endl;
}else{
std::cout << A(A.rows()-1,A.cols()-1) << "]]" << std::endl;
std::cout << A(i, j);
if (j + 1 < A.cols()){
std::cout << ", ";
}
}
std::cout << "]";
if (i + 1 < A.rows())
std::cout << "," << std::endl;
}
std::cout << "]" << std::endl;
return true;
}
+2 -2
View File
@@ -51,8 +51,8 @@
* Small vectors and matrices are kept serial because the overhead of starting
* worker threads can be larger than the work itself.
*/
static const panic::types::uint_t add_omp_min_work = 10000;
//static const panic::types::uint_t add_omp_min_work = 10000;
static const panic::types::uint_t add_omp_min_work = 0;
//---------------------------------------------------------------------------------------------------------------------------
// INPLEMENTATION
//---------------------------------------------------------------------------------------------------------------------------
+122
View File
@@ -0,0 +1,122 @@
/**++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
*
* PANIC
* Portable Algorithms and Numerics In C++
*
* Scientific computing from scratch, with feeling.
*
* Copyright (c) 2026 Michelle Bausager
*
* This file is part of PANIC.
*
* PANIC is free software licensed under the GNU General Public License v3.0 or later.
* You may redistribute and/or modify it under the terms of the GPL.
*
* PANIC is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY;
* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
* See the LICENSE file for the full license text.
*
* SPDX-License-Identifier: GPL-3.0-or-later
*
*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
*
* Project Name: PANIC
* Module Name: random
* File Name: seed.cpp
* Revision: 0.1.0
* Date: 28-06-2026
* Author: Michelle Bausager
*
* Description:
* Defines seed for use in other functions in random/
*
*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
//---------------------------------------------------------------------------------------------------------------------------
// INCLUDE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
#include <random/seed.hpp>
#include <config/types.hpp>
//---------------------------------------------------------------------------------------------------------------------------
// FUNCTION PROTOTYPE
//---------------------------------------------------------------------------------------------------------------------------
namespace panic{
namespace random{
//--------------------------------------------------------------------------------------------------------------------------
// Constructor Name : panic::random::seed::seed_t
//
// Description:
// Creates an seed object.
//--------------------------------------------------------------------------------------------------------------------------
seed_t::seed_t(){
value = 1;
}
//--------------------------------------------------------------------------------------------------------------------------
// Constructor Name : panic::random::seed::seed_t
//
// Description:
// Creates an seed object with a specific seed.
//--------------------------------------------------------------------------------------------------------------------------
seed_t::seed_t(panic::types::uint_t seed){
set(seed);
}
//--------------------------------------------------------------------------------------------------------------------------
// Constructor Name : panic::random::seed::set_seed
//
// Description:
// Sets the seed to a none-zero value
//--------------------------------------------------------------------------------------------------------------------------
bool seed_t::set(panic::types::uint_t seed){
if (seed == 0){
value = 1;
}else{
value = seed;
}
return 1;
}
//--------------------------------------------------------------------------------------------------------------------------
// Constructor Name : panic::random::seed::get
//
// Description:
// Returns the seed.
//--------------------------------------------------------------------------------------------------------------------------
panic::types::uint_t seed_t::get(){
return value;
}
//--------------------------------------------------------------------------------------------------------------------------
// Constructor Name : panic::random::seed::state_at
//
// Description:
// Returns a pshodo-random number based on the input index and the seed value.
//--------------------------------------------------------------------------------------------------------------------------
panic::types::uint_t seed_t::state_at(panic::types::uint_t index) const{
panic::types::uint_t state = value + index;
// Mix it a few times so close indexes do not start too similarly.
//next_state(state);
//next_state(state);
//next_state(state);
state ^= state >> 16;
state *= static_cast<panic::types::uint_t>(0x7feb352d);
state ^= state >> 15;
state *= static_cast<panic::types::uint_t>(0x846ca68b);
state ^= state >> 16;
return state;
}
} // namespace tensor
} // namespace panic
+260
View File
@@ -0,0 +1,260 @@
/**++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
*
* PANIC
* Portable Algorithms and Numerics In C++
*
* Scientific computing from scratch, with feeling.
*
* Copyright (c) 2026 Michelle Bausager
*
* This file is part of PANIC.
*
* PANIC is free software licensed under the GNU General Public License v3.0 or later.
* You may redistribute and/or modify it under the terms of the GPL.
*
* PANIC is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY;
* without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
* See the LICENSE file for the full license text.
*
* SPDX-License-Identifier: GPL-3.0-or-later
*
*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
*
* Project Name: PANIC
* Module Name: random
* File Name: uniform.cpp
* Revision: 0.1.0
* Date: 29-06-2026
* Author: Michelle Bausager
*
* Description:
* Defines the dense layers used in neural network
*
*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++*/
//---------------------------------------------------------------------------------------------------------------------------
// INCLUDE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
#include <random/uniform.hpp>
#include <config/types.hpp>
#include <config/omp.hpp>
#include <random/seed.hpp>
#include <tensor/vector.hpp>
//---------------------------------------------------------------------------------------------------------------------------
// VARIABLE DESCRIPTION
//---------------------------------------------------------------------------------------------------------------------------
static const panic::types::uint_t uniform_omp_min_work = 5000;
static const panic::types::uint_t max_uint_t = ~static_cast<panic::types::uint_t>(0);
static const panic::types::real_t real_max_unit = static_cast<panic::types::real_t>(max_uint_t);
static panic::types::uint_t uniform_state = static_cast<panic::types::uint_t>(1);
//---------------------------------------------------------------------------------------------------------------------------
// FUNCTION PROTOTYPE
//---------------------------------------------------------------------------------------------------------------------------
namespace panic{
namespace random{
//--------------------------------------------------------------------------------------------------------------------------
// Function Name : panic::random::uniform
//
// Description:
// returns a value from a uniform distribution bewteen 0 and 1
//--------------------------------------------------------------------------------------------------------------------------
panic::types::real_t uniform(){
panic::random::seed_t seed;
uniform_state += static_cast<panic::types::uint_t>(1);
seed.set(uniform_state);
return (static_cast<panic::types::real_t>(seed.state_at(seed.get()) / real_max_unit));
}
//--------------------------------------------------------------------------------------------------------------------------
// Function Name : panic::random::uniform
//
// Description:
// returns a value from a uniform distribution bewteen min and max
//--------------------------------------------------------------------------------------------------------------------------
template <typename T>
T uniform(const T min, const T max){
return (min + static_cast<T>(static_cast<panic::types::real_t>((max - min)) * uniform()));
}
//--------------------------------------------------------------------------------------------------------------------------
// EXPLICIT TEMPLATE INSTANTIATION
//
// The implementation is in this .cpp file.
// Build the overload for the official PANIC numeric types.
//--------------------------------------------------------------------------------------------------------------------------
template panic::types::uint_t
uniform<panic::types::uint_t>(const panic::types::uint_t min,
const panic::types::uint_t mix
);
template panic::types::int_t
uniform<panic::types::int_t>(const panic::types::int_t min,
const panic::types::int_t max
);
template panic::types::real_t
uniform<panic::types::real_t>(const panic::types::real_t min,
const panic::types::real_t max
);
//--------------------------------------------------------------------------------------------------------------------------
// Function Name : panic::random::uniform
//
// Description:
// Fills a vector with a uniform distribution
//--------------------------------------------------------------------------------------------------------------------------
bool uniform(panic::tensor::real_vector& a){
panic::random::seed_t seed;
panic::types::uint_t work = a.size();
uniform_state += work;
seed.set(uniform_state);
PANIC_OMP_PARALLEL_FOR_IF(work > uniform_omp_min_work)
for (panic::types::uint_t i = 0; i < work; ++i){
a[i] = (static_cast<panic::types::real_t>(seed.state_at(i)) / real_max_unit);
}
return true;
}
//--------------------------------------------------------------------------------------------------------------------------
// Function Name : panic::random::uniform
//
// Description:
// Fills a vector with a uniform distribution with limits
//--------------------------------------------------------------------------------------------------------------------------
template <typename T>
bool uniform(panic::tensor::vector<T>& a, const T min, const T max){
panic::random::seed_t seed;
panic::types::uint_t work = a.size();
uniform_state += work;
seed.set(uniform_state);
panic::types::real_t temp = static_cast<panic::types::real_t>(max - min);
PANIC_OMP_PARALLEL_FOR_IF(work > uniform_omp_min_work)
for (panic::types::uint_t i = 0; i < work; ++i){
a[i] = min + static_cast<T>((temp*(static_cast<panic::types::real_t>(seed.state_at(i)) / real_max_unit)));
}
return true;
}
//--------------------------------------------------------------------------------------------------------------------------
// EXPLICIT TEMPLATE INSTANTIATION
//
// The implementation is in this .cpp file.
// Build the overload for the official PANIC numeric types.
//--------------------------------------------------------------------------------------------------------------------------
template bool uniform<panic::types::uint_t>(panic::tensor::vector<panic::types::uint_t>& a,
const panic::types::uint_t min,
const panic::types::uint_t mix
);
template bool uniform<panic::types::int_t>(panic::tensor::vector<panic::types::int_t>& a,
const panic::types::int_t min,
const panic::types::int_t max
);
template bool uniform<panic::types::real_t>(panic::tensor::vector<panic::types::real_t>& a,
const panic::types::real_t min,
const panic::types::real_t max
);
//--------------------------------------------------------------------------------------------------------------------------
// Function Name : panic::random::uniform
//
// Description:
// Fills a matrix with a uniform distribution
//--------------------------------------------------------------------------------------------------------------------------
bool uniform(panic::tensor::real_matrix& A){
panic::random::seed_t seed;
panic::types::uint_t rows = A.rows();
panic::types::uint_t cols = A.cols();
panic::types::uint_t work = rows*cols;
uniform_state += work;
seed.set(uniform_state);
PANIC_OMP_PARALLEL_FOR_IF(work > uniform_omp_min_work)
for (panic::types::uint_t i = 0; i < rows; ++i){
for (panic::types::uint_t j = 0; j < cols; ++j){
panic::types::uint_t index = i * cols + j;
A(i,j) = (static_cast<panic::types::real_t>(seed.state_at(index)) / real_max_unit);
}
}
return true;
}
//--------------------------------------------------------------------------------------------------------------------------
// Function Name : panic::random::uniform
//
// Description:
// Fills a matrix with a uniform distribution with limits
//--------------------------------------------------------------------------------------------------------------------------
template <typename T>
bool uniform(panic::tensor::matrix<T>& A, const T min, const T max){
panic::random::seed_t seed;
panic::types::uint_t rows = A.rows();
panic::types::uint_t cols = A.cols();
panic::types::uint_t work = rows*cols;
uniform_state += work;
seed.set(uniform_state);
panic::types::real_t temp = static_cast<panic::types::real_t>(max - min);
PANIC_OMP_PARALLEL_FOR_IF(work > uniform_omp_min_work)
for (panic::types::uint_t i = 0; i < rows; ++i){
for (panic::types::uint_t j = 0; j < cols; ++j){
panic::types::uint_t index = i * cols + j;
A(i,j) = min + static_cast<T>((temp*(static_cast<panic::types::real_t>(seed.state_at(index)) / real_max_unit)));
}
}
return true;
}
//--------------------------------------------------------------------------------------------------------------------------
// EXPLICIT TEMPLATE INSTANTIATION
//
// The implementation is in this .cpp file.
// Build the overload for the official PANIC numeric types.
//--------------------------------------------------------------------------------------------------------------------------
template bool uniform<panic::types::uint_t>(panic::tensor::matrix<panic::types::uint_t>& a,
const panic::types::uint_t min,
const panic::types::uint_t mix
);
template bool uniform<panic::types::int_t>(panic::tensor::matrix<panic::types::int_t>& a,
const panic::types::int_t min,
const panic::types::int_t max
);
template bool uniform<panic::types::real_t>(panic::tensor::matrix<panic::types::real_t>& a,
const panic::types::real_t min,
const panic::types::real_t max
);
} // namespace random
} // namespace panic