wolfhece.drowning_victims.drowning_functions ============================================ .. py:module:: wolfhece.drowning_victims.drowning_functions Module Contents --------------- .. py:data:: EPS :value: 0.15 .. py:data:: G .. py:data:: P_ATM :value: 101325 .. py:data:: RHO_F :value: 1000 .. py:function:: beta_find(x, *data) Iterative function that finds the parameters alpha and beta defining a beta distribution and scales it in a range. The function fits the parameters based on two given percentiles of the data. :param x: Array containing initial guesses for alpha and beta. :param data: Tuple containing (down, up, mini, maxi, perc1, perc2). :return: Array containing the differences between the calculated and given percentiles. .. py:function:: Body_motion(a_RK, batch_turb, CFL, Delta, epsilon, H_pre, H_post, Human, k, NbX, NbY, Pos_bp, resurface, sinking, time_b, t_Wolf_pre, t_Wolf_post, time_goal, U_bp, Ux_pre, Ux_post, Uy_pre, Uy_post, Z_param) Function calculating the motion of the body at each time step using a Runge-Kutta method. From body position, velocity and flow environment, the function determines the flow velocity at the body position and calculates its new velocities, checking for collisions with walls. :param a_RK: Runge-Kutta coefficient. :param batch_turb: Batch turbulence. :param CFL: Courant-Friedrichs-Lewy number. :param Delta: Array of delta values. :param epsilon: Epsilon value. :param H_pre: Pre-update water depth. :param H_post: Post-update water depth. :param Human: Human parameters array. :param k: Turbulence kinetic energy. :param NbX: Number of cells in X direction. :param NbY: Number of cells in Y direction. :param Pos_bp: Body position array for calculations. :param resurface: Resurface array. :param sinking: Sinking array. :param time_b: Time array for savings. :param t_Wolf_pre: Pre-update Wolf time. :param t_Wolf_post: Post-update Wolf time. :param time_goal: Time goal. :param U_bp: Body velocity array for calculations. :param Ux_pre: Pre-update X velocity. :param Ux_post: Post-update X velocity. :param Uy_pre: Pre-update Y velocity. :param Uy_post: Post-update Y velocity. :param Z_param: Z parameters array. :return: Updated Human, Pos_b, resurface, sinking, time_b, U_b arrays. .. py:function:: Body_temperature(BSA, mass, T_w, t) Gives the body temperature at time t considering ONLY the conduction. :param BSA: Body surface area [m²]. :param mass: Body mass [kg]. :param T_w: Water temperature [°C]. :param t: Time from the beginning. :return: ADD (Accumulated Degree Days) and body temperature. .. py:function:: Body_volume_variation(alpha_1, BSA, FRC, H_f, m_b, time_b, T_w, TLC, V_b0, V_clothes1, V_clothes2, z) Function calculating the body volume variation due to the putrefaction gases. The method is described in Delhez et al. 2025, "Predicting the buoyancy and Postmortem Submersion Interval of victims of river drowning". :param alpha_1: Alpha parameter. :param BSA: Body surface area. :param FRC: Functional residual capacity. :param H_f: Water depth. :param m_b: Body mass. :param time_b: Time array for savings. :param T_w: Water temperature. :param TLC: Total lung capacity. :param V_b0: Initial body volume. :param V_clothes1: Clothes volume for temperature < 15°C. :param V_clothes2: Clothes volume for temperature >= 15°C. :param z: Depth array. :return: Body volume. .. py:function:: Collision(Delta, Pos, Pos_p, U, walls) Correct the bodies that hit a wall by correcting its velocity and replace it on the wall (+marge). :param Delta: Array of delta values. :param Pos: Body position array for calculations. :param Pos_p: Previous body position array. :param U: Body velocity array for calculations. :param walls: Walls array. :return: Corrected position and velocity arrays. .. py:function:: Flow_time_t(batch_turb, dt, epsilon, H_f, k, turb_type, U_bp, U_x_vec, U_y_vec, z, z_0) Calculates the difference between the body velocity and the flow velocity (-> relative body velocity) and its sign before and after evaluating the turbulences. :param batch_turb: Batch turbulence. :param dt: Time step. :param epsilon: Epsilon value. :param H_f: Water depth. :param k: Turbulence kinetic energy. :param turb_type: Turbulence type. :param U_bp: Body velocity array for calculations. :param U_x_vec: X velocity vector. :param U_y_vec: Y velocity vector. :param z: Depth array. :param z_0: Roughness length. :return: U_x_dif, U_x_sign, U_y_dif, U_y_sign arrays. .. py:function:: interp_mailles_mat(Delta, epsilon, H_0, H_1, index_b, k, NbX, NbY, Pos_b, t_WOLF_perc, t_Wolf_perc_insta, Ux_0, Ux_1, Uy_0, Uy_1) Determines the flow velocity and height at the body position based on the value of the cells in which the body is and the next cells in x, y and xy. It is a spatial bi-linear and temporal linear interpolation. Method described in Delhez et al. 2025, Lagrangian modelling of the drift of a victim of river drowning. :param Delta: Array of delta values. :param epsilon: Epsilon value. :param H_0: Pre-update water depth. :param H_1: Post-update water depth. :param index_b: Body position indices. :param k: Turbulence kinetic energy. :param NbX: Number of cells in X direction. :param NbY: Number of cells in Y direction. :param Pos_b: Body position array for savings. :param t_WOLF_perc: Percentage of Wolf time. :param t_Wolf_perc_insta: Instantaneous percentage of Wolf time. :param Ux_0: Pre-update X velocity. :param Ux_1: Post-update X velocity. :param Uy_0: Pre-update Y velocity. :param Uy_1: Post-update Y velocity. :return: du_insta, epsilon_v, ind_walls, H_v, k_v, Ux_v, Uy_v, walls arrays. .. py:function:: known_1(g, mini, maxi, down, up, perc1, perc2) Used in the fit of the beta parameters alpha and beta by iteration. :param g: Guess array. :param mini: Minimum value. :param maxi: Maximum value. :param down: Lower percentile value. :param up: Upper percentile value. :param perc1: Lower percentile. :param perc2: Upper percentile. :return: Tuple containing alpha and beta. .. py:function:: Loading(Path_loading, Pos_b, time_b, U_b) Loads the results of a previous simulation and returns the data needed to start from these results. :param Path_loading: Path to the loading file. :param Pos_b: Body position array for savings. :param time_b: Time array for savings. :param U_b: Body velocity array for savings. :return: count_initial, Human, n_loaded, Pos_b, time_b, U_b, Z_param arrays. .. py:function:: Loop_management(progress_queue, process_id, a_RK, BC_cells, count, count_Wolf, CFL, Delta, Human_np, i_initial, n_b, n_saved, n_t, NbX, NbY, Path_Saving, Path_Wolf, Pos, Pos_b, resurface, sinking, time, time_b, time_goal, U, U_b, wanted_time, wanted_Wolf, Z_param_np) Main loop of the code. Calculates the motion of each body at each time in the loop and updates the flow when needed. Everything is based on the array "still" which contains the index of all the bodies that need more calculations, as we work with variable time step for each body. If a body is out of the domain, it is out of still; if a body reaches a time for which we need a save, it is out of still; if a body reaches a time for which the flow needs to be updated, it is out of still. :param progress_queue: Queue for progress updates. :param process_id: Process ID. :param a_RK: Runge-Kutta coefficient. :param BC_cells: Boundary condition cells. :param count: Count of saved states. :param count_Wolf: Count of Wolf states. :param CFL: Courant-Friedrichs-Lewy number. :param Delta: Array of delta values. :param Human_np: Human parameters array. :param i_initial: Initial index. :param n_b: Number of bodies. :param n_saved: Number of saved states. :param n_t: Number of time steps. :param NbX: Number of cells in X direction. :param NbY: Number of cells in Y direction. :param Path_Saving: Path to save results. :param Path_Wolf: Path to Wolf results. :param Pos: Body position array for calculations. :param Pos_b: Body position array for savings. :param resurface: Resurface array. :param sinking: Sinking array. :param time: Time array for calculations. :param time_b: Time array for savings. :param time_goal: Time goal. :param U: Body velocity array for calculations. :param U_b: Body velocity array for savings. :param wanted_time: Array of wanted times. :param wanted_Wolf: Array of wanted Wolf times. :param Z_param_np: Z parameters array. :return: Updated Pos_b, resurface, sinking, time_b, U_b arrays. .. py:function:: Motion_equations(CAM, CDA, CLA, CSA, dt, du_insta, m_b, mu, U_bp, U_x_dif, U_x_sign, U_y_dif, U_y_sign, V_b, vertical) Calculates the body acceleration and velocity based on the motion equation with the Flow to the body forces (Drag, Side and Lift), Gravity and Buoyancy, Friction with the bottom and added mass effect. Equations described in Delhez et al., 2025 "Lagrangian modelling of the drift of a victim of river drowning" :param CAM: Added mass coefficient. :param CDA: Drag coefficient. :param CLA: Lift coefficient. :param CSA: Side force coefficient. :param dt: Time step. :param du_insta: Instantaneous velocity difference. :param m_b: Body mass. :param mu: Friction coefficient. :param U_bp: Body velocity array for calculations. :param U_x_dif: X velocity difference. :param U_x_sign: X velocity sign. :param U_y_dif: Y velocity difference. :param U_y_sign: Y velocity sign. :param V_b: Body volume. :param vertical: Vertical position. :return: Body acceleration and velocity arrays. .. py:function:: Parallel_loop(args) Used to run the code in Multiprocess. :param args: Arguments for the Loop_management function. :return: Result of the Loop_management function. .. py:function:: Preparation_parallelisation(progress_queue, a_RK, BC_cells, count, count_Wolf, CFL, Delta, Human_np, i_initial, n_b, n_saved, n_parallel, n_t, NbX, NbY, Path_saving, Path_Wolf, Pos, Pos_b, resurface, sinking, time, time_b, time_goal, U, U_b, wanted_time, wanted_Wolf, Z_param_np) Splits the arrays in the number we want to be ran in MultiProcess. :param progress_queue: Queue for progress updates. :param a_RK: Runge-Kutta coefficient. :param BC_cells: Boundary condition cells. :param count: Count of saved states. :param count_Wolf: Count of Wolf states. :param CFL: Courant-Friedrichs-Lewy number. :param Delta: Array of delta values. :param Human_np: Human parameters array. :param i_initial: Initial index. :param n_b: Number of bodies. :param n_saved: Number of saved states. :param n_parallel: Number of parallel processes. :param n_t: Number of time steps. :param NbX: Number of cells in X direction. :param NbY: Number of cells in Y direction. :param Path_saving: Path to save results. :param Path_Wolf: Path to Wolf results. :param Pos: Body position array for calculations. :param Pos_b: Body position array for savings. :param resurface: Resurface array. :param sinking: Sinking array. :param time: Time array for calculations. :param time_b: Time array for savings. :param time_goal: Time goal. :param U: Body velocity array for calculations. :param U_b: Body velocity array for savings. :param wanted_time: Array of wanted times. :param wanted_Wolf: Array of wanted Wolf times. :param Z_param_np: Z parameters array. :return: List of tasks for parallel processing. .. py:function:: Read_Wolf_GPU_mat(i_Wolf, Path_Wolf) Reads the results of WolfGPU at a particular time and returns the water depth and velocities' matrix. :param i_Wolf: Index of the WolfGPU result. :param Path_Wolf: Path to the WolfGPU results. :return: Water depth and velocities' matrix. .. py:function:: Read_Wolf_GPU_metadata(Path_Wolf) Reads the parameters of the WolfGPU simulation and returns the relevant ones. :param Path_Wolf: Path to the WolfGPU results. :return: Boundary condition cells, Wolf time step, grid spacing, water depth, number of cells in X and Y directions, and total simulation time. .. py:function:: Save_wanted_time(Path, Pos_b, U_b, count) Saves the body positions and velocities at a given time. :param Path: Path to save the results. :param Pos_b: Body position array for savings. :param U_b: Body velocity array for savings. :param count: Count of saved states. .. py:function:: Skinfold(n_b, known, Human) Determines the body density based on its age and BMI using Siri's equation. An error percentage body fat is also set, according to Meeuwsen et al., 2010. :param n_b: Number of bodies. :param known: Known parameter. :param Human: Human parameters array. :return: Updated Human parameters and error percentage body fat. .. py:function:: state_of_run(progress_queue, frame, interval) Monitoring function that displays progress every `interval` seconds. :param progress_queue: Queue for progress updates. :param frame: Frame object for GUI updates. :param interval: Time interval for updates. .. py:function:: U_turbulence(batch, dt, epsilon, H, k, turb_type, U_x, U_x_dif, U_y, U_y_dif, z, z_0) Adjust the flow velocity with turbulence. 4 different evaluations are proposed but each uses a log law of the wall to go from depth-averaged velocity to velocity related to the body vertical position. :param batch: Batch turbulence. :param dt: Time step. :param epsilon: Epsilon value. :param H: Water depth. :param k: Turbulence kinetic energy. :param turb_type: Turbulence type. :param U_x: X velocity. :param U_x_dif: X velocity difference. :param U_y: Y velocity. :param U_y_dif: Y velocity difference. :param z: Depth array. :param z_0: Roughness length. :return: Adjusted X and Y velocities.