8.494e-02 0016 6.05e-12
1.199e+01 0026 5.72e-10
1.198e+03 0002 1.20e+03
- 7.030e+03 0001 nan
- 7.030e+03 0001 nan
- 7.030e+03 0001 nan
- 7.030e+03 0001 nan
- 7.030e+03 0001 nan
+ 7.030e+03 0001 nan
+ 7.030e+03 0001 nan
+ 7.030e+03 0001 nan
+ 7.030e+03 0001 nan
+ 7.030e+03 0001 nan
7.030e+03 0001 nan
----------------------------------------------------
-Exception on processing:
+Exception on processing:
--------------------------------------------------------
An error occurred in line <2476> of file <\step-33.cc> in function
void Step33::ConservationLaw<dim>::run() [with int dim = 2]
-The violated condition was:
+The violated condition was:
nonlin_iter <= 10
The name and call sequence of the exception was:
ExcMessage ("No convergence in nonlinear solver")
-Additional Information:
+Additional Information:
No convergence in nonlinear solver
--------------------------------------------------------
<img src="https://www.dealii.org/images/steps/developer/step-33.oscillation.gif " alt="" height="300">
-As we see, when the heavy mass of fluid hits the left bottom corner,
+As we see, when the heavy mass of fluid hits the left bottom corner,
some oscillation occurs and lead to the divergence of the iteration. A lazy solution to
this issue is add more viscosity. If we set the diffusion power $\eta = 1.5$ instead of $2.0$,
the simulation would be able to survive this crisis. Then, the result looks like this:
The heavy mass of fluid is drawn down the slope by gravity, where
it collides with the ski lodge and is flung into the air! Hopefully everyone
-escapes! And also, we can see the boundary between heavy mass and light mass blur quickly
+escapes! And also, we can see the boundary between heavy mass and light mass blur quickly
due to the artificial viscosity.
We can also visualize the evolution of the adaptively refined grid:
<h4>Cache the explicit part of residual</h4>
The residual calulated in ConservationLaw::assemble_cell_term function
-read
+read
$R_i = \left(\frac{\mathbf{w}^{k}_{n+1} - \mathbf{w}_n}{\delta t}
- , \mathbf{z}_i \right)_K +
+ , \mathbf{z}_i \right)_K +
\theta \mathbf{B}({\mathbf{w}^{k}_{n+1})(\mathbf{z}_i)_K +
- (1-\theta) \mathbf{B}({\mathbf{w}_{n}) (\mathbf{z}_i)_K $
-means that we calculate the spacial residual twice at one Newton
-iteration step: once respect to the current solution ${\mathbf{w}^{k}_{n+1}$
-and another respect to the last time step solution $\mathbf{w}_{n}$ which
-remains the same during all Newton interations through one timestep.
+ (1-\theta) \mathbf{B}({\mathbf{w}_{n}) (\mathbf{z}_i)_K $
+means that we calculate the spacial residual twice at one Newton
+iteration step: once respect to the current solution ${\mathbf{w}^{k}_{n+1}$
+and another respect to the last time step solution $\mathbf{w}_{n}$ which
+remains the same during all Newton interations through one timestep.
Cache up the explicit part of residual
$ \mathbf{B}({\mathbf{w}_{n}) (\mathbf{z}_i)_K}$
during Newton iteration will save lots of labor.
// residual, adding its negative to the right hand side vector, and adding
// its derivative with respect to the local variables to the Jacobian
// (i.e. the Newton matrix). Recall that the cell contributions to the
- // residual read
+ // residual read
// $R_i = \left(\frac{\mathbf{w}^{k}_{n+1} - \mathbf{w}_n}{\delta t} ,
- // \mathbf{z}_i \right)_K $ $ +
+ // \mathbf{z}_i \right)_K $ $ +
// \theta \mathbf{B}({\mathbf{w}^{k}_{n+1})(\mathbf{z}_i)_K $ $ +
// (1-\theta) \mathbf{B}({\mathbf{w}_{n}) (\mathbf{z}_i)_K $ where
- // $\mathbf{B}({\mathbf{w})(\mathbf{z}_i)_K =
+ // $\mathbf{B}({\mathbf{w})(\mathbf{z}_i)_K =
// - \left(\mathbf{F}(\mathbf{w}),\nabla\mathbf{z}_i\right)_K $ $
// + h^{\eta}(\nabla \mathbf{w} , \nabla \mathbf{z}_i)_K $ $
- // - (\mathbf{G}(\mathbf {w}), \mathbf{z}_i)_K $ for both
+ // - (\mathbf{G}(\mathbf {w}), \mathbf{z}_i)_K $ for both
// ${\mathbf{w} = \mathbf{w}^k_{n+1}$ and ${\mathbf{w} = \mathbf{w}_{n}}$ ,
- // $\mathbf{z}_i$ is the $i$th vector valued test function.
+ // $\mathbf{z}_i$ is the $i$th vector valued test function.
// Furthermore, the scalar product
// $\left(\mathbf{F}(\mathbf{w}), \nabla\mathbf{z}_i\right)_K$ is
// understood as $\int_K \sum_{c=1}^{\text{n\_components}}
// Next, in order to compute the cell contributions, we need to evaluate
- // $F({\mathbf w}^k_{n+1})$, $G({\mathbf w}^k_{n+1})$ and
+ // $F({\mathbf w}^k_{n+1})$, $G({\mathbf w}^k_{n+1})$ and
// $F({\mathbf w}_n)$, $G({\mathbf w}_n)$ at all quadrature
// points. To store these, we also need to allocate a bit of memory. Note
// that we compute the flux matrices and right hand sides in terms of
// single nonzero component (more on this topic can be found in the @ref
// vector_valued module). It will be represented by the variable
// <code>component_i</code> below. With this, the residual term can be
- // re-written as
+ // re-written as
// @f{eqnarray*}
- // R_i &=&
+ // R_i &=&
// \left(\frac{(\mathbf{w}_{n+1} -
// \mathbf{w}_n)_{\text{component\_i}}}{\delta
// t},(\mathbf{z}_i)_{\text{component\_i}}\right)_K \\
// x_d}\right)_K \\
// &+& \sum_{d=1}^{\text{dim}} h^{\eta} \left( \theta \frac{\partial
// \mathbf{w^k_{n+1}}_{\text{component\_i}}}{\partial x_d} + (1-\theta)
- // \frac{\partial \mathbf{w_n}_{\text{component\_i}}}{\partial x_d} ,
+ // \frac{\partial \mathbf{w_n}_{\text{component\_i}}}{\partial x_d} ,
// \frac{\partial (\mathbf{z}_i)_{\text{component\_i}}}{\partial x_d} \right)_K\\
// &-& \left( \theta\mathbf{G}({\mathbf{w}^k_n+1} )_{\text{component\_i}} +
// (1-\theta)\mathbf{G}({\mathbf{w}_n} )_{\text{component\_i}} ,
// Next, we need to define the values of the conservative variables
// ${\mathbf W}$ on this side of the face ($ {\mathbf W}^+$)
- // and on the opposite side (${\mathbf W}^-$), for both ${\mathbf W} =
- // {\mathbf W}^k_{n+1}$ and ${\mathbf W} = {\mathbf W}_n$.
+ // and on the opposite side (${\mathbf W}^-$), for both ${\mathbf W} =
+ // {\mathbf W}^k_{n+1}$ and ${\mathbf W} = {\mathbf W}_n$.
// The "this side" values can be
// computed in exactly the same way as in the previous function, but note
// that the <code>fe_v</code> variable now is of type FEFaceValues or
boundary_values[q],
Wminus[q]);
// Here we assume that boundary type, boundary normal vector and boundary data values
- // maintain the same during time advancing.
+ // maintain the same during time advancing.
EulerEquations<dim>::compute_Wminus (parameters.boundary_conditions[boundary_id].kind,
fe_v.normal_vector(q),
Wplus_old[q],
projects / dealii.git / commitdiff