The expiration time for new job artifacts in CI/CD pipelines is now 30 days (GitLab default). Previously generated artifacts in already completed jobs will not be affected by the change. The latest artifacts for all jobs in the latest successful pipelines will be kept. More information: https://gitlab.lrz.de/help/user/admin_area/settings/continuous_integration.html#default-artifacts-expiration

Commit 7e38204f authored by Benedikt Zoennchen's avatar Benedikt Zoennchen
Browse files

DistMesh on the GPU

parent 2c365466
kernel void computeForces(
__global float2* vertices,
__global int4* edges,
__global float2* lengths,
__global float* scalingFactor,
__global float2* forces)
{
int i = get_global_id(0);
float2 p1 = vertices[edges[i].s0];
float2 p2 = vertices[edges[i].s1];
float2 v = normalize(p1-p2);
// F= ...
float len = lengths[i].s0;
float desiredLen = lengths[i].s1 * (*scalingFactor) * 1.2f;
float lenDiff = (desiredLen - len);
lenDiff = lenDiff > 0 ? lenDiff : 0;
float2 partialForce = v * lenDiff; // TODO;
forces[edges[i].s0] = forces[edges[i].s0] + partialForce; // TODO sync
//forces[edges[i].s0] = (float2) (1.0f, 1.0f);
//forces[edges[i].s1] = forces[edges[i].s1] - partialForce; // TODO sync
}
kernel void moveVertices(__global float2* vertices, __global float2* forces, const float delta) {
int i = get_global_id(0);
//float2 force = (float2)(1.0, 1.0);
float2 force = forces[i];
vertices[i] = vertices[i] + (force * delta);
}
// computation of the scaling factor:
kernel void computeLengths(
__global float2* vertices,
__global int4* edges,
__global float2* lengths,
__global float2* qLengths)
{
int i = get_global_id(0);
float2 p1 = vertices[edges[i].s0];
float2 p2 = vertices[edges[i].s1];
float2 v = p1-p2;
//TODO: desiredLenfunction required
float desiredLen = 1.0f;
float2 len = (float2) (length(v), desiredLen);
lengths[i] = len;
qLengths[i] = (float2) (v.s0*v.s0 + v.s1*v.s1, desiredLen*desiredLen);
}
// kernel for multiple work-groups
kernel void computePartialSF(__const int size, __global float2* qlengths, __local float2* partialSums, __global float2* output) {
int gid = get_global_id(0);
if(gid < size){
int global_index = gid;
float2 accumulator = (float2)(0, 0);
// Loop sequentially over chunks of input vector
while (global_index < size) {
float2 element = qlengths[global_index];
accumulator += element;
global_index += get_global_size(0);
}
int lid = get_local_id(0);
int group_size = get_local_size(0);
float2 len = qlengths[gid];
partialSums[lid] = len;
barrier(CLK_LOCAL_MEM_FENCE);
// group_size has to be a power of 2!
for(int i = group_size/2; i > 0; i>>=1){
if(lid < i && lid + i < size) {
partialSums[lid] += partialSums[lid + i];
}
barrier(CLK_LOCAL_MEM_FENCE);
}
if(lid == 0){
output[get_group_id(0)] = partialSums[0];
}
}
}
// kernel for 1 work-group
kernel void computeCompleteSF(__const int size, __global float2* qlengths, __local float2* partialSums, __global float* scaleFactor) {
int gid = get_global_id(0);
int lid = get_local_id(0);
if(lid < size){
int group_size = get_local_size(0);
partialSums[lid] = qlengths[lid];
barrier(CLK_LOCAL_MEM_FENCE);
// group_size has to be a power of 2!
for(int i = group_size/2; i > 0; i>>=1){
if(lid < i && lid + i < size) {
partialSums[lid] += partialSums[lid + i];
}
barrier(CLK_LOCAL_MEM_FENCE);
}
if(lid == 0) {
*scaleFactor = sqrt(partialSums[0].s0 / partialSums[0].s1);
}
}
}
\ No newline at end of file
#ifdef DOUBLE_FP
#ifdef AMD_FP
#pragma OPENCL EXTENSION cl_amd_fp64 : enable
#else
#ifndef CL_VERSION_1_2
#pragma OPENCL EXTENSION cl_khr_fp64 : enable
#endif
#endif
#define varfloat double
#define _255 255.0
#else
#define varfloat float
#define _255 255.0f
#endif
#ifdef USE_TEXTURE
#define OUTPUT_TYPE __write_only image2d_t
#else
#define OUTPUT_TYPE global uint *
#endif
/**
* For a description of this algorithm please refer to
* http://en.wikipedia.org/wiki/Mandelbrot_set
* @author Michael Bien
*/
kernel void mandelbrot(
const int width, const int height,
const varfloat x0, const varfloat y0,
const varfloat rangeX, const varfloat rangeY,
OUTPUT_TYPE output, global uint *colorMap,
const int colorMapSize, const int maxIterations
) {
unsigned int ix = get_global_id(0);
unsigned int iy = get_global_id(1);
varfloat r = x0 + ix * rangeX / width;
varfloat i = y0 + iy * rangeY / height;
varfloat x = 0;
varfloat y = 0;
varfloat magnitudeSquared = 0;
int iteration = 0;
while ( magnitudeSquared < 4 && iteration < maxIterations ) {
varfloat x2 = x*x;
varfloat y2 = y*y;
y = 2 * x * y + i;
x = x2 - y2 + r;
magnitudeSquared = x2+y2;
iteration++;
}
if ( iteration == maxIterations ) {
#ifdef USE_TEXTURE
write_imagef(output, (int2)(ix, iy), (float4)0);
#else
output[iy * width + ix] = 0;
#endif
} else {
float alpha = (float)iteration / maxIterations;
int colorIndex = (int)(alpha * colorMapSize);
#ifdef USE_TEXTURE
// We could have changed colorMap to a texture + sampler, but the
// unpacking below has minimal overhead and it's kinda interesting.
uint c = colorMap[colorIndex];
write_imageui(output, (int2)(ix, iy), (uint4)(
(c >> 0) & 0xFF,
(c >> 8) & 0xFF,
(c >> 16) & 0xFF,
0xFF
));
#else
output[iy * width + ix] = colorMap[colorIndex];
#endif
// monochrom
//output[iy * width + ix] = 255*iteration/maxIterations;
}
}
\ No newline at end of file
kernel void add(__global float2* a, __global float2* b, __global float2* c) {
int i = get_global_id(0);
c[i] = a[i] * b[i];
}
\ No newline at end of file
package org.vadere.util.geometry.mesh.gen;
import org.jetbrains.annotations.NotNull;
import org.vadere.util.geometry.mesh.inter.IFace;
import org.vadere.util.geometry.shapes.IPoint;
/**
* Created by bzoennchen on 06.09.17.
*/
public class AFace<P extends IPoint> implements IFace<P> {
/**
* One of the half-edges bordering this face.
*/
private int edge;
private int id;
private boolean border;
private boolean destroyed = false;
/**
* Default constructor. To construct a face where you have already some half-edges
* bordering this face.
*
* @param edge one of the half-edges bordering this face.
*/
AFace(@NotNull final int id, @NotNull final int edge) {
this(id, edge, false);
}
AFace(@NotNull final int id,@NotNull final int edge, boolean border) {
this.border = border;
this.edge = edge;
this.id = id;
}
/**
* This constructor can be used for constructing a new face without having
* constructed the bordering half-edges jet.
*/
AFace(@NotNull final int id, @NotNull boolean border) {
this.border = border;
this.edge = -1;
this.id = id;
}
AFace() {
this.border = false;
}
boolean isBorder() {
return border;
}
void destroy() {
setEdge(-1);
destroyed = true;
}
/**
* Sets one of the half-edges bordering this face.
*
* @param edge half-edge bordering this face
*/
void setEdge(final int edge) {
this.edge = edge;
}
int getEdge() {
return edge;
}
int getId() {
return id;
}
boolean isDestroyed() {
return destroyed;
}
}
package org.vadere.util.geometry.mesh.gen;
import org.jetbrains.annotations.NotNull;
import org.vadere.util.geometry.mesh.inter.IHalfEdge;
import org.vadere.util.geometry.shapes.IPoint;
/**
* Created by bzoennchen on 06.09.17.
*/
public class AHalfEdge<P extends IPoint> implements IHalfEdge<P> {
private int id;
/**
* point at the end of the half edge.
*/
private int end;
/**
* next half-edge around the face.
*/
private int next;
/**
* previous half-edge around the face.
*/
private int previous;
private int twin;
/**
* the face the half-edge borders.
*/
private int face;
private boolean destroyed;
protected AHalfEdge(@NotNull final int id, @NotNull final int end, @NotNull final int face) {
this.id = id;
this.end = end;
this.face = face;
this.destroyed = false;
}
protected AHalfEdge(@NotNull final int id, @NotNull final int end) {
this.id = id;
this.end = end;
this.face = -1;
this.destroyed = false;
}
int getFace() {
return face;
}
void setFace(final int face) {
this.face = face;
}
int getEnd() {
return end;
}
boolean hasNext() {
return next != -1;
}
int getNext() {
return next;
}
int getPrevious() {
return previous;
}
int getId() {
return id;
}
int getTwin() { return twin; }
/**
* removes the cyclic pointer structure such that the GC can deleteBoundaryFace these objects.
*/
void destroy() {
setNext(-1);
setPrevious(-1);
setFace(-1);
destroyed = true;
}
boolean isValid() {
return next != -1 && previous != -1 && face != -1;
}
void setPrevious(final int previous) {
this.previous = previous;
}
void setNext(final int next) {
this.next = next;
}
void setTwin(final int twin) {
this.twin = twin;
}
void setEnd(final int end) {
this.end = end;
}
public boolean isDestroyed() {
return destroyed;
}
@Override
public String toString() {
return ""+id;
}
/*
* A half-edge is defined by its end vertex and its face. In a geometry there can not be more than
* one half-edge part of face and ending at end.
*/
/*@Override
public boolean equals(Object o) {
if (this == o) return true;
if (o == null || getClass() != o.getClass()) return false;
PHalfEdge<?> halfEdge = (PHalfEdge<?>) o;
if (!end.equals(halfEdge.end)) return false;
return face != null ? face.equals(halfEdge.face) : halfEdge.face == null;
}
@Override
public int hashCode() {
int result = end.hashCode();
result = 31 * result + (face != null ? face.hashCode() : 0);
return result;
}*/
}
package org.vadere.util.geometry.mesh.gen;
import org.jetbrains.annotations.NotNull;
import org.vadere.util.geometry.mesh.inter.IMesh;
import org.vadere.util.geometry.shapes.IPoint;
import org.vadere.util.triangulation.IPointConstructor;
import java.util.*;
import java.util.function.Predicate;
import java.util.stream.Collectors;
import java.util.stream.Stream;
/**
* Created by bzoennchen on 06.09.17.
*/
public class AMesh<P extends IPoint> implements IMesh<P, AVertex<P>, AHalfEdge<P>, AFace<P>> {
private List<AFace<P>> faces;
private boolean elementRemoved;
private int numberOfVertices;
private int numberOfEdges;
private int numberOfFaces;
//private List<PFace<P>> borderFaces;
private AFace<P> boundary;
private List<AHalfEdge<P>> edges;
private IPointConstructor<P> pointConstructor;
private List<AVertex<P>> vertices;
public AMesh(final IPointConstructor<P> pointConstructor) {
this.faces = new ArrayList<>();
//this.borderFaces = new ArrayList<>();
this.edges = new ArrayList<>();
this.vertices = new ArrayList<>();
this.boundary = new AFace<>(-1, true);
//this.faces.add(boundary);
this.pointConstructor = pointConstructor;
this.elementRemoved = false;
this.numberOfFaces = 0;
this.numberOfEdges = 0;
this.numberOfVertices = 0;
}
@Override
public AHalfEdge<P> getNext(@NotNull AHalfEdge<P> halfEdge) {
return edges.get(halfEdge.getNext());
}
@Override
public AHalfEdge<P> getPrev(@NotNull AHalfEdge<P> halfEdge) {
return edges.get(halfEdge.getPrevious());
}
@Override
public AHalfEdge<P> getTwin(@NotNull AHalfEdge<P> halfEdge) {
return edges.get(halfEdge.getTwin());
}
@Override
public AFace<P> getFace(@NotNull AHalfEdge<P> halfEdge) {
int edgeId = halfEdge.getFace();
if(edgeId == -1) {
if(halfEdge.isDestroyed()) {
throw new IllegalArgumentException(halfEdge + " is already destroyed.");
}
return boundary;
}
else {
return faces.get(halfEdge.getFace());
}
}
@Override
public AHalfEdge<P> getEdge(@NotNull AVertex<P> vertex) {
return edges.get(vertex.getEdge());
}
@Override
public AHalfEdge<P> getEdge(@NotNull AFace<P> face) {
return edges.get(face.getEdge());
}
@Override
public P getPoint(@NotNull AHalfEdge<P> halfEdge) {
return getVertex(halfEdge).getPoint();
}
@Override
public AVertex<P> getVertex(@NotNull AHalfEdge<P> halfEdge) {
return vertices.get(halfEdge.getEnd());
}
// the vertex should not be contained in vertices, only the up/down
@Override
public AVertex<P> getDown(@NotNull AVertex<P> vertex) {
return vertices.get(vertex.getDown());
}
// the vertex should not be contained in vertices, only the up/down
@Override
public void setDown(@NotNull AVertex<P> up, @NotNull AVertex<P> down) {
up.setDown(down.getId());
}
@Override
public P getPoint(@NotNull AVertex<P> vertex) {
return vertex.getPoint();
}
@Override
public AFace<P> getFace() {
return faces.stream().filter(f -> !isDestroyed(f)).findAny().get();
}
@Override
public boolean isBoundary(@NotNull AFace<P> face) {
return face == boundary;
}
@Override
public boolean isBoundary(@NotNull AHalfEdge<P> halfEdge) {
return halfEdge.getFace() == boundary.getId();
}
@Override
public boolean isDestroyed(@NotNull AFace<P> face) {
return face.isDestroyed();
}
@Override
public boolean isDestroyed(@NotNull AHalfEdge<P> edge) {
return edge.isDestroyed();
}
@Override
public boolean isDestroyed(@NotNull AVertex<P> vertex) {
return vertex.isDestroyed();
}
@Override
public void setTwin(@NotNull AHalfEdge<P> halfEdge, @NotNull AHalfEdge<P> twin) {
halfEdge.setTwin(twin.getId());
twin.setTwin(halfEdge.getId());
}
@Override
public void setNext(@NotNull AHalfEdge<P> halfEdge, @NotNull AHalfEdge<P> next) {
halfEdge.setNext(next.getId());
next.setPrevious(halfEdge.getId());
}
@Override
public void setPrev(@NotNull AHalfEdge<P> halfEdge, @NotNull AHalfEdge<P> prev) {
halfEdge.setPrevious(prev.getId());
prev.setNext(halfEdge.getId());
}
@Override
public void setFace(@NotNull AHalfEdge<P> halfEdge, @NotNull AFace<P> face) {
halfEdge.setFace(face.getId());
}
@Override
public void setEdge(@NotNull AFace<P> face, @NotNull AHalfEdge<