Init project after finished it and start of the report

This commit is contained in:
2026-10-12 19:45:30 +02:00
commit b2dde5453d
13 changed files with 781 additions and 0 deletions
+387
View File
@@ -0,0 +1,387 @@
//
// Created by dukantic on 29/09/2026.
//
#ifndef RABBIT_HOLE_H
#define RABBIT_HOLE_H
#include "rabbit_hole.h"
#include <stdio.h>
#include <math.h>
#include "mt19937ar.h"
__uint128_t i,j;
double x_gauss[2];
void init_rabbits(rabbit_hole_t *rabbits) {
rabbits->female_litter=0;
for (i = 0; i < N_MATURE; i++) {
rabbits->male_mature[i] =0;
rabbits->female_mature[i] =0;
}
for (i = 0; i < N_KITTEN; i++) {
rabbits->kitten_male[i] =0;
rabbits->kitten_female[i] =0;
}
}
double sk_sigma_m, sk_mu_m, sk_sigma_f, sk_mu_f;
double sk_draw;
__uint128_t sk_temp_male, sk_temp_female;
void shift_kitten(rabbit_hole_t *rabbits, __uint128_t *new_male, __uint128_t *new_female) {
(*new_male) = rabbits->kitten_male[N_KITTEN-1];
(*new_female) = rabbits->kitten_female[N_KITTEN-1];
for (i = N_KITTEN - 3; i < N_KITTEN-1; i++) {
sk_temp_male = 0;
sk_temp_female = 0;
if (rabbits->kitten_male[i]<100 || rabbits->kitten_female[i]<100) {
for ( j = 0; j < rabbits->kitten_male[i]; j++) {
sk_draw = genrand_real1();
if (sk_draw <= 0.66)
sk_temp_male ++;
}
for (j = 0; j < rabbits->kitten_female[i]; j++) {
sk_draw = genrand_real1();
if (sk_draw <= 0.66)
sk_temp_female ++;
}
rabbits->kitten_male[i]-= sk_temp_male;
rabbits->kitten_female[i]-= sk_temp_female;
(*new_male) += sk_temp_male;
(*new_female) += sk_temp_female;
}
else
{
gauss(x_gauss);
sk_mu_m = (double)rabbits->kitten_male[i] * 0.66;
sk_sigma_m = sqrt(sk_mu_m * 0.33);
sk_mu_f = (double)rabbits->kitten_female[i] * 0.66;
sk_sigma_f = sqrt(sk_mu_f * 0.33);
x_gauss[0] = CLAMP_BTW(x_gauss[0] * sk_sigma_m + sk_mu_m, 0.0, rabbits->kitten_male[i]);
x_gauss[1] = CLAMP_BTW( x_gauss[1] * sk_sigma_f + sk_mu_f, 0.0, rabbits->kitten_female[i]);
(*new_male) += (__uint128_t)roundl(x_gauss[0]);
(*new_female) += (__uint128_t)roundl(x_gauss[1]);
rabbits->kitten_male[i] -= (__uint128_t)roundl(x_gauss[0]);
rabbits->kitten_female[i] -= (__uint128_t)roundl(x_gauss[1]);
}
}
for (i = N_KITTEN-1; i >=1; i--) {
rabbits->kitten_male[i] = rabbits->kitten_male[i-1];
rabbits->kitten_female[i] = rabbits->kitten_female[i-1];
}
rabbits->kitten_male[0] = 0;
rabbits->kitten_female[0] = 0;
}
__uint128_t *males;
__uint128_t *females;
double draw;
void shift_mature(rabbit_hole_t* rabbits) {
males = rabbits->male_mature;
females = rabbits->female_mature;
for (i = N_MATURE-1; i >=1; i--) {
males[i] = males[i-1];
females[i] = females[i-1];
}
males[0] = 0;
females[0] = 0;
}
double dp_threshold_death, dp_draw_death;
double sigma_m,sigma_f, mu_m, mu_f;
__uint128_t temp_male, temp_female;
void death_process(rabbit_hole_t *rabbits) {
for (i =0; i < N_MATURE; i++) {
dp_threshold_death = 1 - (0.4 + ((i >= 9)? (double)(i-9) * 0.1: 0)) / 12;
if (rabbits->male_mature[i] < 100 || rabbits->female_mature[i] < 100) {
temp_male = 0;
temp_female = 0;
for (j = 0; j < rabbits->male_mature[i]; j++) {
dp_draw_death = genrand_real1();
if (dp_draw_death <= dp_threshold_death)
temp_male ++;
}
for (j = 0; j < rabbits->female_mature[i]; j++) {
dp_draw_death = genrand_real1();
if (dp_draw_death <= dp_threshold_death)
temp_female ++;
}
rabbits->male_mature[i] = temp_male;
rabbits->female_mature[i] = temp_female;
}else {
mu_f = (double)rabbits->female_mature[i] * dp_threshold_death;
sigma_f = sqrt(mu_f * (1.0-dp_threshold_death));
mu_m = (double)rabbits->male_mature[i] * dp_threshold_death;
sigma_m = sqrt(mu_m * (1.0-dp_threshold_death));
gauss(x_gauss);
x_gauss[0] = CLAMP_BTW(x_gauss[0] * sigma_f + mu_f, 0.0, (double)rabbits->female_mature[i]);
x_gauss[1] = CLAMP_BTW(x_gauss[1] * sigma_m + mu_m, 0.0, (double)rabbits->male_mature[i]);
rabbits->female_mature[i] = (__uint128_t)roundl(x_gauss[0]);
rabbits->male_mature[i] = (__uint128_t)roundl(x_gauss[1]);
}
}
for (i = 0; i < N_KITTEN; i++) {
if (rabbits->kitten_male[i] < 100 || rabbits->kitten_female[i] < 100) {
temp_male = 0;
temp_female = 0;
for (j = 0; j < rabbits->kitten_male[i]; j++) {
dp_draw_death = genrand_real1();
if (dp_draw_death > 0.05) //0.35/7
temp_male ++;
}
for (j = 0; j < rabbits->kitten_female[i]; j++) {
dp_draw_death = genrand_real1();
if (dp_draw_death > 0.05) //0.35/7
temp_female ++;
}
rabbits->kitten_male[i] = temp_male;
rabbits->kitten_female[i] = temp_female;
}else {
mu_f = (double)rabbits->kitten_female[i] * 0.95;
sigma_f = sqrt(mu_f * 0.05);
mu_m = (double)rabbits->kitten_male[i] * 0.95;
sigma_m = sqrt(mu_m * 0.05);
gauss(x_gauss);
x_gauss[0] = x_gauss[0] * sigma_f + mu_f;
x_gauss[1] = x_gauss[1] * sigma_m + mu_m;
rabbits->kitten_female[i] = (__uint128_t)roundl(x_gauss[0]);
rabbits->kitten_male[i] = (__uint128_t)roundl(x_gauss[1]);
}
}
}
__uint128_t nb_females;
int taille1;
__uint128_t new_litter;
double draw_litter;
double ul_mu_k,ul_sigma_k, ul_mu_sex, ul_sigma_sex;
double ul_mu_litter, ul_sigma_litter;
void update_litter(rabbit_hole_t *rabbits, __uint128_t *new_kitten_male, __uint128_t *new_kitten_female) {
(*new_kitten_female) = 0;
(*new_kitten_male) = 0;
// sigma par femelle carre = ((b-a+1)caree- 1)/12 = 4
// sigma par femelle = 2
//sigma total = sqrt(sigma par femelle carre * N)
//mu = ((9+3)/2 ) * N = 6 * N
if (rabbits->female_litter < 100) {
for (i = 0; i < rabbits->female_litter; i++) {
taille1 = 3 + (int)(4.0 * genrand_real2()); //TODO
taille1 = 3 + (int)(4.0 * genrand_real2()); //TODO
for (int k = 0; k < taille1; k++) {
if (genrand_real1() > 0.5) (*new_kitten_male)++;
else (*new_kitten_female)++;
}
}
}else {
ul_mu_k = 6.0*(double)rabbits->female_litter;
ul_sigma_k = sqrt(4.0 * (double)rabbits->female_litter);
gauss(x_gauss);
x_gauss[0] = CLAMP_BTW(x_gauss[0] * ul_sigma_k+ ul_mu_k, 0.0 , rabbits->female_litter * 9.0);
ul_mu_sex = x_gauss[0] * 0.5;
ul_sigma_sex = sqrt(ul_mu_sex * 0.5);
x_gauss[1] = CLAMP_BTW(x_gauss[1] * ul_sigma_sex + ul_mu_sex , 0.0, rabbits->female_litter * 9.0);
(*new_kitten_female) = (__uint128_t)roundl(x_gauss[1]) ;
(*new_kitten_male) = (__uint128_t)roundl(x_gauss[0] - x_gauss[1]);
}
nb_females = 0;
for (i = 0; i < N_MATURE; i++) {
nb_females += rabbits->female_mature[i];
}
if (nb_females < 100) {
new_litter = 0;
for (i = 0; i < nb_females; i++) {
draw_litter = genrand_real1();
if (draw_litter > 0.5)
new_litter++;
}
rabbits->female_litter = new_litter;
}else {
ul_mu_litter = (double)nb_females * 0.5;
ul_sigma_litter = sqrt(ul_mu_litter * 0.5);
gauss(x_gauss);
x_gauss[0] = CLAMP_BTW(x_gauss[0]* ul_sigma_litter + ul_mu_litter ,0.0,nb_females);
rabbits->female_litter = (__uint128_t)roundl(x_gauss[0]);
}
}
char total_rabbits(rabbit_hole_t *rabbits, __uint128_t *total) {
*total = 0;
__uint128_t total_temp;
for (i = 0; i < N_MATURE; i++) {
total_temp = *total;
*total += rabbits->male_mature[i];
*total += rabbits->female_mature[i];
if (*total < total_temp) //overflow
{
return 0;
}
}
for (i = 0; i < N_KITTEN; i++) {
total_temp = *total;
*total += rabbits->kitten_male[i];
*total += rabbits->kitten_female[i];
if (*total < total_temp) //overflow
{
return 0;
}
}
return 1;
}
static void u128_to_str(__uint128_t v, char *buf)
{
char tmp[45];
int i = 0, j = 0;
if (v == 0) {
buf[0] = '0';
buf[1] = '\0';
return;
}
while (v > 0) {
tmp[i++] = '0' + (char)(v % 10);
v /= 10;
}
while (i > 0)
buf[j++] = tmp[--i];
buf[j] = '\0';
}
static void print_u128(__uint128_t v)
{
char buf[45];
u128_to_str(v, buf);
printf("%s", buf);
}
static void print_row_u128(const char *label, const __uint128_t *arr, size_t n)
{
fputs(label, stdout);
fputs("| ", stdout);
for (size_t i = 0; i < n; i++) {
print_u128(arr[i]);
fputs(" | ", stdout);
}
putchar('\n');
}
char print_rabbits(rabbit_hole_t *rabbits, __uint128_t *month)
{
printf("During time : ");
print_u128(*month / 12);
printf(" years, ");
print_u128(*month % 12);
printf(" months\n");
print_row_u128("Male : ", rabbits->male_mature, N_MATURE);
print_row_u128("Female : ", rabbits->female_mature, N_MATURE);
print_row_u128("Male kitten : ", rabbits->kitten_male, N_KITTEN);
print_row_u128("Female kitten : ", rabbits->kitten_female, N_KITTEN);
printf("\nLitter : | ");
print_u128(rabbits->female_litter);
printf("\n=> ");
__uint128_t total;
if ( ! total_rabbits(rabbits,&total)) {
puts("======================OVERFLOW===================");
return 0;
}
print_u128(total);
printf("\n");
puts("########################################");
return 1;
}
void csv_header() {
fprintf(stderr, "Month, Total Rabbit, Litters,");
for ( i = 0; i < N_MATURE; i++) {
fprintf(stderr,"Male (age %d),",i+1);
}
for (i = 0; i < N_MATURE; i++) {
fprintf(stderr,"Female (age %d),",i+1);
}
for (i = 0; i < N_KITTEN; i++) {
fprintf(stderr,"Kitten Male (age %d),",i+1);
}
for (i = 0; i < N_KITTEN; i++) {
fprintf(stderr,"Kitten Female (age %d),",i+1);
}
fprintf(stderr, "\n");
}
char csv_rabbits(rabbit_hole_t *rabbits, const __uint128_t *month) {
__uint128_t total;
if (!total_rabbits(rabbits,&total))
return 0;
char number_string[512];
u128_to_str(*month, number_string);
fprintf(stderr, "%s,", number_string);
u128_to_str(total,number_string);
fprintf(stderr, "%s,", number_string);
u128_to_str(rabbits->female_litter,number_string);
fprintf(stderr, "%s,", number_string);
for ( i = 0; i < N_MATURE; i++) {
u128_to_str(rabbits->male_mature[i],number_string);
fprintf(stderr,"%s,", number_string);
}
for (i = 0; i < N_MATURE; i++) {
u128_to_str(rabbits->female_mature[i],number_string);
fprintf(stderr,"%s,", number_string);
}
for (i = 0; i < N_KITTEN; i++) {
u128_to_str(rabbits->kitten_male[i],number_string);
fprintf(stderr,"%s,", number_string);
}
for (i = 0; i < N_KITTEN; i++) {
u128_to_str(rabbits->kitten_female[i],number_string);
fprintf(stderr,"%s,", number_string);
}
fprintf(stderr, "\n");
return 1;
}
__uint128_t new_male;
__uint128_t new_female;
__uint128_t new_kitten_male;
__uint128_t new_kitten_female;
char process_month(rabbit_hole_t *rabbits,__uint128_t *month) {
shift_kitten(rabbits, &new_male, &new_female);
if (*month % 12 == 0) {
shift_mature(rabbits);
}
rabbits->male_mature[0] += new_male;
rabbits->female_mature[0] += new_female;
death_process(rabbits);
update_litter(rabbits, &new_kitten_male, &new_kitten_female);
rabbits->kitten_male[0] += new_kitten_male;
rabbits->kitten_female[0] += new_kitten_female;
(*month)++;
if (!csv_rabbits(rabbits,month))
return 0;
return 1;
}
#endif