My laptop battery dies quickly these days. Certain apps (cough, cough, Slack), have very high idle CPU-usage. You can *pause* these programs with

```
killall -STOP Slack
```

And later you can *resume* the application with

```
killall -CONT Slack
```

In a few words, explain what this blog is about…or else

My laptop battery dies quickly these days. Certain apps (cough, cough, Slack), have very high idle CPU-usage. You can *pause* these programs with

```
killall -STOP Slack
```

And later you can *resume* the application with

```
killall -CONT Slack
```

The SIGGRAPH Latex style changed to the Libertine font. Here’re the steps to convince Latexit to use the new stylesheet and then to convince Illustrator to use the libertine font for drag and drop math.

```
mkdir ~/Library/texmf/tex/latex/local/acmart.cls/
cp ~/Dropbox/boundary/Paper/acmart.cls ~/Library/texmf/tex/latex/local/acmart.cls
```

In Latexit, open up Preferences, add a new SIGGRAPH “Template” containing:

```
\documentclass[sigconf, review]{acmart}
\pagenumbering{gobble}
```

If you try to drag and drop these into illustrator you’ll see that illustrator has replaced the nice math font with Myriad or something silly.

Download the TTF libertine font pack

Drag this into FontBook.app

```
cp /usr/local/texlive/2015/texmf-dist/fonts/type1/public/libertine/*.pfb ~/Library/Application\ Support/Adobe/Fonts/
```

**Update:** I also had to issue:

```
cp /usr/local/texlive/2015/texmf-dist/fonts/type1/public/txfonts/*.pfb ~/Library/Application\ Support/Adobe/Fonts/
cp /usr/local/texlive/2015/texmf-dist/fonts/type1/public//newtx/*.pfb ~/Library/Application\ Support/Adobe/Fonts/
```

If you see boxes with X’s replacing symbols after dragging and dropping from LaTeXit, then drag into Finder instead (to create a .pdf file), then open this directly and Illustrator will give a warning and tell you which font it’s (still) missing.

Restart Adobe Illustrator

**Update:**

I also added:

```
cp /usr/local/texlive/2017/texmf-dist/fonts/type1/public/amsfonts/cm/*.pfb ~/Library/Application\ Support/Adobe/Fonts/
```

I wrote a small mex wrapper for the *graph segmentation* part of the “Graph Based Image Segmentation” code. Most of the previously matlab implementations/wrappers worked on images. I want to apply this to geometry so I needed access to the graph segmentation directly. Here’s the wrapper (soon to be part of gptoolbox):

