slm: OpenCL code base  0.1
trajectory.cl
Go to the documentation of this file.
1 
13 
14 #ifdef KERNEL_INTEGRATE_TRAJECTORY
15 static inline void trajectory_record( __global const float2 *uv_array,
41  __global const bool *mask_array,
42  __global ushort *traj_nsteps_array,
43  __global float *traj_length_array,
44  __global char2 *trajectory_vec,
45  const uint global_id,
46  const uint seed_idx,
47  const float2 current_seed_point_vec )
48 {
49  // Private variables - non-constant within this kernel instance
50  __private uint idx, prev_idx, n_steps=0u;
51  __private float l_trajectory=0.0f, dl=0.0f, dt=DT_MAX;
52  __private float2 uv1_vec, uv2_vec, dxy1_vec, dxy2_vec,
53  vec=current_seed_point_vec, prev_vec=vec, next_vec;
54  // Start
55  prev_vec = vec;
56  idx = get_array_idx(vec);
57  // Loop downstream until the pixel is masked, i.e., we've exited the basin or grid,
58  // or if the streamline is too long (l or n)
59  while (!mask_array[idx] && (l_trajectory<MAX_LENGTH && n_steps<(MAX_N_STEPS-1))) {
60  compute_step_vec(dt, uv_array, &dxy1_vec, &dxy2_vec, &uv1_vec, &uv2_vec,
61  vec, &next_vec, &idx);
62  if (!mask_array[idx]) {
63  if (runge_kutta_step_record(&dt, &dl, &l_trajectory, &dxy1_vec, &dxy2_vec,
64  &vec, &prev_vec, &next_vec, &n_steps, &idx, &prev_idx, trajectory_vec))
65  break;
66  } else {
67  euler_step_record(&dt, &dl, &l_trajectory, uv1_vec,
68  &vec, prev_vec, &n_steps, trajectory_vec);
69  break;
70  }
71  }
72  // Record this final trajectory point and return
73  finalize_trajectory(global_id, n_steps, l_trajectory,
74  traj_nsteps_array, traj_length_array);
75  return;
76 }
77 #endif
78 
79 
80 #ifdef KERNEL_INTEGRATE_TRAJECTORY
81 static inline void trajectory_jittered( __global const float2 *uv_array,
108  __global const bool *mask_array,
109  __global uint *slc_array,
110  __global uint *slt_array,
111  const uint global_id,
112  const uint seed_idx,
113  const float2 current_seed_point_vec,
114  const uint initial_rng_state )
115 {
116  // Private variables - non-constant within this kernel instance
117  __private uint idx, prev_idx, n_steps=0u, rng_state=initial_rng_state;
118  __private float l_trajectory=0.0f, dl=0.0f, dt=DT_MAX;
119  __private float2 uv1_vec, uv2_vec, dxy1_vec, dxy2_vec,
120  vec=current_seed_point_vec, prev_vec, next_vec;
121  // Start by recording the seed point
122  prev_vec = vec;
123  idx = get_array_idx(vec);
124  if (mask_array[idx]) return;
125  atomic_write_sl_data(&slt_array[idx], &slc_array[idx], l_trajectory);
126  // Loop downstream until the pixel is masked, i.e., we've exited the basin or grid,
127  // or if the streamline is too long (in l_trajectory or n_steps)
128  while (!mask_array[idx] && (l_trajectory<MAX_LENGTH && n_steps<(MAX_N_STEPS-1))) {
129  compute_step_vec_jittered(dt, uv_array, &rng_state, &dxy1_vec, &dxy2_vec,
130  &uv1_vec, &uv2_vec, vec, &next_vec, &idx);
131  if (!mask_array[idx]) {
132  if (runge_kutta_step_write_sl_data(&dt, &dl, &l_trajectory,
133  &dxy1_vec, &dxy2_vec,
134  &vec, &prev_vec, next_vec,
135  &n_steps, &idx, &prev_idx,
136  mask_array, slt_array, slc_array))
137  break;
138  } else {
139  euler_step_write_sl_data(&dt, &dl, &l_trajectory, uv1_vec,
140  &vec, prev_vec, &n_steps, &idx, &prev_idx,
141  mask_array, slt_array, slc_array);
142  break;
143  }
144 
145  }
146  return;
147 }
148 #endif
149 
#define DT_MAX
Definition: info.h:48
static void euler_step_record(float *dt, float *dl, float *l_trajectory, const float2 uv_vec, float2 *vec, const float2 prev_vec, uint *n_steps, __global char2 *trajectory_vec)
Compute a single Euler integration step of of a streamline.
static void atomic_write_sl_data(__global uint *slt, __global uint *slc, const float l_trajectory)
Add the current streamline length (l_trajectory) to the current pixel of the slt accumulation array...
Definition: writearray.cl:31
static void trajectory_jittered(__global const float2 *uv_array, __global const bool *mask_array, __global uint *slc_array, __global uint *slt_array, const uint global_id, const uint seed_idx, const float2 current_seed_point_vec, const uint initial_rng_state)
Integrate a jittered flow path downstream or upstream.
Definition: trajectory.cl:107
static bool runge_kutta_step_record(float *dt, float *dl, float *l_trajectory, float2 *dxy1_vec, float2 *dxy2_vec, float2 *vec, float2 *prev_vec, float2 *next_vec, uint *n_steps, uint *idx, uint *prev_idx, __global char2 *trajectory_vec)
Compute a single step of 2nd-order Runge-Kutta numerical integration of a streamline given precompute...
static bool runge_kutta_step_write_sl_data(float *dt, float *dl, float *l_trajectory, float2 *dxy1_vec, float2 *dxy2_vec, float2 *vec, float2 *prev_vec, const float2 next_vec, uint *n_steps, uint *idx, uint *prev_idx, __global const bool *mask_array, __global uint *slt_array, __global uint *slc_array)
Compute a single step of 2nd-order Runge-Kutta numerical integration of a streamline given precompute...
static void compute_step_vec(const float dt, const __global float2 *uv_array, float2 *dxy1_vec, float2 *dxy2_vec, float2 *uv1_vec, float2 *uv2_vec, const float2 vec, float2 *next_vec, uint *idx)
Compute a 2nd-order Runge-Kutta integration step along a streamline.
Definition: computestep.cl:41
static void finalize_trajectory(const uint global_id, uint n_steps, float l_trajectory, __global ushort *traj_nsteps_array, __global float *traj_length_array)
Record the (final) trajectory length and count of integration steps to global arrays traj_length_arra...
static uint get_array_idx(float2 vec)
Compute the array index of the padded grid pixel pointed to by a float2 grid position vector (choice ...
Definition: essentials.cl:62
#define MAX_LENGTH
Definition: info.h:44
static void euler_step_write_sl_data(float *dt, float *dl, float *l_trajectory, const float2 uv_vec, float2 *vec, const float2 prev_vec, uint *n_steps, uint *idx, uint *prev_idx, __global const bool *mask_array, __global uint *slt_array, __global uint *slc_array)
Compute a single Euler integration step of of a streamline.
#define MAX_N_STEPS
Definition: info.h:49
static void compute_step_vec_jittered(const float dt, const __global float2 *uv_array, uint *rng_state, float2 *dxy1_vec, float2 *dxy2_vec, float2 *uv1_vec, float2 *uv2_vec, const float2 vec, float2 *next_vec, uint *idx)
Compute a jittered 2nd-order Runge-Kutta integration step along a streamline.
Definition: computestep.cl:85
static void trajectory_record(__global const float2 *uv_array, __global const bool *mask_array, __global ushort *traj_nsteps_array, __global float *traj_length_array, __global char2 *trajectory_vec, const uint global_id, const uint seed_idx, const float2 current_seed_point_vec)
Integrate a streamline downstream or upstream; record the trajectory.
Definition: trajectory.cl:40