This example illustrates how to compute a gradient of a second order difference and how to vizualize the result using an export to Asymptote. If we are given three points on the manifold , the 2sphere, the second order difference can be interpreted as (twice) the distance of the midpoint to , since in Euclidean space . The first figure illustrates this situation on the sphere for points on the equator and near the north pole, to be precise
1
2
3
4
5
M = Sn(2);
x = [1;0;0];
z = [0;1;0];
y = M.geopoint([0;0;1], M.geopoint(x,z,0.5),0.05);
c = M.midPoint(x,z);
The gradient with respect to of this second order absolute difference is derived in [1] and implemented within MVIRT in gradSecondOrderMidpointModel
or can be computed using DxGeo
and DyGeo
of the manifold
class as
1
2
3
4
5
% we take the inner differential (for x and z)
inner = M.log(c,y)./M.dist(c,y);
gradX = M.AdjDxGeo(x,z,0.5,inner);
gradY =  M.log(y,c)./M.dist(c,y);
gradZ = M.AdjDyGeo(x,z,0.5,inner);
To illustrate the vectors and the points in Matlab, we can use quiver3
and plotS2
together as
1
2
3
4
5
6
7
8
9
10
11
12
pts = cat(2,x,y,z);
vecs = cat(2,gradX,gradY,gradZ);
quiver3(pts(1,:),pts(2,:),pts(3,:),vecs(1,:),vecs(2,:),vecs(3,:));
vecs1 = permute(cat(3,pts,vecs),[1,3,2]);
hold on
plotS2(cat(2,x,y,z,c),'o','MarkerSize',3);
pts1 = cat(2,x,y,z,c);
plotS2(M.geodesic(x,z,'t',0:1/100:1),'');
geo1 = M.geodesic(x,z,'t',0:1/100:1);
plotS2(M.geodesic(y,c,'t',0:1/100:1),'b');
geo2 = M.geodesic(y,c,'t',0:1/100:1);
hold off
plotting the decent direction tangent vectors , all 4 points, and two geodesics. The corresponding image in Asymptote is shown in the second Figure
We perform a step towards the negative gradient direction with step size , i.,e.
1
2
3
xN = M.exp(x,gradX);
yN = M.exp(y,gradY);
zN = M.exp(z,gradZ);
yields finally the third figure
and a reduction of the second order absolute difference from about (old dark green segment) to (light green segment).
Creating the additional arrays of points
1
2
3
pts2 = cat(2,xN,yN,zN,M.midPoint(xN,zN));
geo3 = M.geodesic(xN,zN,'t',0:1/100:1);
geo4 = M.geodesic(yN,M.midPoint(xN,zN),'t',0:1/100:1);
we can use the exportSpherePCT2Asy
export
to export points, curves (connected points) and tangent vectors to Asymptote with
1
2
3
4
5
6
7
exportSpherePCT2Asy({pts1,pts2},{geo1,geo2,geo3,geo4},{vecs1},
{ [0;0;.66],[.33;0;.66],...
[.6;.6;.6],[0;0.5;0],[0.5;.6;.5],[0;1;0],...
[0;0.5;1]},
'File','gradTV2onS2/grad.asy',...
'DotSize',[5,5,1,1,1,1.2,2],...
'Camera',[1,.8,.6],'ArrowHead',10);
where the second and third line specify the colors of all 7 sets.
See also
References

Bačák, M, Bergmann, R, Steidl, G and Weinmann, A (2016). A second order nonsmooth variational model for restoring manifoldvalued images. SIAM Journal on Scientific Computing. 38 A567–A597
We introduce a new nonsmooth variational model for the restoration of manifoldvalued data which includes second order differences in the regularization term. While such models were successfully applied for realvalued images, we introduce the second order difference and the corresponding variational models for manifold data, which up to now only existed for cyclic data. The approach requires a combination of techniques from numerical analysis, convex optimization and differential geometry. First, we establish a suitable definition of absolute second order differences for signals and images with values in a manifold. Employing this definition, we introduce a variational denoising model based on first and second order differences in the manifold setup. In order to minimize the corresponding functional, we develop an algorithm using an inexact cyclic proximal point algorithm. We propose an efficient strategy for the computation of the corresponding proximal mappings in symmetric spaces utilizing the machinery of Jacobi fields. For the \(n\)sphere and the manifold of symmetric positive definite matrices, we demonstrate the performance of our algorithm in practice. We prove the convergence of the proposed exact and inexact variant of the cyclic proximal point algorithm in Hadamard spaces. These results which are of interest on its own include, e.g., the manifold of symmetric positive definite matrices.@article{BBSW16, author = {Ba{\v{c}}{\'a}k, Miroslav and Bergmann, Ronny and Steidl, Gabriele and Weinmann, Andreas}, title = {A second order nonsmooth variational model for restoring manifoldvalued images}, journal = {SIAM Journal on Scientific Computing}, year = {2016}, volume = {38}, number = {1}, pages = {A567A597}, doi = {10.1137/15M101988X}, eprint = {1506.02409}, eprinttype = {arXiv} }