```
// mexopts = gptoolbox_mexopts('Static',false,'Debug',true);
// mex('segment_graph.cpp',mexopts{:});
#ifdef MEX
# include <mex.h>
# include <igl/C_STR.h>
# include <igl/matlab/mexErrMsgTxt.h>
# undef assert
# define assert( isOK ) ( (isOK) ? (void)0 : (void) ::mexErrMsgTxt(C_STR(__FILE__<<":"<<__LINE__<<": failed assertion `"<<#isOK<<"'"<<std::endl) ) )
#endif
#include "segment-graph.h"
#include <igl/matlab/mexErrMsgTxt.h>
#include <igl/matlab/parse_rhs.h>
#include <igl/unique.h>
#include <igl/matlab/prepare_lhs.h>
#include <igl/matlab/requires_arg.h>
#include <igl/matlab/validate_arg.h>
#include <igl/matlab/MexStream.h>
#include <Eigen/Sparse>
void mexFunction(
int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
using namespace igl::matlab;
using namespace Eigen;
using namespace std;
igl::matlab::MexStream mout;
std::streambuf *outbuf = std::cout.rdbuf(&mout);
mexErrMsgTxt(nrhs>0,"Too few inputs");
mexErrMsgTxt(mxIsSparse(prhs[0]),"Matrix should be sparse");
const mxArray * mx_data = prhs[0];
const int m = mxGetM(mx_data);
const int n = mxGetN(mx_data);
mexErrMsgTxt(n == mxGetM(prhs[0]), "Matrix should be square");
assert(mxIsSparse(mx_data));
assert(mxGetNumberOfDimensions(mx_data) == 2);
// TODO: It should be possible to directly load the data into the sparse
// matrix without going through the triplets
// Copy data immediately
double * pr = mxGetPr(mx_data);
mwIndex * ir = mxGetIr(mx_data);
mwIndex * jc = mxGetJc(mx_data);
const int num_edges = mxGetNzmax(mx_data);
edge * edges = new edge[num_edges];
int k = 0;
for(int j=0; j<n;j++)
{
// Iterate over inside
while(k<(int)jc[j+1])
{
//cout<<ir[k]<<" "<<j<<" "<<pr[k]<<endl;
assert((int)ir[k]<m);
assert((int)j<n);
edges[k].a = ir[k];
edges[k].b = j;
edges[k].w = pr[k];
k++;
}
}
// defaults
int min_size = 0;
// Threshold
int c = sqrt((double)n);
{
int i = 1;
while(i<nrhs)
{
mexErrMsgTxt(mxIsChar(prhs[i]),"Parameter names should be strings");
// Cast to char
const char * name = mxArrayToString(prhs[i]);
if(strcmp("Threshold",name) == 0)
{
requires_arg(i,nrhs,name);
validate_arg_scalar(i,nrhs,prhs,name);
validate_arg_double(i,nrhs,prhs,name);
c = (double)*mxGetPr(prhs[++i]);
}else if(strcmp("MinSize",name) == 0)
{
requires_arg(i,nrhs,name);
validate_arg_scalar(i,nrhs,prhs,name);
validate_arg_double(i,nrhs,prhs,name);
min_size = (int)((double)*mxGetPr(prhs[++i]));
}
i++;
}
}
universe *u = segment_graph(n, num_edges, edges, c);
// post process small components
for (int i = 0; i < num_edges; i++) {
int a = u->find(edges[i].a);
int b = u->find(edges[i].b);
if ((a != b) && ((u->size(a) < min_size) || (u->size(b) < min_size)))
u->join(a, b);
}
switch(nlhs)
{
case 1:
{
plhs[0] = mxCreateDoubleMatrix(m,1, mxREAL);
Eigen::VectorXi C(m);
for(int i = 0;i<m;i++)
{
C(i) = u->find(i);
}
Eigen::VectorXi uC,I,J;
igl::unique(C,uC,I,J);
prepare_lhs_index(J,plhs);
}
default: break;
}
delete[] edges;
delete u;
std::cout.rdbuf(outbuf);
}
```

It takes the graph as a *sparse matrix* and outputs the component ids:

```
C = segment_graph(A);
```

A while ago, I tried to get glfw playing nicely with matlab’s mex files. I didn’t quite succeed.

My current project requires “GP-GPU” programming. I’m currently using glfw’s “background window” feature to create an OpenGL context. My first tests show that matlab’s mex will play nicely **if**`glfwTerminate()`

is *never called*:

```
#include <igl/opengl/glfw/background_window.h>
#include <mex.h>
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
GLFWwindow * window;
igl::opengl::glfw::background_window(window);
glfwDestroyWindow(window);
//glfwTerminate();
}
```

You can compile this (on mac os x) with:

```
mex( ...
'-largeArrayDims','-DMEX','CXXFLAGS=$CXXFLAGS -std=c++11', ...
'LDFLAGS=\$LDFLAGS -framework Foundation -framework AppKit -framework Accelerate -framework OpenGL -framework AppKit -framework Carbon -framework QuartzCore -framework IOKit ', ...
'-I/usr/local/igl/libigl//external/nanogui/ext/eigen/', ...
'-I/usr/local/igl/libigl/include', ...
'-I/usr/local/igl/libigl/external/nanogui/ext/glfw/include/', ...
'-L/usr/local/igl/libigl/lib/','-lglfw3', ...
'background_glfw.cpp');
```

I don’t know *why* including `glfwTerminate()`

causes matlab to *sometimes* crash. The error report is impossible to read and seems to have something to do with Quartz.

Here’s a little demonstration of how to use gptoolbox and MATLAB to generate an offset surfaces of a triangle mesh. This takes a mesh in `V`

,`F`

and creates a mesh `SV`

,`SF`

of the isosurface at signed distance `iso`

:

```
% Extract offset at minus 3% of bounind box diagonal length
iso = -0.03;
% Resolution grid → resolution of extracted offset surface
side = 60;
% Amount of smoothing to apply to distance field
sigma = 1.4;
bbd = norm(max(V)-min(V));
% Pad generously for positive iso values
[BC,side,r] = voxel_grid([V;max(V)+iso*1;min(V)-iso*1],side);
D = signed_distance(BC,V,F);
D = reshape(D,side([2 1 3]));
% Smooth signed distance field
D = imfilter(D,fspecial('gaussian',9,sigma),'replicate');
BC3 = reshape(BC,[side([2 1 3]) 3]);
% Use matlab's built-in marching cubes iso-surface mesher (works most of the time)
surf = isosurface(BC3(:,:,:,1),BC3(:,:,:,2),BC3(:,:,:,3),D,iso*bbd);
SV = surf.vertices;
SF = surf.faces;
```

Here’s a blue bunny with a positive offset surface, an orange “cage”:

Here’s a blue bunny with a negative offset surface. This is useful for hollowing out objects to 3d print:

Because the iso-surface extraction will over tesselate low curvature patches of the output surface, it would make a lot of sense to remesh/decimate this mesh.

(to create these fancy renderings:)

```
clf;
hold on;
t = tsurf(F,V,'EdgeColor','none',fsoft, 'FaceVertexCData',repmat(blue,size(V,1),1),'FaceAlpha',1+(iso<0)*(0.35-1),fphong);
ts = tsurf(SF,SV,'EdgeAlpha',0.2+(iso<0)*(0-0.2),fsoft,'FaceVertexCData',repmat(orange,size(SV,1),1),fphong,'FaceAlpha',1+(iso>0)*(0.2-1));
apply_ambient_occlusion(ts);
hold off;
axis equal;
view(-20,20)
camlight;
t.SpecularStrength = 0.04;
l = light('Position',[5 -5 10],'Style','local');
add_shadow(t,l,'Color',0.8*[1 1 1],'Fade','local','Ground',[0 0 -1 min([V(:,3);SV(:,3)])]);
set(gca,'pos',[0 0 1 1])
set(gca,'Visible','off');
set(gcf,'Color','w');
drawnow;
```

```
git branch | grep -v "^\\*" | xargs git branch -D
```

This morning I played around with a very simple way of quad meshing a triangle mesh. This is more or less the “ideal pipeline” for quad meshing via parameterization. It will likely only work for very simple meshes. I’m pleasantly surprised that it works at all.

```
% Load a mesh with boundary and no topological (donut) holes
[V,F] = load_mesh('~/Dropbox/models/beetle.off');
V = V*axisangle2matrix([1 0 0],-pi/2);
% Construct the LSCM matrix
[~,Q] = lscm(V,F,[],[]);
M = massmatrix(V,F);
% Use Fiedler vector as parameterization
[EV,ED] = eigs(Q,repdiag(M,2),5,'sm');
U = reshape(EV(:,end-3),[],2);
% try to find "canonical" rotation
[sU,sS,sV] = svd(U'*U);
U = U*sU;
% Lay down a grid in 2D
h = ((max(U(:,1))-min(U(:,1)))/(160-2));
[X,Y] = meshgrid(min(U(:,1))-h:h:max(U(:,1))+h,min(U(:,2))-h:h:max(U(:,2))+h);
% Find all quads that are strictly inside the parameterized mesh
[QQ,QU] = surf2patch(X,Y,0*Y);QU = QU(:,1:2);
[sqrD,I,C] = signed_distance([QU QU(:,1)*0],[U U(:,1)*0],F);
QQ = QQ(all(abs(sqrD(QQ))<1e-10,2),:);
[QU,I] = remove_unreferenced(QU,QQ);
QQ = fliplr(I(QQ));
% Snap boundary vertices of quad mesh to boundary of mesh, and
% solve little Dirichlet problem to diffuse the displacements smoothly
QF = [QQ(:,[1 2 3]);QQ(:,[1 3 4])];
L = cotmatrix(QU,QF);
b = unique(outline(QQ));
[~,~,QU(b,:)] = point_mesh_squared_distance(QU(b,:),U,outline(F));
QU = min_quad_with_fixed(-L,[],b,QU(b,:));
% Locate each quad mesh vertex in the parameterized mesh
[sqrD,I,C] = signed_distance([QU QU(:,1)*0],[U U(:,1)*0],F);
B = barycentric_coordinates(C(:,1:2),U(F(I,1),:),U(F(I,2),:),U(F(I,3),:));
% Interpolate 3D positions at quad mesh vertex locations
QV = V(F(I,1),:).*B(:,1) + V(F(I,2),:).*B(:,2) + V(F(I,3),:).*B(:,3);
```

To create the image above use:

```
clf;
hold on;
tq = tsurf(QQ,QV,'FaceVertexCData',repmat(blue,size(QV,1),1),'EdgeColor','none',fsoft,fphong);
tt = tsurf(F,V-[0.8 0 0],'FaceVertexCData',repmat(1-0.8*(1-blue),size(V,1),1),'EdgeAlpha',0.5,fsoft,fphong);
axis equal;
camproj('persp');
view(-38,17);
drawn;
camlight;
plot_edges(QV,[QQ(:,2:3);QQ(:,[4 1])],'Color',orange);
plot_edges(QV,[QQ(:,1:2);QQ(:,3:4)],'w');hold off;
l = light('Position',[-10 0 10],'Style','infinite');
sq = add_shadow(tq,l,'Fade','infinite','Color',[0.8 0.8 0.8]);
st = add_shadow(tt,l,'Fade','infinite','Color',[0.8 0.8 0.8]);
set(gca,'Visible','off');
set(gcf,'Color','w');
```

There’s a lot to improve about the quad mesh: the boundaries are not aligned with the meshing direction, the meshing directions are not always aligned with principle curvature, etc. Not bad for 30 mins of coding using tools that are just lying around.

Here’s my tested plan for using git and github to manage student assignments for my upcoming course.

We’re going to create a student repo and an instructor fork for *each assignment*. The motivation for separating assignments is to ensure that the diffs that students are submitting for their solutions are tidy. In my case, there’s not much shared code between assignments

Create **public*** repo that students will see: `test-gp-hw0-student`

.

If you try to fork this on Github, it will tell you that you already have a fork. So instead trick github into forking it via the github importer, import `test-gp-hw0-student`

as a **private** `test-gp-hw0-instructor`

. Now you have two repos. We’ll set up a branch on `test-gp-hw0-instructor`

to track `test-gp-hw0-student`

for easy merging:

```
INSTRUCTOR_REPO=test-gp-hw0-instructor
STUDENT_REPO=test-gp-hw0-student
git clone https://github.com/alecjacobson/$INSTRUCTOR_REPO.git
cd $INSTRUCTOR_REPO
git remote add student https://github.com/alecjacobson/$STUDENT_REPO.git
git fetch student
git checkout -b student-tracker --track student/master
git checkout master
cd ..
```

To test this out, let’s create a file visible to everyone (students and instructors) in the student repo:

```
git clone https://github.com/alecjacobson/$STUDENT_REPO.git
cd $STUDENT_REPO
echo "visible to everyone" >> hello-everyone.md
git add hello-everyone.md
git commit -m "add file visible to everyone" hello-everyone.md
git push
cd ../
```

Now, let’s merge change into the instructor repo

```
cd $INSTRUCTOR_REPO
git checkout master
git pull student master
cat hello-everyone.md
cd ..
```

Just to be thorough, let’s try the reverse: creating a file just visible to instructors:

```
cd $INSTRUCTOR_REPO
echo "visible to just instructors" >> hello-instructors.md
git add hello-instructors.md
git commit -m "add file visible to just instructors" hello-instructors.md
git push
cd ..
```

Open an incognito/private window and create a new github account for a fake student (pro tip: gmail lets you put dots anywhere in your email address so that alec.jacobson@gmail.com redirects to alecjacobson@gmail.com, but most websites don’t recognize that these are the same). My fake students name is `stewartdent`

.

Clone the repository and try to access the instructor repo (assuming the url was leaked):

```
git clone https://stewartdent@github.com/stewartdent/$STUDENT_REPO.git stewartdent
cd stewartdent
git remote add instructor https://stewartdent@github.com/alecjacobson/$INSTRUCTOR_REPO.git
git fetch instructor
```

You should see an error:

```
remote: Repository not found.
fatal: repository 'https://stewartdent@github.com/alecjacobson/test-gp-hw0-instructor.git/' not found
```

Now, let’s act as the student to submit the homework

```
# cd stewartdent
echo "stewartdent was here" >> assignment-01.md
git add assignment-01.md
git commit -m "Stewart Dent's Assignment 01 submission" assignment-01.md
git push
cd ..
```

The student has pushed their submission onto their own repo. They can continue to work on it using git in all its glory. When they’re down they “hand in” their homework via a pull request onto the public `test-gp-hw0-student`

repo. The student does this by visiting the Github page for their own fork of `test-gp-hw0-student`

and clicking “New Pull Request”.

The instructor can now grade this. I’m imagining that the `test-gp-hw0-instructor`

repo has a solution that can be compared (in some way or another) to the students work. To pull the students work into the instructor copy (but never merge):

```
cd $INSTRUCTOR_REPO
git checkout -b stewartdent-master master
git pull https://github.com/stewartdent/test-gp-hw0-student.git master
cat empty-form.md
# Record grade.
git checkout master
git branch -D stewartdent-master
```

Then the instructor can close the pull request on `test-gp-hw0-student`

.

Questions or discussions related to the homework can be conducted on the “issues” page of `test-gp-hw0-student`

* The fact that this repo is public should immediately indicate that student’s work submitted via pull request will also be public and importantly will be visible to other students. Unlike the most common form of cheating, this means that Student A can copy the work of Student B without the explicit consent of Student B. Student B’s only possible defense against this would be to submit their work just before the deadline. As far as I know, I can’t hide pull requests from the public (that would be a nice solution) and I prefer not to use any method that relies on student’s having to create a private fork (since they cost money; although, other class require very over-priced textbooks, so a 6-month github subscription is not so crazy, but it does go against the open-source spirit; there’s also the Big Tobacco style get’em while their young https://education.github.com/pack).

The `hypot`

function in matlab purports to be more numerically stable at computing the hypotenuse of a (right-)triangle in 2D. The example the `help hypot`

gives is:

```
a = 3*[1e300 1e-300];
b = 4*[1e300 1e-300];
c1 = sqrt(a.^2 + b.^2)
c2 = hypot(a,b)
```

where you see as output:

```
c1 =
Inf 0
```

and

```
c2 =
5e+300 5e-300
```

this is a *compiled* built-in function so you can’t just `open hypot`

to find out what its doing. It might just be *pre-dividing* by the maximum absolute value. There’s probably a better reference for this, but I found it in: “Vector length and normalization difficulties” by Mike Day.

Continuing this example:

```
m = max(abs([a;b]))
c3 = m.*sqrt((a./m).^2 + (b./m).^2)
```

produces

```
c3 =
5e+300 5e-300
```

While matlab’s hypot only accepts two inputs, pre-dividing by the maximum obviously extends to any dimension.

My typical go-to formula for the area of a triangle *with vertices in nD* is Heron’s formula. Actually this is a formula for the area of a triangle *given its side lengths*. Presumably it’s easy to compute side lengths in any dimension: e.g., a = sqrt((B-C).(B-C)).

Kahan showed that the vanilla Heron’s formula is numerically poor if your (valid) triangle side lengths are given as floating point numbers. He improved this formula by sorting the side-lengths and carefully rearranging the terms in the formula to reduce roundoff.

However, if you’re *vertices* don’t actually live in R^n but rather F^n where F is the fixed precision floating-point grid, then computing side lengths will in general incur some roundoff error. And it may so happen (and it does so happen) that this roundoff error *invalidates* the triangle side lengths. Here, *invalid* means that the side lengths don’t obey the triangle inequality. This might happen for nearly degenerate (zero area) triangles: one (floating point) side length ends up longer than the sum of the others. Kahan provides an assertion to check for this, but there’s no proposed remedy that I could find. One could just declare that these edge-lengths must of come from a zero-area triangle. But returning zero might let the user happily work with very nonsensical triangle side lengths in other contexts *not coming from an embedded triangle*. You could have two versions: one for “known embeddings” (with assertions off and returning 0) and one for “intrinsic/metric only” (with assertions on and returning NaN). Or you could try to fudge in an epsilon a nd hope you choose it above the roundoff error threshold but below the tolerance for perceived “nonsense”. These are messy solutions.

An open question (open in the personal sense of “I don’t know the answer”) is what is the most robust way to compute triangle area in nD *with floating point vertex positions*. A possible answer might be: compute the darn floating point edge lengths, use Kahan’s algorithm and replace NaNs with zeros. That seems unlikely to me because (as far as I can tell) there are 4 sqrt calls, major sources of floating point error.

Re-reading Shewchuk’s robustness notes, I saw his formula for triangle area for points A,B,C in 3D:

Area3D(A,B,C) = sqrt( Area2D(Axy,Bxy,Cxy)^2 + Area2D(Ayz,Byz,Cyz)^2 + Area2D(Azx,Bzx,Czx)^2 ) / 2

where Area2D(A,B,C) is computed as the determinant of [[A;B;C] 1]. This formula reduces the problem of computing 3D triangle area into computing 2D triangle area on the “shadows” (projections) of the triangle on each axis-aligned plane. This lead me to think of a natural generalization to nD:

AreaND(A,B,C) = sqrt( ∑_{i<j} ( Area2D(Aij,Bij,Cij)^2 ) / 2

This formula computes the area of the shadow on all (N choose 2) axis-aligned planes. Since the sqrt receives a sum of squares as in argument there’s no risk of getting a NaN. There’s a clear drawback that this formula is O(N^2) vs Heron’s/Kahan’s O(N).