Skip to content

Use FMA and reduce repeat computation #23

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
JulianKnodt wants to merge 1 commit into
base: main
Choose a base branch
from
Open

Use FMA and reduce repeat computation #23

Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
36 changes: 31 additions & 5 deletions src/cumc.cu
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -353,7 +353,9 @@ namespace cumc
Scalar d0 = data[v0];
Scalar d1 = data[v1];

Scalar t = (d1 != d0) ? clamp((iso - d0) / (d1 - d0), Scalar(0.0), Scalar(1.0)) : Scalar(0.5);
const Scalar t = (d1 != d0)
? clamp((iso - d0) / (d1 - d0), Scalar(0.0), Scalar(1.0)) \
: Scalar(0.5);

Vertex<Scalar> p0 = {Scalar(x), Scalar(y), Scalar(z)};
Vertex<Scalar> p1 = {Scalar(x + offset[0]), Scalar(y + offset[1]), Scalar(z + offset[2])};
Expand All @@ -364,9 +366,31 @@ namespace cumc
p1 += deform[v1];
}

return p0 + (p1 - p0) * t;
return fma_vertex(t, p1, fma_vertex(-t, p0, p0));
}

__device__ Vertex<double> fma_vertex(
const double t, const Vertex<double> & a, const Vertex<double> & b
) {
Vertex<double> out;
out.x = fma(t, a.x, b.x);
out.y = fma(t, a.y, b.y);
out.z = fma(t, a.z, b.z);
return out;
}

__device__ Vertex<float> fma_vertex(
const float t, const Vertex<float> & a, const Vertex<float> & b
) {
Vertex<float> out;
out.x = fmaf(t, a.x, b.x);
out.y = fmaf(t, a.y, b.y);
out.z = fmaf(t, a.z, b.z);
return out;
}



template <typename Scalar, typename IndexType>
__global__ void create_cell_mc_verts_kernel(CuMC<Scalar, IndexType> mc,
Scalar const *__restrict__ data,
Expand Down Expand Up @@ -421,7 +445,8 @@ namespace cumc
IndexType v1 = mc.gA(x + offset[0], y + offset[1], z + offset[2]);
Scalar d0 = data[v0];
Scalar d1 = data[v1];
Scalar t = (d1 != d0) ? clamp((iso - d0) / (d1 - d0), Scalar(0.0), Scalar(1.0)) : Scalar(0.5);
const Scalar delta = d1 - d0;
const Scalar t = (d1 != d0) ? clamp((iso - d0) / delta, Scalar(0.0), Scalar(1.0)) : Scalar(0.5);

Vertex<Scalar> p0 = {Scalar(x), Scalar(y), Scalar(z)};
Vertex<Scalar> p1 = {Scalar(x + offset[0]), Scalar(y + offset[1]), Scalar(z + offset[2])};
Expand All @@ -445,8 +470,9 @@ namespace cumc
atomicAdd(&adj_deform[v1].z, adj_p1.z);
}

Scalar adj_d0 = (iso - d1) / ((d1 - d0) * (d1 - d0)) * adj_t;
Scalar adj_d1 = (d0 - iso) / ((d1 - d0) * (d1 - d0)) * adj_t;
const auto delta_sqr = delta * delta;
Scalar adj_d0 = (iso - d1) / delta_sqr * adj_t;
Scalar adj_d1 = (d0 - iso) / delta_sqr * adj_t;

atomicAdd(&adj_data[v0], adj_d0);
atomicAdd(&adj_data[v1], adj_d1);
Expand Down