[SCM] Gerris Flow Solver branch, upstream, updated. e8f73a07832050124d2b8bf6c6f35b33180e65a8
Stephane Popinet
popinet at users.sf.net
Tue Nov 24 12:25:13 UTC 2009
The following commit has been merged in the upstream branch:
commit 4fd1ca04bdc0db3f6751c82b567c8c287e496891
Author: Stephane Popinet <popinet at users.sf.net>
Date: Tue Sep 29 12:25:59 2009 +1000
General Orthogonal Coordinates
darcs-hash:20090929022559-d4795-d08e3becf85fc7fcf7f29be0f0241bb1ad7f8719.gz
diff --git a/doc/figures/lonlat.tm b/doc/figures/lonlat.tm
index d476007..a4d290e 100644
--- a/doc/figures/lonlat.tm
+++ b/doc/figures/lonlat.tm
@@ -232,21 +232,21 @@
\;
<\equation*>
- m<rsub|x>*m<rsub|y>*d*x*d*y*\<partial\><rsub|t>\<b-U\>+<big|sum><rsub|f<rsub|x>>m<rsub|y>*d*y*\<b-F\><rsub|x>(\<b-U\>)+<big|sum><rsub|f<rsub|y>>m<rsub|x>*d*x*\<b-F\><rsub|y>(\<b-U\>)=0
+ h<rsub|\<lambda\>>*h<rsub|\<theta\>>*d\<lambda\>*d\<theta\>*\<partial\><rsub|t>\<b-U\>+<big|sum><rsub|f<rsub|\<lambda\>>>h<rsub|\<theta\>>*d\<theta\>*\<b-F\><rsub|\<lambda\>>(\<b-U\>)+<big|sum><rsub|f<rsub|\<theta\>>>h<rsub|\<lambda\>>*d\<lambda\>*\<b-F\><rsub|\<theta\>>(\<b-U\>)=0
</equation*>
can be rewritten
<\equation*>
- m<rsub|x>*m<rsub|y>*\<partial\><rsub|t>\<b-U\>+<frac|1|d*x>*<big|sum><rsub|f<rsub|x>>m<rsub|y>*\<b-F\><rsub|x>(\<b-U\>)+<frac|1|d*y>*<big|sum><rsub|f<rsub|y>>m<rsub|x>*\<b-F\><rsub|y>(\<b-U\>)=0.
+ h<rsub|\<lambda\>>*h<rsub|\<theta\>>*\<partial\><rsub|t>\<b-U\>+<frac|1|d\<lambda\>>*<big|sum><rsub|f<rsub|\<lambda\>>>h<rsub|\<theta\>>*\<b-F\><rsub|\<lambda\>>(\<b-U\>)+<frac|1|d\<theta\>>*<big|sum><rsub|f<rsub|\<theta\>>>h<rsub|\<lambda\>>*\<b-F\><rsub|\<theta\>>(\<b-U\>)=0.
</equation*>
- The differential form can be recovered by making <math|d*x> and <math|d*y>
- tend to zero
+ The differential form can be recovered by making <math|d\<lambda\>> and
+ <math|d\<theta\>> tend to zero
<\with|mode|math>
<\eqnarray*>
- <tformat|<table|<row|<cell|>|<cell|m<rsub|x>*m<rsub|y>*\<partial\><rsub|t>\<phi\>+\<partial\><rsub|x>(m<rsub|y>*\<phi\>*u)+\<partial\><rsub|y>(m<rsub|x>*\<phi\>*v)=0,>|<cell|>>|<row|<cell|>|<cell|m<rsub|x>*m<rsub|y>*\<partial\><rsub|t>(\<phi\>*u)+\<partial\><rsub|x><left|[>m<rsub|y>*<left|(>\<phi\>*u<rsup|2>+g*<frac|\<phi\><rsup|2>|2><right|)><right|]>+\<partial\><rsub|y><left|(>m<rsub|x>*\<phi\>*u*v<right|)>=0,>|<cell|>>|<row|<cell|>|<cell|m<rsub|x>*m<rsub|y>*\<partial\><rsub|t>(\<phi\>*v)+\<partial\><rsub|x>(m<rsub|y>*\<phi\>*u*v)+\<partial\><rsub|y><left|[>m<rsub|x>*<left|(>\<phi\>*v<rsup|2>+g*<frac|\<phi\><rsup|2>|2><right|)><right|]>=0.>|<cell|>>>>
+ <tformat|<table|<row|<cell|>|<cell|h<rsub|\<lambda\>>*h<rsub|\<theta\>>*\<partial\><rsub|t>\<phi\>+\<partial\><rsub|\<lambda\>>(h<rsub|\<theta\>>*\<phi\>*u)+\<partial\><rsub|\<theta\>>(h<rsub|\<lambda\>>*\<phi\>*v)=0,>|<cell|>>|<row|<cell|>|<cell|h<rsub|\<lambda\>>*h<rsub|\<theta\>>*\<partial\><rsub|t>(\<phi\>*u)+\<partial\><rsub|\<lambda\>><left|[>h<rsub|\<theta\>>*<left|(>\<phi\>*u<rsup|2>+g*<frac|\<phi\><rsup|2>|2><right|)><right|]>+\<partial\><rsub|\<theta\>><left|(>h<rsub|\<lambda\>>*\<phi\>*u*v<right|)>=0,>|<cell|>>|<row|<cell|>|<cell|h<rsub|\<lambda\>>*h<rsub|\<theta\>>*\<partial\><rsub|t>(\<phi\>*v)+\<partial\><rsub|\<lambda\>>(h<rsub|\<theta\>>*\<phi\>*u*v)+\<partial\><rsub|\<theta\>><left|[>h<rsub|\<lambda\>>*<left|(>\<phi\>*v<rsup|2>+g*<frac|\<phi\><rsup|2>|2><right|)><right|]>=0.>|<cell|>>>>
</eqnarray*>
</with>
@@ -254,21 +254,21 @@
<\with|mode|math>
<\eqnarray*>
- <tformat|<table|<row|<cell|>|<cell|\<partial\><rsub|t>\<phi\>+<frac|u*|m<rsub|x>>*\<partial\><rsub|x>\<phi\>*+<frac|v*|m<rsub|y>>*\<partial\><rsub|y>\<phi\>+<frac|\<phi\>*|m<rsub|x>*m<rsub|y>>*<left|[>\<partial\><rsub|x>(m<rsub|y>*u)+\<partial\><rsub|y>(m<rsub|x>*v)<right|]>=0,>|<cell|>>|<row|<cell|>|<cell|\<partial\><rsub|t>u*+<frac|u*|m<rsub|x>>*\<partial\><rsub|x>u+<frac|v*|m<rsub|y>>*\<partial\><rsub|y>u*+<frac|g|m<rsub|x>>*\<partial\><rsub|x>\<phi\>+g*<frac|\<phi\>|2*m<rsub|x>*m<rsub|y>>*\<partial\><rsub|x>m<rsub|y>=0,>|<cell|>>|<row|<cell|>|<cell|\<partial\><rsub|t>v+<frac|u|m<rsub|x>>*\<partial\><rsub|x>v+<frac|v|m<rsub|y>>*\<partial\><rsub|y>v+<frac|g|m<rsub|y>>*\<partial\><rsub|y>\<phi\>+g*<frac|\<phi\>|2*m<rsub|x>*m<rsub|y>>*\<partial\><rsub|y>m<rsub|x>=0.>|<cell|>>>>
+ <tformat|<table|<row|<cell|>|<cell|\<partial\><rsub|t>\<phi\>+<frac|u*|h<rsub|\<lambda\>>>*\<partial\><rsub|\<lambda\>>\<phi\>*+<frac|v*|h<rsub|\<theta\>>>*\<partial\><rsub|\<theta\>>\<phi\>+<frac|\<phi\>*|h<rsub|\<lambda\>>*h<rsub|\<theta\>>>*<left|[>\<partial\><rsub|\<lambda\>>(h<rsub|\<theta\>>*u)+\<partial\><rsub|\<theta\>>(h<rsub|\<lambda\>>*v)<right|]>=0,>|<cell|>>|<row|<cell|>|<cell|\<partial\><rsub|t>u*+<frac|u*|h<rsub|\<lambda\>>>*\<partial\><rsub|\<lambda\>>u+<frac|v*|h<rsub|\<theta\>>>*\<partial\><rsub|\<theta\>>u*+<frac|g|h<rsub|\<lambda\>>>*\<partial\><rsub|\<lambda\>>\<phi\>+g*<frac|\<phi\>|2*h<rsub|\<lambda\>>*h<rsub|\<theta\>>>*\<partial\><rsub|\<lambda\>>h<rsub|\<theta\>>=0,>|<cell|>>|<row|<cell|>|<cell|\<partial\><rsub|t>v+<frac|u|h<rsub|\<lambda\>>>*\<partial\><rsub|\<lambda\>>v+<frac|v|h<rsub|\<theta\>>>*\<partial\><rsub|\<theta\>>v+<frac|g|h<rsub|\<theta\>>>*\<partial\><rsub|\<theta\>>\<phi\>+g*<frac|\<phi\>|2*h<rsub|\<lambda\>>*h<rsub|\<theta\>>>*\<partial\><rsub|\<theta\>>h<rsub|\<lambda\>>=0.>|<cell|>>>>
</eqnarray*>
</with>
Introducing the notation
<\equation*>
- d<rsub|t>=\<partial\><rsub|t>+<frac|u*|m<rsub|x>>*\<partial\><rsub|x>+<frac|v*|m<rsub|y>>*\<partial\><rsub|y>,
+ d<rsub|t>=\<partial\><rsub|t>+<frac|u*|h<rsub|\<lambda\>>>*\<partial\><rsub|\<lambda\>>+<frac|v*|h<rsub|\<theta\>>>*\<partial\><rsub|\<theta\>>,
</equation*>
this can be rewritten
<\with|mode|math>
<\eqnarray*>
- <tformat|<table|<row|<cell|>|<cell|d<rsub|t>\<phi\>+<frac|\<phi\>*|m<rsub|x>*m<rsub|y>>*<left|[>\<partial\><rsub|x>(m<rsub|y>*u)+\<partial\><rsub|y>(m<rsub|x>*v)<right|]>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>u*<with|color|blue|<with|color|green|<with|color|green|-f<rsub|G>*v>>>+<frac|g|m<rsub|x>>*\<partial\><rsub|x>\<phi\><with|color|red|+<frac|g*\<phi\>|2*m<rsub|x>*m<rsub|y>>*\<partial\><rsub|x>m<rsub|y>>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>v<with|color|green|<with|color|blue|<with|color|green|+f<rsub|G>*u>>>+<frac|g|m<rsub|y>>*\<partial\><rsub|y>\<phi\><with|color|red|+*<frac|g*\<phi\>|2*m<rsub|x>*m<rsub|y>>*\<partial\><rsub|y>m<rsub|x>>=0.>|<cell|>>>>
+ <tformat|<table|<row|<cell|>|<cell|d<rsub|t>\<phi\>+<frac|\<phi\>*|h<rsub|\<lambda\>>*h<rsub|\<theta\>>>*<left|[>\<partial\><rsub|\<lambda\>>(h<rsub|\<theta\>>*u)+\<partial\><rsub|\<theta\>>(h<rsub|\<lambda\>>*v)<right|]>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>u*<with|color|blue|<with|color|green|<with|color|green|-f<rsub|G>*v>>>+<frac|g|h<rsub|\<lambda\>>>*\<partial\><rsub|\<lambda\>>\<phi\><with|color|red|+<frac|g*\<phi\>|2*h<rsub|\<lambda\>>*h<rsub|\<theta\>>>*\<partial\><rsub|\<lambda\>>h<rsub|\<theta\>>>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>v<with|color|green|<with|color|blue|<with|color|green|+f<rsub|G>*u>>>+<frac|g|h<rsub|\<theta\>>>*\<partial\><rsub|\<theta\>>\<phi\><with|color|red|+*<frac|g*\<phi\>|2*h<rsub|\<lambda\>>*h<rsub|\<theta\>>>*\<partial\><rsub|\<theta\>>h<rsub|\<lambda\>>>=0.>|<cell|>>>>
</eqnarray*>
</with>
@@ -276,7 +276,7 @@
of Williamson et al, 1991) are indicated in red and green respectively and
<\equation*>
- f<rsub|G>\<equiv\><frac|v*\<partial\><rsub|x>m<rsub|y>-u*\<partial\><rsub|y>m<rsub|x>|m<rsub|x>*m<rsub|y>>*.*
+ f<rsub|G>\<equiv\><frac|v*\<partial\><rsub|\<lambda\>>h<rsub|\<theta\>>-u*\<partial\><rsub|\<theta\>>h<rsub|\<lambda\>>|h<rsub|\<lambda\>>*h<rsub|\<theta\>>>*.*
</equation*>
<subsubsection|Application to spherical coordinates>
@@ -284,8 +284,8 @@
For spherical coordinates
<\equation*>
- m<rsub|x>=r*cos y,<hspace|10pt>m<rsub|y>=r,<hspace|10pt>\<partial\><rsub|x>m<rsub|y>=0,<hspace|10pt>\<partial\><rsub|y>m<rsub|x>=-r*sin
- y,
+ h<rsub|\<lambda\>>=r*cos \<theta\>,<hspace|10pt>h<rsub|\<theta\>>=r,<hspace|10pt>\<partial\><rsub|\<lambda\>>h<rsub|\<theta\>>=0,<hspace|10pt>\<partial\><rsub|\<theta\>>h<rsub|\<lambda\>>=-r*sin
+ \<theta\>,
</equation*>
which gives
@@ -293,11 +293,11 @@
<\with|mode|math>
<\eqnarray*>
<tformat|<table|<row|<cell|>|<cell|d<rsub|t>\<phi\>+<frac|\<phi\>*|r*cos
- y>*<left|[>\<partial\><rsub|x>u+\<partial\><rsub|y>(v*cos
- y)<right|]>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>u*<with|color|blue|<with|color|green|<with|color|green|-<frac|u*v|r>*tan
- y>>>+<frac|g|r*cos y>*\<partial\><rsub|x>\<phi\>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>v<with|color|green|<with|color|blue|<with|color|green|<with|color|green|+<frac|u<rsup|2>|r>*tan
- y>>>>+<frac|g|m<rsub|y>>*\<partial\><rsub|y>\<phi\><with|color|red|-<frac|g*\<phi\>|2*r>*tan
- y>=0.>|<cell|>>>>
+ \<theta\>>*<left|[>\<partial\><rsub|\<lambda\>>u+\<partial\><rsub|\<theta\>>(v*cos
+ \<theta\>)<right|]>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>u*<with|color|blue|<with|color|green|<with|color|green|-<frac|u*v|r>*tan
+ \<theta\>>>>+<frac|g|r*cos \<theta\>>*\<partial\><rsub|\<lambda\>>\<phi\>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>v<with|color|green|<with|color|blue|<with|color|green|<with|color|green|+<frac|u<rsup|2>|r>*tan
+ \<theta\>>>>>+<frac|g|r>*\<partial\><rsub|\<theta\>>\<phi\><with|color|red|-<frac|g*\<phi\>|2*r>*tan
+ \<theta\>>=0.>|<cell|>>>>
</eqnarray*>
</with>
@@ -306,14 +306,14 @@
For polar coordinates
<\equation*>
- m<rsub|x>=1,<hspace|10pt>m<rsub|y>=x,<hspace|10pt>\<partial\><rsub|x>m<rsub|y>=1,<hspace|10pt>\<partial\><rsub|y>m<rsub|x>=0,
+ h<rsub|\<lambda\>>=1,<hspace|10pt>h<rsub|\<theta\>>=\<lambda\>,<hspace|10pt>\<partial\><rsub|\<lambda\>>h<rsub|\<theta\>>=1,<hspace|10pt>\<partial\><rsub|\<theta\>>h<rsub|\<lambda\>>=0,
</equation*>
which gives
<\with|mode|math>
<\eqnarray*>
- <tformat|<table|<row|<cell|>|<cell|d<rsub|t>\<phi\>+<frac|\<phi\>*|x>*<left|[>\<partial\><rsub|x>(x*u)+\<partial\><rsub|y>v<right|]>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>u*<with|color|blue|<with|color|green|<with|color|green|-<frac|v<rsup|2>|x>>>>+g*\<partial\><rsub|x>\<phi\><with|color|red|+g*<frac|\<phi\>|2*x>>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>v<with|color|green|<with|color|blue|<with|color|green|<with|color|green|+<frac|u*v|x>>>>>+<frac|g|x>*\<partial\><rsub|y>\<phi\>=0.>|<cell|>>>>
+ <tformat|<table|<row|<cell|>|<cell|d<rsub|t>\<phi\>+<frac|\<phi\>*|\<lambda\>>*<left|[>\<partial\><rsub|\<lambda\>>(\<lambda\>*u)+\<partial\><rsub|\<theta\>>v<right|]>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>u*<with|color|blue|<with|color|green|<with|color|green|-<frac|v<rsup|2>|\<lambda\>>>>>+g*\<partial\><rsub|\<lambda\>>\<phi\><with|color|red|+g*<frac|\<phi\>|2*\<lambda\>>>=0,>|<cell|>>|<row|<cell|>|<cell|d<rsub|t>v<with|color|green|<with|color|blue|<with|color|green|<with|color|green|+<frac|u*v|\<lambda\>>>>>>+<frac|g|\<lambda\>>*\<partial\><rsub|\<theta\>>\<phi\>=0.>|<cell|>>>>
</eqnarray*>
</with>
</body>
@@ -350,12 +350,16 @@
<no-break><pageref|auto-3>>
<with|par-left|<quote|1.5fn>|3<space|2spc>The Saint-Venant equations in
- orthogonal curvilinear coordinates <datoms|<macro|x|<repeat|<arg|x>|<with|font-series|medium|<with|font-size|1|<space|0.2fn>.<space|0.2fn>>>>>|<htab|5mm>>
+ general orthogonal coordinates <datoms|<macro|x|<repeat|<arg|x>|<with|font-series|medium|<with|font-size|1|<space|0.2fn>.<space|0.2fn>>>>>|<htab|5mm>>
<no-break><pageref|auto-4>>
<with|par-left|<quote|3fn>|3.1<space|2spc>Application to spherical
coordinates <datoms|<macro|x|<repeat|<arg|x>|<with|font-series|medium|<with|font-size|1|<space|0.2fn>.<space|0.2fn>>>>>|<htab|5mm>>
<no-break><pageref|auto-5>>
+
+ <with|par-left|<quote|3fn>|3.2<space|2spc>Application to polar
+ coordinates <datoms|<macro|x|<repeat|<arg|x>|<with|font-series|medium|<with|font-size|1|<space|0.2fn>.<space|0.2fn>>>>>|<htab|5mm>>
+ <no-break><pageref|auto-6>>
</associate>
</collection>
</auxiliary>
\ No newline at end of file
diff --git a/src/river.c b/src/river.c
index e18d571..6a1a0a7 100644
--- a/src/river.c
+++ b/src/river.c
@@ -194,23 +194,26 @@ static void sources (FttCell * cell, GfsRiver * r)
GFS_VALUE (cell, r->v[2]) += r->dt*r->g/2.*(etaL + etaR)*(zbL - zbR)/delta;
- /* longitude-latitude geometric source terms */
-#if 0
- FttVector p;
- ftt_cell_pos (cell, &p);
- gdouble radius = 3./(2.*M_PI), g = 1.;
- gdouble theta = p.y/radius, tantheta = tan (theta);
- gdouble
- phiu = GFS_VALUE (cell, r->v[1]),
- phiv = GFS_VALUE (cell, r->v[2]),
- phi = GFS_VALUE (cell, r->v[0]);
-#if 1
- GFS_VALUE (cell, r->v[1]) += r->dt*tantheta/radius*phiu*phiv/phi;
- GFS_VALUE (cell, r->v[2]) += r->dt*tantheta/radius*(- phiu*phiu/phi - g*phi*phi/2.);
-#else
- GFS_VALUE (cell, r->v[2]) += r->dt*tantheta/radius*(- g*phi*phi/2.);
-#endif
-#endif
+ /* Geometric source terms (see doc/figures/lonlat.tm) */
+ if (GFS_DOMAIN (r)->cell_metric) {
+ GfsDomain * domain = GFS_DOMAIN (r);
+ FttCellFace face = { cell };
+ gdouble f[FTT_NEIGHBORS];
+ for (face.d = 0; face.d < FTT_NEIGHBORS; face.d++)
+ f[face.d] = (* domain->face_metric) (domain, &face, domain->metric_data);
+ gdouble dh_dl = f[FTT_RIGHT] - f[FTT_LEFT];
+ gdouble dh_dt = f[FTT_TOP] - f[FTT_BOTTOM];
+ gdouble dldh = (* domain->cell_metric) (domain, cell, domain->metric_data)*
+ delta*GFS_SIMULATION (r)->physical_params.L;
+ gdouble
+ phi = GFS_VALUE (cell, r->v1[0]),
+ phiu = GFS_VALUE (cell, r->v1[1]),
+ phiv = GFS_VALUE (cell, r->v1[2]);
+ gdouble fG = phiv*dh_dl - phiu*dh_dt;
+ gdouble g = GFS_SIMULATION (r)->physical_params.g;
+ GFS_VALUE (cell, r->v[1]) += r->dt*(g*phi*phi/2.*dh_dl + fG*phiv)/dldh;
+ GFS_VALUE (cell, r->v[2]) += r->dt*(g*phi*phi/2.*dh_dt - fG*phiu)/dldh;
+ }
}
static void advance (GfsRiver * r, gdouble dt)
@@ -483,8 +486,8 @@ static void river_init (GfsRiver * r)
r->v1[2] = gfs_domain_add_variable (domain, NULL, NULL);
gfs_variable_set_vector (&r->v1[1], 2);
- r->dv[0][0] = gfs_domain_add_variable (domain, "Px", "x-component of the thickness gradient");
- r->dv[1][0] = gfs_domain_add_variable (domain, "Py", "y-component of the thickness gradien");
+ r->dv[0][0] = gfs_domain_add_variable (domain, "Px", "x-component of the depth gradient");
+ r->dv[1][0] = gfs_domain_add_variable (domain, "Py", "y-component of the depth gradien");
r->dv[0][1] = gfs_domain_add_variable (domain, "Ux", "x-component of the flux gradient");
r->dv[1][1] = gfs_domain_add_variable (domain, "Uy", "y-component of the flux gradient");
r->dv[0][2] = gfs_domain_add_variable (domain, "Vx", "x-component of the flux gradient");
diff --git a/test/Makefile.am b/test/Makefile.am
index 3b7dbce..e5d94c1 100644
--- a/test/Makefile.am
+++ b/test/Makefile.am
@@ -27,7 +27,8 @@ TESTDIRS = \
geo \
waves \
nz \
- parabola
+ parabola \
+ lonlat
EXTRA_DIST = \
template.tex \
diff --git a/test/lonlat/coriolis/coriolis.gfs b/test/lonlat/coriolis/coriolis.gfs
new file mode 100644
index 0000000..c101df2
--- /dev/null
+++ b/test/lonlat/coriolis/coriolis.gfs
@@ -0,0 +1,57 @@
+# Title: Circular dam break on a rotating sphere
+#
+# Description:
+#
+# Similar test case but with rotation. See also test case of \cite{rossmanith2004}, Figure 7.
+#
+# \begin{figure}[htbp]
+# \caption{Solution to the rotating shallow-water equations computed
+# on a longitude-latitude grid in the domain
+# $[-75^\circ,75^\circ]\times[-75^\circ,75^\circ]$ with $256\times
+# 256$ points. The Coriolis parameter is set to $f=10$. The solution
+# is shown at times (a) $t=0.4$, (b) $t=0.8$, and (c) $t=1.2$.
+# }
+# \begin{center}
+# \begin{tabular}{cc}
+# (a) \includegraphics[width=0.45\hsize]{isolines-0.4.eps} &
+# (b) \includegraphics[width=0.45\hsize]{isolines-0.8.eps} \\
+# \multicolumn{2}{c}{(c) \includegraphics[width=0.45\hsize]{isolines-1.2.eps}}
+# \end{tabular}
+# \end{center}
+# \end{figure}
+#
+# Author: St\'ephane Popinet
+# Command: gerris2D -m coriolis.gfs
+# Version: 090924
+# Required files: isolines.gfv
+# Running time: 5 minutes
+# Generated files: isolines-0.4.eps isolines-0.8.eps isolines-1.2.eps
+#
+Define LENGTH (150./180.*M_PI)
+
+1 0 GfsRiver GfsBox GfsGEdge {} {
+ PhysicalParams { L = LENGTH }
+ MetricLonLat 1.
+ AdvectionParams { cfl = 0.4 }
+ Refine 8
+ InitFraction P (0.2 - acos(cos(x)*cos(y)))
+ Init {} { P = 0.2 + 1.8*P/LENGTH }
+ Time { end = 1.4 }
+ SourceCoriolis 10.*sin(y)
+ OutputTime { istep = 10 } stderr
+ OutputSimulation { step = 0.4 } sim-%g.gfs
+ OutputSimulation { step = 0.4 } sim-%g.txt { variables = U,V,P format = text }
+# OutputSimulation { istep = 10 } stdout
+ EventScript { start = end } {
+ for i in 0.4 0.8 1.2; do
+ echo "Save stdout { format = EPS line_width = 0.5 }" | \
+ gfsview-batch2D sim-$i.gfs isolines.gfv > isolines-$i.eps
+ done
+ }
+}
+GfsBox {
+ right = Boundary { BcNeumann U 0 }
+ left = Boundary { BcNeumann U 0 }
+ top = Boundary { BcNeumann V 0 }
+ bottom = Boundary { BcNeumann V 0 }
+}
diff --git a/doc/examples/hump/isolines.gfv b/test/lonlat/coriolis/isolines.gfv
similarity index 86%
copy from doc/examples/hump/isolines.gfv
copy to test/lonlat/coriolis/isolines.gfv
index 7708f1e..2da075e 100644
--- a/doc/examples/hump/isolines.gfv
+++ b/test/lonlat/coriolis/isolines.gfv
@@ -1,9 +1,9 @@
# GfsView 2D
View {
- tx = -1.10977 ty = -0.473248
+ tx = 0 ty = 0
sx = 1 sy = 1 sz = 1
q0 = 0 q1 = 0 q2 = 0 q3 = 1
- fov = 17.1233
+ fov = 25.764
r = 0.3 g = 0.4 b = 0.6
res = 1
lc = 0.001
@@ -16,7 +16,7 @@ Isoline {
} {
n.x = 0 n.y = 0 n.z = 1
pos = 0
-} P+Zb {
+} P {
amin = 1
amax = 1
cmap = Jet
@@ -24,7 +24,7 @@ Isoline {
reversed = 0
use_scalar = 1
} {
- n = 30
+ n = 20
}
Boundaries {
r = 0 g = 0 b = 0
diff --git a/doc/examples/hump/isolines.gfv b/test/lonlat/isolines.gfv
similarity index 86%
copy from doc/examples/hump/isolines.gfv
copy to test/lonlat/isolines.gfv
index 7708f1e..2da075e 100644
--- a/doc/examples/hump/isolines.gfv
+++ b/test/lonlat/isolines.gfv
@@ -1,9 +1,9 @@
# GfsView 2D
View {
- tx = -1.10977 ty = -0.473248
+ tx = 0 ty = 0
sx = 1 sy = 1 sz = 1
q0 = 0 q1 = 0 q2 = 0 q3 = 1
- fov = 17.1233
+ fov = 25.764
r = 0.3 g = 0.4 b = 0.6
res = 1
lc = 0.001
@@ -16,7 +16,7 @@ Isoline {
} {
n.x = 0 n.y = 0 n.z = 1
pos = 0
-} P+Zb {
+} P {
amin = 1
amax = 1
cmap = Jet
@@ -24,7 +24,7 @@ Isoline {
reversed = 0
use_scalar = 1
} {
- n = 30
+ n = 20
}
Boundaries {
r = 0 g = 0 b = 0
diff --git a/test/lonlat/lonlat.gfs b/test/lonlat/lonlat.gfs
new file mode 100644
index 0000000..920e726
--- /dev/null
+++ b/test/lonlat/lonlat.gfs
@@ -0,0 +1,117 @@
+# Title: Circular dam break on a sphere
+#
+# Description:
+#
+# An initial circular cylinder collapses and creates shock and
+# rarefaction waves. The initial condition are radially-symmetric and
+# should remain so. The problem is discretised using
+# longitude-latitude spherical coordinates. Deviations from radial
+# symmetry are a measure of the accuracy of treatment of geometric
+# source terms.
+#
+# This test case was proposed by \cite{rossmanith2004}, Figures 5 and 6.
+#
+# \begin{figure}[htbp]
+# \caption{\label{height}Solution to the shallow-water equations computed on a
+# longitude-latitude grid in the domain
+# $[-75^\circ,75^\circ]\times[-75^\circ,75^\circ]$ with $256\times
+# 256$ points. The solution is shown at times (a) $t=0.3$, (b)
+# $t=0.6$, and (c) $t=0.9$. The contours do not appear circular
+# because the solution has been projected down to a plane.}
+# \begin{center}
+# \begin{tabular}{cc}
+# (a) \includegraphics[width=0.45\hsize]{isolines-0.3.eps} &
+# (b) \includegraphics[width=0.45\hsize]{isolines-0.6.eps} \\
+# \multicolumn{2}{c}{(c) \includegraphics[width=0.45\hsize]{isolines-0.9.eps}}
+# \end{tabular}
+# \end{center}
+# \end{figure}
+#
+# \begin{figure}[htbp]
+# \caption{Scatter plot of the (radial) solution shown in Figure
+# \ref{height}. The green line is the average solution. The solution
+# is shown at times (a) $t=0.3$, (b) $t=0.6$, and (c) $t=0.9$.}
+# \begin{center}
+# \begin{tabular}{c}
+# (a) \includegraphics[width=0.8\hsize]{p-0.3.eps} \\
+# (b) \includegraphics[width=0.8\hsize]{p-0.6.eps} \\
+# (c) \includegraphics[width=0.8\hsize]{p-0.9.eps}
+# \end{tabular}
+# \end{center}
+# \end{figure}
+#
+# Author: St\'ephane Popinet
+# Command: gerris2D -m lonlat.gfs
+# Version: 090924
+# Required files: isolines.gfv
+# Running time: 5 minutes
+# Generated files: isolines-0.3.eps isolines-0.6.eps isolines-0.9.eps p-0.3.eps p-0.6.eps p-0.9.eps
+#
+Define LENGTH (150./180.*M_PI)
+
+1 0 GfsRiver GfsBox GfsGEdge {} {
+ PhysicalParams { L = LENGTH }
+ MetricLonLat 1.
+ AdvectionParams { cfl = 0.4 }
+ Refine 8
+ InitFraction P (0.2 - acos(cos(x)*cos(y)))
+ Init {} { P = 0.2 + 1.8*P/LENGTH }
+ Time { end = 0.9 }
+ OutputTime { istep = 10 } stderr
+ OutputSimulation { step = 0.3 } sim-%g.gfs
+ OutputSimulation { step = 0.3 } sim-%g.txt { variables = U,V,P format = text }
+# OutputSimulation { istep = 10 } stdout
+ EventScript { start = end } {
+ for i in 0.3 0.6 0.9; do
+ if awk 'BEGIN { n1 = 0 }{
+ c = cos($1)*cos($2);
+ d = atan2(sqrt(1. - c*c),c)*180./3.14159265359
+ i = int(d*2.)
+ x[i] += d
+ y[i] += $6
+ n[i]++
+ x1[n1] = d;
+ y1[n1++] = $6;
+ }END {
+ for (i = 0; i <= 180; i++)
+ if (n[i] > 0)
+ print x[i]/n[i], y[i]/n[i], n[i];
+ sum = 0.;
+ for (i = 0; i < n1; i++) {
+ j = int(x1[i]*2.)
+ d = y1[i] - y[j]/n[j];
+ sum += d*d;
+ }
+ scatter = sqrt(sum/n1);
+ if (scatter > 0.012) {
+ print scatter > "/dev/stderr";
+ exit (1);
+ }
+ }' < sim-$i.txt > prof-$i.txt ; then :
+ else
+ exit $GFS_STOP;
+ fi
+
+ cat <<EOF | gnuplot
+rdist(x,y)=acos(cos(x)*cos(y))*180./pi
+cdist(x,y)=sqrt(x*x+y*y)*180./pi
+set xlabel 'Angular distance (degree)'
+set ylabel 'Surface height'
+set xtics 0,22.5,90
+set ytics 0,0.25,0.75
+set term postscript eps color lw 2 solid 20
+set output 'p-$i.eps'
+plot [0:90][0:0.75]'sim-$i.txt' u (rdist(\$1,\$2)):6 ps 0.25 pt 6 t '', 'prof-$i.txt' w l lw 2 t ''
+EOF
+
+ echo "Save stdout { format = EPS line_width = 0.5 }" | \
+ gfsview-batch2D sim-$i.gfs isolines.gfv > isolines-$i.eps
+ done
+ }
+}
+GfsBox {
+ right = Boundary { BcNeumann U 0 }
+ left = Boundary { BcNeumann U 0 }
+ top = Boundary { BcNeumann V 0 }
+ bottom = Boundary { BcNeumann V 0 }
+}
diff --git a/test/template.tex b/test/template.tex
index d447a6d..9b1e119 100644
--- a/test/template.tex
+++ b/test/template.tex
@@ -108,6 +108,11 @@ branch only.
\input{parabola/parabola.tex}
+\section{General Orthogonal Coordinates}
+
+\input{lonlat/lonlat.tex}
+\input{lonlat/coriolis/coriolis.tex}
+
\bibliographystyle{plain}
\bibliography{tests}
diff --git a/test/tests.bib b/test/tests.bib
index c29b3d0..26dcf5f 100644
--- a/test/tests.bib
+++ b/test/tests.bib
@@ -197,6 +197,15 @@
pages= {81-85}
}
+ at Article{rossmanith2004,
+ author = {J. A, Rossmanith and D. S. Bale and R. J. LeVeque},
+ title = {A wave propagation algorithm for hyperbolic systems on curved manifolds},
+ journal = {Journal of Computational Physics},
+ year = 2004,
+ volume = 199,
+ pages = {631-662}
+}
+
@article{rudman97,
title={Volume-tracking methods for interfacial flow calculations},
author={Rudman, M.},
--
Gerris Flow Solver
More information about the debian-science-commits
mailing list