# Spherical coordinates

The default radius of the sphere is set to one.

``double Radius = 1.;``

Rather than the classical mathematical convention, we use the geographic convention: x is the longitude within $\left[-180:180\right]$ degrees and y is the latitude within $\left[-90:90\right]$ degrees. $\Delta$ is the characteristic cell length expressed in the same units as Radius.

``````map {
x = x < -180. ? x + 360. : x > 180. ? x - 360. : x;
Delta_x = Delta_y = Δ;
}``````

On trees we need to define consistent refinement functions. The default restriction functions (averages) are already consistent (i.e. they preserve volume and length integrals).

Mathematically we have $cm=fm.y=\mathrm{cos}\left(\varphi \right)$ $fm.x=1$ with $\varphi$ the latitude in radians. This is the definition we use for the length scale factors fm.

For the volume scale factor cm, more care needs to be taken to guarantee discrete volume conservation i.e. $\sum _{children}c{m}_{child}=4c{m}_{parent}$ This is achieved by computing the exact area of a cell i.e. ${\Delta }^{2}{\int }_{\varphi -d\varphi /2}^{\varphi +d\varphi /2}\mathrm{cos}\varphi d\varphi ={\Delta }^{2}\left(\mathrm{sin}\left(\varphi +d\varphi /2\right)-\mathrm{sin}\left(\varphi -d\varphi /2\right)\right)={\Delta }^{2}cm$

``````#if TREE
static void refine_cm_lonlat (Point point, scalar cm)
{
double φ = y*π/180., dphi = Δ/(2.*Radius);
double sinphi = sin(φ);
fine(cm,0,0) = fine(cm,1,0) = (sinphi - sin(φ - dphi))/dphi;
fine(cm,0,1) = fine(cm,1,1) = (sin(φ + dphi) - sinphi)/dphi;
}

static void refine_face_x_lonlat (Point point, scalar fm)
{
if (!is_refined(neighbor(-1)))
fine(fm,0,0) = fine(fm,0,1) = 1.;
if (!is_refined(neighbor(1)) && neighbor(1).neighbors)
fine(fm,2,0) = fine(fm,2,1) = 1.;
fine(fm,1,0) = fine(fm,1,1) = 1.;
}

static void refine_face_y_lonlat (Point point, scalar fm)
{
double φ = y*π/180., dphi = Δ/(2.*Radius);
if (!is_refined(neighbor(0,-1)))
fine(fm,0,0) = fine(fm,1,0) = cos(φ - dphi);
if (!is_refined(neighbor(0,1)) && neighbor(0,1).neighbors)
fine(fm,0,2) = fine(fm,1,2) = cos(φ + dphi);
fine(fm,0,1) = fine(fm,1,1) = cos(φ);
}
#endif

event metric (i = 0) {``````

We initialise the scale factors, taking care to first allocate the fields if they are still constant.

``````  if (is_constant(cm)) {
scalar * l = list_copy (all);
cm = new scalar;
free (all);
all = list_concat ({cm}, l);
free (l);
}
scalar cmv = cm;
foreach() {
double φ = y*π/180., dphi = Δ/(2.*Radius);
cmv[] = (sin(φ + dphi) - sin(φ - dphi))/(2.*dphi);
}

if (is_constant(fm.x)) {
scalar * l = list_copy (all);
fm = new face vector;
free (all);
all = list_concat ((scalar *){fm}, l);
free (l);
}
face vector fmv = fm;
foreach_face(x)
fmv.x[] = 1.;
foreach_face(y)
fmv.y[] = cos(y*π/180.);``````

We set our refinement/prolongation functions on trees.

``````#if TREE
cm.refine = cm.prolongation = refine_cm_lonlat;
fm.x.prolongation = refine_face_x_lonlat;
fm.y.prolongation = refine_face_y_lonlat;
#endif

boundary ({cm, fm});
}``````