fronts.solve_from_guess

fronts.solve_from_guess(D, Si, Sb, o_guess, S_guess, radial=False, max_nodes=1000, verbose=0)

Solve an instance of the general problem starting from a guess of the solution.

Given a positive function D, scalars \(S_i\), \(S_b\) and \(o_b\), and coordinate unit vector \(\mathbf{\hat{r}}\), finds a function S of r and t such that:

\[\begin{split}\begin{cases} \dfrac{\partial S}{\partial t} = \nabla\cdot\left[D\left(S\right)\dfrac{\partial S}{\partial r} \mathbf{\hat{r}}\right ] & r>r_b(t),t>0\\ S(r, 0) = S_i & r>0 \\ S(r_b(t), t) = S_b & t>0 \\ r_b(t) = o_b\sqrt t \end{cases}\end{split}\]

Alternative to the main solve function. This function requires a starting mesh and guess of the solution. It is significantly less robust than solve, and will fail to converge in many cases that the latter can easily handle (whether it converges will usually depend heavily on the problem, the starting mesh and the guess of the solution; it will raise a RuntimeError on failure). However, when it converges it is usually faster than solve, which may be an advantage for some use cases. You should nonetheless prefer solve unless you have a particular use case for which you have found this function to be better.

Possible use cases include refining a solution (note that solve can do that too), optimization runs in which known solutions make good first approximations of solutions with similar parameters and every second of computing time counts, and in the implementation of other solving algorithms. In all these cases, solve should probably be used as a fallback for when this function fails.

Parameters

• D (callable) – Twice-differentiable function that maps the range of S to positive values. It can be called as D(S) to evaluate it at S. It can also be called as D(S, n) with n equal to 1 or 2, in which case the first n derivatives of the function evaluated at the same S are included (in order) as additional return values. While mathematically a scalar function, D operates in a vectorized fashion with the same semantics when S is a numpy.ndarray.

• Si (float) – \(S_i\), the initial value of S in the domain.

• Sb (float) – \(S_b\), the value of S imposed at the boundary.

• o_guess (numpy.array_like, shape (n_guess,)) – Starting mesh in terms of the Boltzmann variable o. Must be strictly increasing. o_guess[0] is \(o_b\), which determines the behavior of the boundary. If zero, it implies that the boundary always exists at \(r=0\). It must be strictly positive if radial is not False. Be aware that a non-zero value implies a moving boundary. On the other end, o_guess[-1] must be large enough to contain the solution to the semi-infinite problem.

• S_guess (float or numpy.array_like, shape (n_guess,)) – Starting guess of S at the points in o_guess. If a single value, the guess is assumed uniform.

• radial ({False, 'cylindrical', 'spherical'}, optional) –

  Choice of coordinate unit vector \(\mathbf{\hat{r}}\). Must be one of the following:

  • False (default)

    \(\mathbf{\hat{r}}\) is any coordinate unit vector in rectangular (Cartesian) coordinates, or an axial unit vector in a cylindrical coordinate system

  • 'cylindrical'

    \(\mathbf{\hat{r}}\) is the radial unit vector in a cylindrical coordinate system

  • 'spherical'

    \(\mathbf{\hat{r}}\) is the radial unit vector in a spherical coordinate system

• max_nodes (int, optional) – Maximum allowed number of mesh nodes.

• verbose ({0, 1, 2}, optional) –

  Level of algorithm’s verbosity. Must be one of the following:

  • 0 (default) : work silently.

  • 1 : display a termination report.

  • 2 : display progress during iterations.

Returns

solution – See Solution for a description of the solution object. Additional fields specific to this solver are included in the object:

• onumpy.ndarray, shape (n,)

  Final solver mesh, in terms of the Boltzmann variable o.

• niterint

  Number of iterations required to find the solution.

• rms_residualsnumpy.ndarray, shape (n-1,)

  RMS values of the relative residuals over each mesh interval.

Return type

Solution

See also

solve()

Notes

Given that the location of the boundary is expressed in terms of the Boltzmann variable, a fixed boundary can be had only if o_guess[0] is 0. Any other o_guess[0] implies a moving boundary. This restriction affects radial problems in particular.

This function works by transforming the partial differential equation with the Boltzmann transformation using ode and then solving the resulting ODE with SciPy’s collocation-based boundary value problem solver scipy.integrate.solve_bvp and a two-point Dirichlet condition that matches the boundary and initial conditions of the problem. Upon that solver’s convergence, it runs a final check on whether the candidate solution also satisfies the semi-infinite condition (which implies \(dS/do\to0\) as \(o\to\infty\)).