Lines Matching +full:petsc +full:- +full:dist
1 #include <petsc/private/dmpleximpl.h>
30 return -0.25 * r5 + 0.5 * r4 + 0.625 * r3 - (5.0 / 3.0) * r2 + 1.0; in GaspariCohn()
33 …return (1.0 / 12.0) * r5 - 0.5 * r4 + 0.625 * r3 + (5.0 / 3.0) * r2 - 5.0 * r + 4.0 - (2.0 / 3.0) … in GaspariCohn()
38 …DMPlexGetLETKFLocalizationMatrix - Compute localization weight matrix for LETKF [move to ml/da/int…
43 + n_obs_vertex - Number of nearest observations to use per vertex (eg, MAX_Q_NUM_LOCAL_OBSERVATIONS…
44 . n_obs_local - Number of local observations
45 . n_dof - Number of degrees of freedom
46 . Vecxyz - Array of vectors containing the coordinates
47 - H - Observation operator matrix
50 . Q - Localization weight matrix (sparse, AIJ format)
55 exactly n_obs_vertex non-zero entries corresponding to the nearest
56 observations, weighted by the Gaspari-Cohn fifth-order piecewise
100 // If n_obs_global < n_obs_vertex, we will pad with -1 indices and 0.0 weights. in DMPlexGetLETKFLocalizationMatrix()
131 …// Host view must match the data layout from VecGetArray (d-major, i-minor implies LayoutLeft for … in DMPlexGetLETKFLocalizationMatrix()
173 /* Temporary storage for top-k per vertex */ in DMPlexGetLETKFLocalizationMatrix()
194 best_idxs_dev(i, k) = -1; in DMPlexGetLETKFLocalizationMatrix()
201 PetscReal diff = v_coords[d] - obs_coords_dev(j, d); in DMPlexGetLETKFLocalizationMatrix()
208 PetscInt pos = n_obs_vertex - 1; in DMPlexGetLETKFLocalizationMatrix()
209 while (pos > 0 && best_dists_dev(i, pos - 1) > dist2) { in DMPlexGetLETKFLocalizationMatrix()
210 best_dists_dev(i, pos) = best_dists_dev(i, pos - 1); in DMPlexGetLETKFLocalizationMatrix()
211 best_idxs_dev(i, pos) = best_idxs_dev(i, pos - 1); in DMPlexGetLETKFLocalizationMatrix()
212 pos--; in DMPlexGetLETKFLocalizationMatrix()
218 current_max_dist = best_dists_dev(i, n_obs_vertex - 1); in DMPlexGetLETKFLocalizationMatrix()
223 PetscReal radius2 = best_dists_dev(i, n_obs_vertex - 1); in DMPlexGetLETKFLocalizationMatrix()
228 if (best_idxs_dev(i, k) != -1) { in DMPlexGetLETKFLocalizationMatrix()
229 PetscReal dist = std::sqrt(best_dists_dev(i, k)); in DMPlexGetLETKFLocalizationMatrix()
231 values_dev(i, k) = GaspariCohn(dist, radius); in DMPlexGetLETKFLocalizationMatrix()
233 indices_dev(i, k) = -1; // Ignore this entry in DMPlexGetLETKFLocalizationMatrix()
240 // Host views must be LayoutRight for MatSetValues (row-major) in DMPlexGetLETKFLocalizationMatrix()